Introduction to Parallel Processing in R using Futureverse - Easier than Ever Before

Henrik Bengtsson

University of California, San Francisco

R Foundation, R Consortium, Bioconductor

@HenrikBengtsson



useR! 2026, Warsaw, Poland (July 6, 2026)

About me

  • Used R since September 2000
  • MSc Computer Science, PhD Mathematical Statistics
  • Applied research in genomics and bioinformatics
  • Always had a need to speed up analyses and share worldwide
  • Creator of the future ecosystem (“futureverse”) - anno 2015
  • Best thing in R: The R community! People nerding out on R

How about you?

Q1. How many R conferences have you attended before?

  • my first!
  • my second
  • three or more

Q2. For how long have you used R?

  • < 1 year
  • 1-2 years
  • 3-5 years
  • > 5 years

Q3. Have you used parallelization in R before?

  • never
  • tried, but failed
  • yes, some
  • yes, a lot

Q4. What R parallelization methods have you used?

Q5. Reason for attending this parallelization tutorial?

Q6. What operating system do you use today?

  • Windows
  • macOS
  • Linux

Q7. How do you run R today?

  • in the terminal
  • in RStudio
  • in Positron
  • in Rgui
  • …?

Why do we parallelize code?

Parallel & distributed processing can be used to:

  • speed up processing (wall time)
  • lower memory footprint (per machine)
  • avoid data transfers (compute where data lives)
  • other reasons, e.g. asynchronous UI

How do we parallelize code?

We may choose to parallelize on:

  • our personal laptop or work desktop computer (single user)
  • a shared powerful computer (multiple users)
  • across many computers, e.g. in the office or in the cloud
  • high-performance compute (HPC) cluster (multiple users) with a job scheduler, e.g. Slurm, Son of Grid Engine (SGE), TORQUE/PBS, …

Your turn: Install and load ‘futureverse’

Step 1: Install packages

install.packages("futureverse")
install.packages("progress")
install.packages("beepr")        # optional

Step 2: Make sure everything is up-to-date

futureverse::futureverse_update()
#> All futureverse packages up-to-date

Step 3: Verify it works

library(futurize)
plan(multisession, workers = 2)
plan(sequential)

Your turn: Measuring processing times

We can measure the processing time of our code using system.time().

Define a slow function

slow_sqrt <- function(x) {
  Sys.sleep(0.5)  # half-second delay per item
  sqrt(x)
}

Running times for different input lengths

# ~1.5 seconds
system.time({ ys <- lapply(1:3, slow_sqrt) })

# ~4.0 seconds
system.time({ ys <- lapply(1:8, slow_sqrt) })

Your turn: Measuring processing times using in-house toc()

We can measure the processing time of our code using tic()/toc() from the futureverse package:

tic <- futureverse:::tic; toc <- futureverse:::toc

tic()
ys <- toc( lapply(1:3, slow_sqrt) )
ys <- toc( lapply(1:8, slow_sqrt) )

Or, equivalently using pipes:

tic()
ys <- lapply(1:3, slow_sqrt) |> toc()
ys <- lapply(1:8, slow_sqrt) |> toc()

Measure running times via progress reporting






All you need to remember is progressify()

Your turn: Introducing progressify

We can get a progress bar by piping our calls to progressify().

library(progressify)
handlers(global = TRUE)  # Must be enabled manually!

ys <- lapply(1:10, slow_sqrt) |> progressify()

Default output (using the built-in utils::txtProgressBar()):

> ys <- lapply(1:10, slow_sqrt) |> progressify()
|=============================        | 80%

Progress bar is erased on completion.

Your turn: Styled progress reporting

We can customize the style of the progress bar using the cli package.

library(progressify)
handlers("cli", global = TRUE)

ys <- lapply(1:10, slow_sqrt) |> progressify()

Output:

> ys <- lapply(1:10, slow_sqrt) |> progressify()
■■■■■■■■■■■■■■■■■■■■■■■■■         80% | ETA:  1s

Progress bar is erased on completion.

Your turn: Make some noise!

We can combine multiple progress handlers:

library(progressify)
handlers("cli", "beepr", global = TRUE)

ys <- lapply(1:10, slow_sqrt) |> progressify()

Output:

> ys <- lapply(1:10, slow_sqrt) |> progressify()
■■■■■■■■■■■■■■■■■■■■■■■■■         80% | ETA:  1s

with ♫♪♬♫ … ♬♫

Progress bar is erased on completion.

Your turn: Keep the progress output when done

By default, progress bars are erased on completion. This feature can be disabled, but to keep things simple, we will use a special “tutorial” handler today:

library(progressify)
handlers("tutorial", global = TRUE)

ys <- lapply(1:10, slow_sqrt) |> progressify()

Output:

> ys <- lapply(1:10, slow_sqrt) |> progressify()
[ 5s] | [============================] 100% [10/10] (ETA:  0s)
>

Progress bar remains beyond completion.

Your turn: Keep all progress updates

By default, the progress bars are cleared at each update, and also at the very end. Because of this, we will use a special “tutorial” handler for this workshop:

library(progressify)
handlers("tutorial", "newline", global = TRUE)

ys <- lapply(1:10, slow_sqrt) |> progressify()

Output:

> ys <- lapply(1:10, slow_sqrt) |> progressify()
[ 0s] - [-----------------------------]   0% [0/10] (ETA:  ?s)
[ 1s] \ [==>--------------------------]  10% [1/10] (ETA:  5s) 
[ 1s] | [=====>-----------------------]  20% [2/10] (ETA:  4s) 
[ 2s] / [========>--------------------]  30% [3/10] (ETA:  4s) 
[ 2s] - [===========>-----------------]  40% [4/10] (ETA:  3s) 
[ 3s] \ [==============>--------------]  50% [5/10] (ETA:  3s) 
[ 3s] | [=================>-----------]  60% [6/10] (ETA:  2s) 
[ 4s] / [====================>--------]  70% [7/10] (ETA:  2s) 
[ 4s] - [=======================>-----]  80% [8/10] (ETA:  1s) 
[ 5s] \ [==========================>--]  90% [9/10] (ETA:  1s) 
[ 5s] | [============================] 100% [10/10] (ETA:  0s)
>

Progress output remains beyond completion.

R: A Functional Programming Language

R is a functional programming language

Mathematics:

\[ \begin{aligned} f(x) &= 2 \cdot x \\ y &= f(1) \\ z &= f(3) \end{aligned} \]

R:

f <- function(x) 2*x
y <- f(1)
z <- f(3)

A function has code, arguments, and a value

f <- function(x) {
  2*x
}

A function has:

  • code, e.g. 2*x
  • arguments (“parameters”), e.g. x
  • a value

When we call the function:

y <- f(4)

R evaluates the code with the value of argument x (= 4), returns the value (= 8), and assigns the value to variable y.

A function can take multiple arguments

f <- function(x, b) {
  b*x
}

A function has:

  • code, e.g. b*x
  • arguments (“parameters”), e.g. x and b
  • a value

When we call the function:

y <- f(4, 2)

R evaluates the code with the values of arguments x (= 4) and b (= 2), returns the value, and assigns the value to variable y.

Arguments can have default values

f <- function(x, b = 2) {
  b*x
}

When we call the function:

y <- f(4)

R evaluates the code with the value of argument x (= 4), the default value of b (= 2), returns the value, and assigns the value to variable y.

Summation

\[ h(\mathbf{x}) = \sum_{i=1}^{n} x_i ; \mathbf{x} = (x_1, x_2, \ldots, x_n) \\ \]

h <- function(x) {
  sum <- 0
  n <- length(x)
  for (i in 1:n) sum <- sum + x[i]
  sum
}

Summation - emulate slowness

\[ h(\mathbf{x}) = \sum_{i=1}^{n} x_i ; \mathbf{x} = (x_1, x_2, \ldots, x_n) \\ \]

h <- function(x) {
  sum <- 0
  for (i in seq_along(x)) {
    Sys.sleep(1.0)   # one-second delay
    sum <- sum + x[i]
  }
  sum
}
h(1:3)  # ~3 seconds
#> [1] 6

Summation can be done in parts

\[ h(\mathbf{x}) = \sum_{i=1}^{n} x_i \]

\[ \phantom{h(\mathbf{x})} = \sum_{i=1}^{k} x_i + \sum_{i=k+1}^{n} x_i \\ \]

\[ \phantom{h(\mathbf{x})} = h(\mathbf{x}_{1:k}) + h(\mathbf{x}_{(k+1):n}) \]

Summation can be done in parts

\[ \begin{aligned} h(\mathbf{x}) &= \sum_{i=1}^{n} x_i \\ \mathbf{x} &= (1, 2, \ldots, 10) \end{aligned} \]

\[ y = h(\mathbf{x}) \]

\[ \begin{aligned} a &= h(\mathbf{x}_{1:5}) \\ b &= h(\mathbf{x}_{6:10}) \\ y &= a + b \end{aligned} \]

\[ y = 55 \]

Summation of parts can be done in any order

\[ ~\\ ~\\ y = h(\mathbf{x}) \]

y <- h(x)

\[ \begin{aligned} a &= h(\mathbf{x}_{1:5}) \\ b &= h(\mathbf{x}_{6:10}) \\ y &= a + b \end{aligned} \]

a <- h(x[1:5])
b <- h(x[6:10])
y <- a + b

\[ \begin{aligned} b &= h(\mathbf{x}_{6:10}) \\ a &= h(\mathbf{x}_{1:5}) \\ y &= a + b \end{aligned} \]

b <- h(x[6:10])
a <- h(x[1:5])
y <- a + b

\[ y = 55 \]

All takes ~10 seconds

Concurrent summation of parts

\[ \text{concurrently} \left\{\begin{matrix} a = h(\mathbf{x}_{1:5}) \\ b = h(\mathbf{x}_{6:10}) \\ \end{matrix}\right. \\ \hspace{10ex} y = a + b \]

How can we calculate \(a\) and \(b\) concurrently in R?

TL;DR

a %<-% h(x[1:5])
b %<-% h(x[6:10])
y <- a + b

a and b take 5 seconds each => everything takes ~5 seconds

Package ‘future’ adds support for concurrent programming in R


library(future)
plan(multisession)

a %<-% h(x[1:5])
b %<-% h(x[6:10])

y <- a + b


  • %<-% is called the future-assignment operator

Your turn: Setup

Step 3: Implement \(h()\) as above, but call it slow_sum()

slow_sum <- function(x) {
  sum <- 0
  for (value in x) {
    Sys.sleep(1.0)    # one-second slowdown per value
    sum <- sum + value
  }
  sum
}

Your turn: Verify

Step 4: Test slow_sum()

  1. Create \(x = (1, 2, ..., 10)\) in R

  2. Call y <- slow_sum(x) and confirm it takes 10 seconds

  3. Call a <- slow_sum(x[1:5]) and confirm it takes 5 seconds

  4. Same for b <- slow_sum(x[6:10])

Your turn: Our very first concurrent evaluation of R code

Step 5: Sum parts concurrently

library(future)
plan(multisession)

x <- 1:10

a %<-% slow_sum(x[1:5])
b %<-% slow_sum(x[6:10])

y <- a + b

  1. Confirm correct result

  2. Confirm that it finishes in 5 seconds - use tic()/toc()

  3. How long does each line take? - use tic()/toc()

  4. Something else?

Note on y <- a + b |> toc()

Note that:

y <- a + b |> toc()

is the same as:

y <- a + toc(b)

which is not what we wanted, i.e.

y <- toc(a + b)

To avoid this, we need to use:

y <- { a + b } |> toc()

which is identical to:

y <- toc({ a + b })

Concurrent Programming in R

Future … what?

  • A future is an abstraction for a value that will be available later
  • The state of a future is either unresolved or resolved
  • The value is the result of an evaluated expression

An R assignment:

v <- expr

A future-value assignment:

f <- future(expr)
v <- value(f)


y <- slow_sum(x)
f <- future({ slow_sum(x) })
y <- value(f)

References

Friedman & Wise (1976, 1977), Hibbard (1976), Baker & Hewitt (1977)

Your turn: future() and value()

library(future)
plan(multisession)

a %<-% slow_sum(x[1:5])
b %<-% slow_sum(x[6:10])

y <- a + b
library(future)
plan(multisession)

fa <- future( slow_sum(x[1:5])  )
fb <- future( slow_sum(x[6:10]) )

y <- value(fa) + value(fb)

  • Verify result

  • Verify total run time

  • Verify run time per line

Functions

\[ \begin{aligned} a &= \text{slow_sum}(\mathbf{x}_{1:5}) \\ b &= \text{slow_sum}(\mathbf{x}_{6:10}) \\ y &= a + b \end{aligned} \]

a <- slow_sum(x[1:5])
b <- slow_sum(x[6:10])
y <- a + b

Functions and Futures

\[ \text{concurrently} \left\{\begin{matrix} a = \text{slow_sum}(\mathbf{x}_{1:5}) \\ b = \text{slow_sum}(\mathbf{x}_{6:10}) \\ \end{matrix}\right. \\ \hspace{10ex} y = a + b \]

fa <- future( slow_sum(x[1:5])  )
fb <- future( slow_sum(x[6:10]) )
y <- value(fa) + value(fb)

Futureverse - A Friendly, Unifying Parallelization Framework in R

  • “Write once, run anywhere” - a simple unified API
  • 100% cross-platform - easy to install
  • Works with any type of parallel backend
  • Very well tested, lots of CPU mileage
    11 years, used by 1,600 packages, top-0.6% downloaded

Low friction:

  • Automatically exports global variables
  • Automatically relays output, messages, and warnings
  • Proper parallel random number generation (RNG)
  • … and much more

Globals automatically found (99.5% worry-free)

Static-code inspection by walking the abstract syntax tree (AST):

x <- rnorm(n = 100)          lobstr::ast( { slow_sum(x) } )
f <- future({ slow_sum(x) })     |  █─`{`
                   └─────────────|  └─█─ slow_sum
                                 |    └─ x

=> globals identified and exported to the worker:

  • slow_sum() - a function (also searched recursively)
  • x - a numeric vector of length 100

Globals can also be specified manually:

f <- future({ slow_sum(x) }, globals = c("slow_sum", "x"))

Apply the same function over and over

For-loop:

xs <- list(1:5, 6:10, 11:15, 16:20)
ys <- list()
for (i in seq_along(xs)) {
  ys[[i]] <- slow_sum(xs[[i]])
}


ys <- list()
for (x in xs) {
  ys <- c(ys, slow_sum(x))
}

Apply the same function over and over

Base R:

xs <- list(1:5, 6:10, 11:15, 16:20)
ys <- lapply(xs, slow_sum)

Tidyverse:

library(purrr)
xs <- list(1:5, 6:10, 11:15, 16:20)
ys <- map(xs, slow_sum)

Apply the same function over and over

Plyr:

library(plyr)
xs <- list(1:5, 6:10, 11:15, 16:20)
ys <- llply(xs, slow_sum)

Foreach:

library(foreach)
xs <- list(1:5, 6:10, 11:15, 16:20)
ys <- foreach(x = xs) %do% slow_sum(x)

Take-home: lapply(), map(), llply(), foreach() all the same


ys <- lapply(xs, slow_sum)

ys <- purrr::map(xs, slow_sum)

ys <- plyr::llply(xs, slow_sum)

ys <- foreach(x = xs) %do% slow_sum(x)

Pipes: Syntactic sugar to spice it up

Base R:

ys <- lapply(xs, slow_sum)

is identical to

ys <- xs |> lapply(slow_sum)

R can’t tell the difference! (e.g. try with lobstr::ast(...))

Tidyverse:

ys <- xs |> map(slow_sum)

Plyr:

ys <- xs |> llply(slow_sum)

Your turn: Draw an Owl

ys <- xs |> map(slow_sqrt)
ys <- xs |> parallel_map(slow_sqrt)

Your turn: Start with a for-loop

Task: Create a function ys <- parallel_map(xs, slow_sqrt)

Combine the idea of:

f <- future( slow_sqrt(x) )
v <- value(f)

with:

ys <- list()
for (i in seq_along(xs)) {
  ys[[i]] <- slow_sqrt(xs[[i]])
}

to make a “concurrent” for-loop.

Hint: You need two for-loops - one for future() and one for value()

Your turn: Concurrent for-loop to concurrent map

# Create futures
fs <- list()
for (i in seq_along(xs)) {
  fs[[i]] <- future( slow_sqrt(xs[[i]]) )
}

# Collect values
ys <- list()
for (i in seq_along(xs)) {
  ys[[i]] <- value( fs[[i]] ) 
}

Coffee Break (30 mins)

Your turn: Create parallel_map()

Building blocks:

# Create futures
fs <- purrr::map(xs, function(x) { future(slow_sqrt(x)) })

# Collect values
ys <- purrr::map(fs, value)

value() is vectorized, i.e. it can operate on lists too

# Create futures
fs <- purrr::map(xs, function(x) { future(slow_sqrt(x)) })

# Collect values
ys <- value(fs)

Your turn: Benchmark map() vs parallel_map()

Compare how fast:

xs <- runif(30)
ys <- purrr::map(xs, slow_sqrt)

is compared to:

plan(multisession, workers = 4)
xs <- runif(30)
ys <- parallel_map(xs, slow_sqrt)
  • Use system.time() or tic()/toc()

Why is value(fs) preferred over lapply(fs, value)?

lapply(fs, value) resolves futures sequentially:

  • No early stopping: If the third future fails, lapply() still blocks waiting for the first and second futures to finish.
  • Orphaned workers: If a future fails, the remaining running futures are left running on the workers, wasting resources.

Vectorized value(fs) is aware of all futures:

  • Early stopping: If any future fails, times out, or its worker crashes, value() immediately throws the error.
  • Automatic cancellation: Outstanding futures are sent a termination signal, instantly freeing up the workers.

ys <- parallel_map(xs, slow_sqrt)

Good job! You’ve learned a lot!



  • The concept of “futures” was invented in 1975 (sic!)
  • Basic building blocks future() and value() are very powerful
  • You can create richer parallelization functions from these, e.g. parallel_map()

Now, let’s simplify our lives …






All you need to remember is futurize()

Futureverse allows you to stick to your favorite coding style

Parallel alternatives to traditional, sequential functions:

ys <- lapply(xs, slow_sqrt)                           # {base}
ys <- lapply(xs, slow_sqrt) |> futurize()             # {futurize}

ys <- map(xs, slow_sqrt)                              # {purrr}
ys <- map(xs, slow_sqrt) |> futurize()                # {futurize}

ys <- foreach(x = xs) %do% slow_sqrt(x)               # {foreach}
ys <- foreach(x = xs) %do% slow_sqrt(x) |> futurize() # {futurize}


future.apply and furrr?

  • futurize() is new since January 2026
  • future.apply, furrr, and doFuture are now only used in the background

Futures can be resolved anywhere!

Our ys <- map(xs, slow_sqrt) |> futurize() can be resolved on lots of different backends - it doesn’t matter where!

On the current machine

plan(sequential)

sequentially (default)

plan(multisession)
plan(multisession, workers=4)

in parallel via background R workers

plan(future.mirai::mirai_multisession)
plan(future.mirai::mirai_multisession, workers=4)

in parallel via future.mirai; low latency

plan(future.callr::callr)
plan(future.callr::callr, workers=4)

in parallel via future.callr; all memory is returned as each future is resolved

Futures can be resolved anywhere!

On other machines

plan(cluster, workers=c("n1", "m2.uni.edu", "vm.cloud.org"))

distributed on other local or remote computers

plan(future.mirai::mirai)

distributed via future.mirai daemon queue

On high-performance compute (HPC) cluster

In parallel on HPC job schedulers for long-running tasks with higher latency.

plan(future.batchtools::batchtools_lsf)
plan(future.batchtools::batchtools_openlava)
plan(future.batchtools::batchtools_sge)
plan(future.batchtools::batchtools_slurm)
plan(future.batchtools::batchtools_torque)

HPC processing: future.batchtools


plan(future.batchtools::batchtools_slurm)                          ‎

fq <- fs::dir_ls(glob = "*.fq")             # 80 files; 200 GB each
bam <- purrr::map(fq, align) |> futurize()  # 1 hour each


{henrik: ~}$ squeue
Job ID   Name            User         Time Use S
-------- --------------- ------------ -------- -
606411   xray            alice        46:22:22 R
606638   fz:purrr-5      henrik       00:52:05 R
606641   python          bob          37:18:30 R
606643   fz:purrr-6      henrik       00:51:55 R
...

Making progress in parallel . . .




 
All you need to remember is futurize() + progressify()

Progress reporting + Parallel processing

You can combine progressify() and futurize() to get progress reporting in parallel:

library(futurize)

## Set up three local workers and one remote
plan(cluster, workers = c("m1", "m2", "m3", "server.uni.edu"))

library(progressify)
handlers("cli", global = TRUE)

xs <- runif(50)
ys <- lapply(xs, slow_sqrt) |> progressify() |> futurize()

Output:

ys <- lapply(xs, slow_sqrt) |> progressify() |> futurize()
■■■■■■■■■■■■■■■■■■■■■■■■■         80% | ETA:  1s

Progress is reported in a near-live fashion.

Your turn: Try different backends

Sequential:

library(futurize)
plan(sequential)  # default

xs <- runif(20)
ys <- purrr::map(xs, slow_sqrt) |> futurize()

Built-in multisession:

plan(multisession)
ys <- purrr::map(xs, slow_sqrt) |> futurize()

Add-on mirai-based multisession:

plan(future.mirai::mirai_multisession)
ys <- purrr::map(xs, slow_sqrt) |> futurize()

Your turn: Benchmark with different number of workers

Run:

library(futurize)
library(progressify)
handlers("tutorial", global = TRUE)

n <- 20
xs <- runif(n)
ys <- purrr::map(xs, slow_sqrt) |> progressify() |> futurize()

Then benchmark with:

plan(multisession, workers = 2)
plan(multisession, workers = 4)
plan(multisession, workers = 6)
plan(multisession, workers = 8)
  • Try to call plan() without arguments

Your turn: Benchmark different backends

Run:

library(futurize)
library(progressify)
handlers("tutorial", global = TRUE)

n <- 20
xs <- runif(n)
ys <- purrr::map(xs, slow_sqrt) |> progressify() |> futurize()

Then benchmark with:

plan(multicore)
plan(multisession)
plan(future.mirai::mirai_multisession)
plan(future.callr::callr)

Messages, Warnings, and Errors




     

Some functions have side effects

g <- function(x, b = 2) {
  if (b <= 0) stop("Argument 'b' must be a positive scalar: ", b)
  sqrt(b * x)
}

What type of results may this function produce?

  • a value
> g(1)
[1] 2
  • an error
> g(1, b = 0)
Error in g(1, b = 0) :
Argument 'b' must be a positive scalar: 0
  • a value and a warning
> g(-1)
[1] NaN
Warning message:
In sqrt(b * x) : NaNs produced

Your turn: Messages, warnings, and errors with futures

g <- function(x, b = 2) {
  message("x = ", x)
  if (b <= 0) stop("Argument 'b' must be a positive scalar: ", b)
  Sys.sleep(3.0)
  sqrt(b * x)
}

f <- future({ g(1) })
v <- value(f)

f <- future({ g(1, b = 0) })
v <- value(f)

f <- future({ g(-1) })
v <- value(f)

Your turn: Errors terminate calls

xs <- c(1, 2, 3, 4)
bs <- c(4, 3, 2, 1)
ys <- purrr::map2_dbl(xs, bs, g)
g <- function(x, b = 2) {
  message("x = ", x)
  if (b <= 0)
    stop("Argument 'b' must be positive: ", b)
  Sys.sleep(0.5)
  sqrt(b * x)
}

xs <- c(1, 2, -3, 4)
bs <- c(4, 3,  2, 1)
ys <- purrr::map2_dbl(xs, bs, g)

xs <- c(1, 2,  3, 4)
bs <- c(4, 3, -2, 1)
ys <- purrr::map2_dbl(xs, bs, g)

Your turn: Errors terminate parallelized calls

library(futurize)
plan(multisession, workers = 4)

xs <- c(1, 2, 3, 4)
bs <- c(4, 3, 2, 1)
ys <- purrr::map2_dbl(xs, bs, g) |> futurize()

xs <- c(1, 2, -3, 4)
bs <- c(4, 3,  2, 1)
ys <- purrr::map2_dbl(xs, bs, g) |> futurize()

xs <- c(1, 2,  3, 4)
bs <- c(4, 3, -2, 1)
ys <- purrr::map2_dbl(xs, bs, g) |> futurize()

Business as usual but in parallel

library(futurize)
plan(multisession, workers = 4)
xs <- c(1, 2, -3, 4)
bs <- c(4, 3, 2, 1)
ys <- purrr::map2_dbl(xs, bs, g) |> futurize()

> ys <- purrr::map2_dbl(xs, bs, g) |> futurize()
x = 1
x = 2
x = -3
x = 4
Warning message:
In sqrt(b * x) : NaNs produced
> ys
[1] 2.00000 2.44949     NaN 2.00000

Things we haven’t talked about

Random number generation (RNG) in parallel processing

Random numbers need extra care when tasks run in parallel.

Using rnorm() without declaring it:

xs <- 1:4
ys <- lapply(xs, function(i) rnorm(1)) |> futurize()
Warning message:
UNRELIABLE VALUE: Future (...) unexpectedly generated random numbers
without specifying argument 'seed'. There is a risk that those random
numbers are not statistically sound and the overall results might be
invalid. To fix this, specify 'seed = TRUE'. ...

Declaring parallel-safe RNG with seed = TRUE:

ys <- lapply(xs, function(i) rnorm(1)) |> futurize(seed = TRUE)

Why it matters:

  • Workers can share or overlap the same random stream => duplicated or correlated draws
  • Results may change with the number of workers and the backend => not reproducible
  • seed = TRUE uses the parallel-safe L’Ecuyer-CMRG generator to give each task its own independent, reproducible stream

When parallelization does not help

Parallelization adds overhead - it is not always a win.

xs <- 1:10000
ys <- purrr::map(xs, sqrt)                 # ~instant
ys <- purrr::map(xs, sqrt) |> futurize()   # much slower!

Parallelization can be slower when:

  • Tasks are too small - per-task overhead (launch, serialize, transfer) exceeds the actual work
  • Data is large - copying inputs and outputs to and from workers dominates
  • The bottleneck is not the CPU - disk-, network-, or memory-bandwidth-bound work does not speed up
  • The code is already vectorized - e.g. sqrt(xs) beats any parallel map

Rules of thumb

  • Parallelize the slow, coarse chunks, not the tiny fast ones
  • Group many small tasks into fewer larger ones (chunking)
  • Measure - use tic()/toc(), don’t assume

Not everything can travel to a worker

Some objects are tied to your current R session and cannot be exported to a worker:

  • Database connections (DBI), open file or socket connections
  • External pointers - xgboost, torch, keras/tensorflow, Rcpp
  • reticulate Python objects
con <- DBI::dbConnect(RSQLite::SQLite(), "my.db")   # tied to this session
ys <- purrr::map(qs, function(q) DBI::dbGetQuery(con, q)) |> futurize()
Error: Detected a non-exportable reference in the globals:
'con' of class 'SQLiteConnection' uses an external pointer

Enable early detection with options(future.globals.onReference = "error").

Fix: create the resource inside each task, on the worker:

ys <- purrr::map(qs, function(q) {
  con <- DBI::dbConnect(RSQLite::SQLite(), "my.db")
  on.exit(DBI::dbDisconnect(con))
  DBI::dbGetQuery(con, q)
}) |> futurize()

Parallelizing select domain-specific calls

library(futurize); plan(multisession, workers = 8)
library(progressify); handlers("tutorial")

# Simulate dataset
n <- 5000
data <- data.frame(y = rbinom(n, 1, 0.5), x1 = rnorm(n), x2 = rnorm(n))

# Fit logistic regression
fit_glm <- function(data, indices) {
  coef(glm(y ~ x1 + x2, data = data[indices, ], family = binomial))
}

# Bootstrap
b <- boot::boot(data, statistic = fit_glm, R = 1000) |> progressify() |> futurize()
  • I’m giving a talk on futurize on Wednesday at 14:00

The same futurize() works across your stats packages

Popular statistics & ML packages support futurize() - no new API to learn:

# Bootstrap ({boot})
boot::boot(data, statistic, R = 1000)         |> futurize()

# Cross-validation ({glmnet})
glmnet::cv.glmnet(x, y)                       |> futurize()

# Model tuning ({caret})
caret::train(y ~ ., data = d, method = "rf")  |> futurize()

# Additive models ({mgcv})
mgcv::bam(y ~ s(x1) + s(x2), data = d)        |> futurize()

# Mixed-model bootstrap ({lme4})
lme4::bootMer(fit, statistic, nsim = 1000)    |> futurize()

Your existing call, plus one pipe - it runs on whatever backend plan() you chose.

Under the hood: futurize() is a transpiler

futurize() doesn’t parallelize your call directly - it rewrites (“transpiles”) it into the equivalent Futureverse call, then evaluates that.

Pass eval = FALSE to see the generated code instead of running it:

lapply(xs, fcn) |> futurize(eval = FALSE)
future.apply::future_lapply(xs, fcn, future.seed = FALSE, future.globals = TRUE,
    future.stdout = TRUE, future.conditions = "condition", future.scheduling = 1,
    future.label = "fz:base::lapply-%d")

The same holds for domain-specific packages - inspect what boot generates:

boot::boot(data, statistic, R = 1000) |> futurize(eval = FALSE)
local({
    cl <- do.call(future::makeClusterFuture, args = list(seed = FALSE,
        globals = TRUE, stdout = TRUE, conditions = "condition",
        label = "fz:boot::boot"))
    oopts <- options(future.ClusterFuture.clusterEvalQ = "ignore")
    on.exit(options(oopts))
    boot::boot(data, statistic, R = 1000, parallel = "snow",
        ncpus = 2L, cl = cl)
})

Each supported package has its own transpiler mapping the sequential call to a parallel backend.

Nested parallelization with some care

  • “Structured concurrency” -
    Futureverse protects against over-parallelization

Sequential
processing

Parallelization
with 4 workers

Nested parallelization
with 5 workers each
running 3 parallel tasks

Cheat sheet (1/2): choose a backend & parallelize

library(futurize)

## 1. Choose WHERE to run (once)
plan(sequential)                        # default - no parallelism
plan(multisession, workers = 4)         # background R workers
plan(cluster, workers = c("n1", "n2"))  # other machines

## 2. Parallelize - just add a pipe
ys <- lapply(xs, fcn)             |> futurize()
ys <- purrr::map(xs, fcn)         |> futurize()
ys <- foreach(x = xs) %do% fcn(x) |> futurize()

Cheat sheet (2/2): progress, RNG & building blocks

## 3. Add progress reporting
library(progressify)
handlers("cli", global = TRUE)
ys <- purrr::map(xs, fcn) |> progressify() |> futurize()

## 4. Reproducible random numbers
ys <- purrr::map(xs, fcn) |> futurize(seed = TRUE)

## Under the hood: the building blocks
f <- future(fcn(x))
v <- value(f)

The future is bright!

Coming: Efficient memory handling

  • Shared memory support is on the horizon:
    • Powered by the upcoming mori package
    • Credit to the old Bioconductor SharedObject package (which is very similar)
  • Shared memory only works on parallel workers running on the same machine
  • I aim to make it fully automatic, so you will not have to do anything
  • When released, all existing code will automatically see an improvement

Coming: Resource management

  • Declarative Resource API:
ys <- lapply(xs, fit) |> 
        resources(memory = 2*GiB, gpu = TRUE, time = 5*hours) |> futurize()

  • Avoid Resource Waste:

    • Conservative: prevent impossible runs - early stopping
    • Scheduling: avoid running out of memory, HPC job scheduling, …

  • Monitor Resources:

    • Track and monitor resource usage (e.g. CPU, GPU, memory, runtime) during execution
  • Predict Resources: (possibly also)
    • Model resource needs over time, e.g. runtime and memory
    • Predict resource needs on the input, e.g. suggest specifications, better progress ETA

Thank you - Enjoy useR! 💜

Appendix

foreach() is not a for loop!

We must never think of foreach() as a for-loop

The name foreach() tricks us into believing it’s a for-loop



fq <- fs::dir_ls(glob = "*.fq")
bam <- list()
for (ii in seq_along(fq)) {
  bam[[ii]] <- align(fq[[ii]])
}



   !=

library(foreach)
fq <- fs::dir_ls(glob = "*.fq")
bam <- list()
foreach(ii = seq_along(fq)) %dopar% {
  bam[[ii]] <- align(fq[[ii]])
}

This does not work because … foreach() is not a for-loop!

> str(bam)
List of 80
 $ : chr "sample_01.bam"
 $ : chr "sample_02.bam"
 ...
 $ : chr "sample_80.bam"
> 
> str(bam)
 list()
>

Repeat after me: foreach() is not a for-loop! 🎶


for (ii in 1:1000) {
  message("foreach() is not a for-loop!")
  Sys.sleep(2)
  beepr::beep("wilhelm")
}

foreach() is just like lapply() …


fq <- fs::dir_ls(glob = "*.fq")
bam <- list()
lapply(seq_along(fq), function(ii) {
  bam[[ii]] <- align(fq[[ii]])
})

This does not work, because the <- assignment is done inside a function

… and just like map() …


fq <- fs::dir_ls(glob = "*.fq")
bam <- list()
purrr::map(seq_along(fq), function(ii) {
  bam[[ii]] <- align(fq[[ii]])
})

This does not work, because the <- assignment is done inside a function

Take-home: lapply(), map(), foreach() return values


fq <- fs::dir_ls(glob = "*.fq")

bam <- lapply(fq, align)

bam <- purrr::map(fq, align)

bam <- foreach(x = fq) %do% align(x)