%\VignetteIndexEntry{1. Introduction to BiocParallel}
%\VignetteKeywords{parallel, Infrastructure}
%\VignettePackage{BiocParallel}
%\VignetteEngine{knitr::knitr}


\documentclass{article}

<<style, eval=TRUE, echo=FALSE, results="asis">>=
BiocStyle::latex()
@

<<setup, echo=FALSE>>=
suppressPackageStartupMessages({
    library(BiocParallel)
    library(VariantAnnotation)
    library(GenomicAlignments)
    library(RNAseqData.HNRNPC.bam.chr14)
    library(TxDb.Hsapiens.UCSC.hg19.knownGene)
})
@

\newcommand{\BiocParallel}{\Biocpkg{BiocParallel}}

\title{Introduction to \BiocParallel}
\author{Valerie Obenchain, Vincent Carey, Michael Lawrence,
  Martin Morgan\footnote{\url{Martin.Morgan@RoswellPark.org}}}
\date{Edited: January 23, 2021; Compiled: \today}

\begin{document}

\maketitle

\tableofcontents

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Introduction}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Numerous approaches are available for parallel computing in \R{}. The
CRAN Task View for high performance and parallel computing provides
useful high-level summaries and package categorization.
\url{http://cran.r-project.org/web/views/HighPerformanceComputing.html}
Most Task View packages cite or identify one or more of \CRANpkg{snow},
\CRANpkg{Rmpi}, \CRANpkg{multicore} or \CRANpkg{foreach} as relevant
parallelization infrastructure. Direct support in \R{} for parallel
computing started with release 2.14.0 with inclusion of the
\CRANpkg{parallel} package which contains modified versions of
\CRANpkg{multicore} and \CRANpkg{snow}.

A basic objective of \BiocParallel{} is to reduce the complexity faced
when developing and using software that performs parallel computations.
With the introduction of the \Rcode{BiocParallelParam} object,
\BiocParallel{} aims to provide a unified interface to existing parallel
infrastructure where code can be easily executed in different
environments. The \Rcode{BiocParallelParam} specifies the environment of
choice as well as computing resources and is invoked by `registration` or
passed as an argument to the \BiocParallel{} functions.

\BiocParallel{} offers the following conveniences over the
`roll your own` approach to parallel programming.

\begin{itemize}
  \setlength{\itemsep}{5pt}
  \item{unified interface:}{
    \Rcode{BiocParallelParam} instances define the method of parallel
    evaluation (multi-core, snow cluster, etc.) and computing resources
    (number of workers, error handling, cleanup, etc.).
  }
  \item{parallel iteration over lists, files and vectorized operations:}{
    \Rcode{bplapply}, \Rcode{bpmapply} and \Rcode{bpvec} provide
    parallel list iteration and vectorized operations. \Rcode{bpiterate}
    iterates through files distributing chunks to parallel workers.
  }
  \item{cluster scheduling:}{
    When the parallel environment is managed by a cluster scheduler
    through \CRANpkg{batchtools}, job management and result retrieval are
    considerably simplified.
  }
  \item{support of \Rcode{foreach}:}{
    The \CRANpkg{foreach} and \CRANpkg{iterators} packages are fully supported.
    Registration of the parallel back end uses \Rcode{BiocParallelParam}
    instances.
  }
\end{itemize}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Quick start}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

The \Rpackage{BiocParallel} package is available at bioconductor.org
and can be downloaded via \Rcode{BiocManager}:
<<BiocManager, eval=FALSE>>=

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("BiocParallel")
@

Load \BiocParallel{}.
<<BiocParallel>>=
library(BiocParallel)
@

The test function simply returns the square root of ``x''.
<<quickstart_FUN>>=
FUN <- function(x) { round(sqrt(x), 4) }
@

Functions in \BiocParallel use the registered back-ends for parallel
evaluation. The default is the top entry of the registry list.
<<quickstart_registry>>=
registered()
@
%%
Configure your R session to always use a particular back-end configure
by setting options named after the back ends in an \Rcode{.Rprofile}
file, e.g.,
<<configure_registry, eval=FALSE>>=
options(MulticoreParam=quote(MulticoreParam(workers=4)))
@

When a \BiocParallel{} function is invoked with no \Rcode{BPPARAM}
argument the default back-end is used.
<<quickstart_bplapply_default, eval=FALSE>>=
bplapply(1:4, FUN)
@

Environment specific back-ends can be defined for any of the registry entries.
This example uses a 2-worker SOCK cluster.
<<quickstart_snow>>=
param <- SnowParam(workers = 2, type = "SOCK")
bplapply(1:4, FUN, BPPARAM = param)
@

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{The \BiocParallel{} Interface}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\subsection{Classes}

\subsubsection{\Rcode{BiocParallelParam}}

\Rcode{BiocParallelParam} instances configure different parallel evaluation
environments. Creating or \Rcode{register()}ing a `\Rcode{Param}`
allows the same code to be used in different parallel environments without
a code re-write. Params listed are supported on all of Unix, Mac and Windows
except \Rcode{MulticoreParam} which is Unix and Mac only.

\begin{itemize}
  \setlength{\itemsep}{5pt}
  \item{\Rcode{SerialParam}: }

  Supported on all platforms.

  Evaluate \BiocParallel-enabled code with parallel evaluation disabled.
  This approach is useful when writing new scripts and trying to debug code.

  \item{\Rcode{MulticoreParam}: }

  Supported on Unix and Mac. On Windows, \Rcode{MulticoreParam} dispatches
  to \Rcode{SerialParam}.

  Evaluate \BiocParallel-enabled code using multiple cores on a single
  computer. When available, this is the most efficient and least
  troublesome way to parallelize code. Windows does not support
  multi-core evaluation (the \Rcode{MulticoreParam} object can be
  used, but evaluation is serial). On other operating systems, the
  default number of workers equals the value of the global option
  \Rcode{mc.cores} (e.g., \Rcode{getOption("mc.cores")}) or, if that
  is not set, the number of cores returned by
  \Rcode{parallel::detectCores() - 2}; when number of cores cannot be
  determined, the default is 1.

  \Rcode{MulticoreParam} uses 'forked' processes with 'copy-on-change'
  semantics -- memory is only copied when it is changed. This makes it
  very efficient to invoke compared to other back-ends.

  There are several important caveats to using
  \Rcode{MulticoreParam}. Forked processes are not available on
  Windows. Some environments, e.g., \emph{RStudio}, do not work well
  with forked processes, assuming that \R{} code evaluation is
  single-threaded. Some external resources, e.g., access to files or
  data bases, maintain state in a way that assumes the resource is
  accessed only by a single thread.  A subtle cost is that \R's
  garbage collector runs periodically, and 'marks' memory as in
  use. This effectively triggers a copy of the marked memory. \R's
  generational garbage collector is triggered at difficult-to-predict
  times; the effect in a long-running forked process is that the
  memory is eventually copied. See
  \href{https://support.bioconductor.org/p/70196/#70509}{this post}
  for additional details.

  \Rcode{MulticoreParam} is based on facilities originally implemented
  in the \CRANpkg{multicore} package and subsequently the
  \CRANpkg{parallel} package in base \R{}.

  \item{\Rcode{SnowParam}: }

  Supported on all platforms.

  Evaluate \BiocParallel-enabled code across several distinct \R{}
  instances, on one or several computers. This is a straightforward approach
  for executing parallel code on one or several computers, and is
  based on facilities originally implemented in the \CRANpkg{snow} package.
  Different types of \CRANpkg{snow} `back-ends' are supported, including
  socket and MPI clusters.

  \item{\Rcode{BatchtoolsParam}: }

  Applicable to clusters with formal schedulers.

  Evaluate \BiocParallel-enabled code by submitting to a cluster scheduler
  like SGE.

  \item{\Rcode{DoparParam}: }

  Supported on all platforms.

  Register a parallel back-end supported by the \CRANpkg{foreach} package
  for use with \BiocParallel.
\end{itemize}

The simplest illustration of creating \Rcode{BiocParallelParam} is
<<BiocParallelParam_SerialParam>>=
serialParam <- SerialParam()
serialParam
@

Most parameters have additional arguments influencing behavior, e.g.,
specifying the number of `cores' to use when creating a
\Rcode{MulticoreParam} instance
<<BiocParallelParam_MulticoreParam>>=
multicoreParam <- MulticoreParam(workers = 8)
multicoreParam
@

Arguments are described on the corresponding help page, e.g.,
\Rcode{?MulticoreParam}.

\subsubsection{\Rcode{register()}ing \Rcode{BiocParallelParam} instances}

The list of registered \Rcode{BiocParallelParam} instances represents
the user's preferences for different types of back-ends. Individual
algorithms may specify a preferred back-end, and different back-ends
maybe chosen when parallel evaluation is nested.

The registry behaves like a `stack' in that the last entry registered
is added to the top of the list and becomes the ``next used``
(i.e., the default).

\Rcode{registered} invoked with no arguments lists all back-ends.
<<register_registered>>=
registered()
@

\Rcode{bpparam} returns the default from the top of the list.
<<register_bpparam>>=
bpparam()
@

Add a specialized instance with \Rcode{register}. When
\Rcode{default} is TRUE, the new instance becomes the default.
<<register_BatchtoolsParam>>=
default <- registered()
register(BatchtoolsParam(workers = 10), default = TRUE)
@

BatchtoolsParam has been moved to the top of the list and
is now the default.
<<register_BatchtoolsParam2>>=
names(registered())
bpparam()
@

Restore the original registry

<<register_restore>>=
for (param in rev(default))
    register(param)
@

\subsection{Functions}

\subsubsection{Parallel looping, vectorized and aggregate operations}

These are used in common functions, implemented as much as possible
for all back-ends. The functions (see the help pages, e.g.,
\Rcode{?bplapply} for a full definition) include

\begin{description}
  \item{\Rcode{bplapply(X, FUN, ...)}: }

    Apply in parallel a function \Rcode{FUN} to each element of \Rcode{X}.
    \Rcode{bplapply} invokes \Rcode{FUN} \Rcode{length(X)} times, each
    time with a single element of \Rcode{X}.

  \item{\Rcode{bpmapply(FUN, ...)}: }

    Apply in parallel a function \Rcode{FUN} to the first, second, etc.,
    elements of each argument in \ldots.

  \item{\Rcode{bpiterate(ITER, FUN, ...)}: }

    Apply in parallel a function \Rcode{FUN} to the output of
    function \Rcode{ITER}. Data chunks are returned by \Rcode{ITER}
    and distributed to parallel workers along with \Rcode{FUN}.
    Intended for iteration though an undefined number of data chunks
    (i.e., records in a file).

  \item{\Rcode{bpvec(X, FUN, ...)}: }

    Apply in parallel a function \Rcode{FUN} to subsets of \Rcode{X}.
    \Rcode{bpvec} invokes function \Rcode{FUN} as many times as there are
    cores or cluster nodes, with \Rcode{FUN} receiving a subset (typically
    more than 1 element, in contrast to \Rcode{bplapply}) of \Rcode{X}.

  \item{\Rcode{bpaggregate(x, data, FUN, ...)}: }

    Use the formula in \Rcode{x} to aggregate \Rcode{data} using \Rcode{FUN}.
\end{description}

\subsubsection{Parallel evaluation environment}

These functions query and control the state of the parallel evaluation
environment.

\begin{description}
  \item{\Rcode{bpisup(x)}: }
    Query a \Rcode{BiocParallelParam} back-end \Rcode{x} for its status.

  \item{\Rcode{bpworkers}; \Rcode{bpnworkers}: }
    Query a \Rcode{BiocParallelParam} back-end for the number of workers
    available for parallel evaluation.

  \item{\Rcode{bptasks}: }
    Divides a job (e.g., single call to *lapply function) into tasks.
    Applicable to \Rcode{MulticoreParam} only; \Rcode{DoparParam} and
    \Rcode{BatchtoolsParam} have their own approach to dividing a job
    among workers.

  \item{\Rcode{bpstart(x)}: }
    Start a parallel back end specified by \Rcode{BiocParallelParam} \Rcode{x},
    if possible.

  \item{\Rcode{bpstop(x)}: }
    Stop a parallel back end specified by \Rcode{BiocParallelParam} \Rcode{x}.
\end{description}

\subsubsection{Error handling and logging}

Logging and advanced error recovery is available in \Rcode{BiocParallel} 1.1.25
and later. For a more details see the vignette titled "Error Handling and
Logging":

<<error-vignette, eval=FALSE>>=
browseVignettes("BiocParallel")
@

\subsubsection{Locks and counters}

Inter-process (i.e., single machine) locks and counters are supported
using \Rcode{ipclock()}, \Rcode{ipcyield()}, and friends. Use these to
synchronize computation, e.g., allowing only a single process to write
to a file at a time.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Use cases}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Sample data are BAM files from a transcription profiling experiment
available in the \Rpackage{RNAseqData.HNRNPC.bam.chr14} package.

<<use_cases_data>>=
library(RNAseqData.HNRNPC.bam.chr14)
fls <- RNAseqData.HNRNPC.bam.chr14_BAMFILES
@

\subsection{Single machine}

Common approaches on a single machine are to use multiple cores in
forked processes, or to use clusters of independent processes.

For purely \R-based computations on non-Windows computers, there are
substantial benefits, such as shared memory, to be had using forked
processes. However, this approach is not portable across platforms,
and fails when code uses functionality, e.g., file or data base
access, that assumes only a single thread is accessing the
resource. While use of forked processes with \Rcode{MulticoreParam} is
an attractive solution for scripts using pure \R{} functionality,
robust and complex code often requires use of independent processes
and \Rcode{SnowParam}.

\subsubsection{Forked processes with \Rcode{MulticoreParam}}

This example counts overlaps between BAM files and a defined set of ranges.
First create a GRanges with regions of interest (in practice this could be
large).
<<forking_gr, message=FALSE>>=
library(GenomicAlignments) ## for GenomicRanges and readGAlignments()
gr <- GRanges("chr14", IRanges((1000:3999)*5000, width=1000))
@

A \Rclass{ScanBamParam} defines regions to extract from the files.
<<forking_param>>=
param <- ScanBamParam(which=range(gr))
@

FUN counts overlaps between the ranges in `gr` and the files.
<<forking_FUN>>=
FUN <- function(fl, param) {
  gal <- readGAlignments(fl, param = param)
  sum(countOverlaps(gr, gal))
}
@

All parameters necessary for running a job in a multi-core environment
are specified in the \Rclass{MulticoreParam} instance.
<<forking_default_multicore>>=
MulticoreParam()
@

The \BiocParallel{} functions, such as \Rfunction{bplapply}, use
information in the \Rclass{MulticoreParam} to set up the appropriate
back-end and pass relevant arguments to low-level functions.
\begin{verbatim}
> bplapply(fls[1:3], FUN, BPPARAM = MulticoreParam(), param = param)
$ERR127306
[1] 1185

$ERR127307
[1] 1123

$ERR127308
[1] 1241
\end{verbatim}

Shared memory environments eliminate the need to pass large data between
workers or load common packages. Note that in this code the GRanges data was
not passed to all workers in \Rcode{bplapply} and FUN did not need to load
\Biocpkg{GenomicAlignments} for access to the \Rcode{readGAlignments} function.

Problems with forked processes occur when code implementating
functionality used by the workers is not written in anticipation of
use by forked processes. One example is the database connection
underlying Bioconductor's \Rcode{org.*} packages. This psudoe-code
<<db_problems, eval = FALSE>>=
library(org.Hs.eg.db)
FUN <- function(x, ...) {
    ...
    mapIds(org.Hs.eg.db, ...)
    ...
}
bplapply(X, FUN, ..., BPPARAM = MulticoreParam())
@
%% 
is likely to fail, because \Rcode{library(org.Hs.eg.db)} opens a
database connection that is accessed by multiple processes. A solution
is to ensure that the database is opened independently in each process
<<db_solution_1, eval = FALSE>>=
FUN <- function(x, ...) {
    library(org.Hs.eg.db)
    ...
    mapIds(org.Hs.eg.db, ...)
    ...
}
bplapply(X, FUN, ..., BPPARAM = MulticoreParam())
@

\subsubsection{Clusters of independent processes with \Rcode{SnowParam}}
Both Windows and non-Windows machines can use the cluster approach to spawn
processes. \BiocParallel{} back-end choices for clusters on a single machine
are \Rclass{SnowParam} for configuring a Snow cluster or the \Rclass{DoparParam}
for use with the \Rpackage{foreach} package.

To re-run the counting example, FUN needs to modified such that `gr` is
passed as a formal argument and required libraries are loaded on each worker.
(In general, this is not necessary for functions defined in a package name
space, see Section~\ref{sec:developers}.)
<<cluster_FUN>>=
FUN <- function(fl, param, gr) {
  suppressPackageStartupMessages({
     library(GenomicAlignments)
  })
  gal <- readGAlignments(fl, param = param)
  sum(countOverlaps(gr, gal))
}
@

Define a 2-worker SOCK Snow cluster.
<<cluster_snow_param>>=
snow <- SnowParam(workers = 2, type = "SOCK")
@

A call to \Rcode{bplapply} with the \Rclass{SnowParam} creates the
cluster and distributes the work.
<<cluster_bplapply>>=
bplapply(fls[1:3], FUN, BPPARAM = snow, param = param, gr = gr)
@

The FUN written for the cluster adds some overhead due to the passing of the
GRanges and the loading of \Biocpkg{GenomicAlignments} on each
worker. This approach, however, has the advantage that it works on
most platforms and does not require a coding change when switching between
windows and non-windows machines.

If several \Rcode{bplapply()} statements are likely to require the
same resource, it often makes sense to create a cluster once using
\Rcode{bpstart()}. The workers are re-used by each call to
\Rcode{bplapply()}, so they do not have to re-load packages, etc.
<<db_solution_2, eval = FALSE>>=
register(SnowParam())  # default evaluation
bpstart()              # start the cluster
...
bplapply(X, FUN1, ...)
...
bplapply(X, FUN2, ...)  # re-use workers
...
bpstop()
@ 

\subsection{\emph{Ad hoc} cluster of multiple machines}

We use the term \emph{ad hoc} cluster to define a group of machines that can
communicate with each other and to which the user has password-less
log-in access. This
example uses a group of compute machines ("the rhinos") on the FHCRC network.

\subsubsection{Sockets}

On Linux and Mac OS X, a socket cluster is created across machines
by supplying machine names as the \Rcode{workers} argument to a
\Rclass{BiocParallelParam} instance instead of a number. Each name
represents an \R{} process; repeat names indicate multiple workers on
the same machine.

Create a \Rclass{SnowParam} with 2 cpus from `rhino01` and 1 from `rhino02`.
<<ad_hoc_sock_snow_param>>=
hosts <- c("rhino01", "rhino01", "rhino02")
param <- SnowParam(workers = hosts, type = "SOCK")
@

Execute FUN 4 times across the workers.
\begin{verbatim}
> FUN <- function(i) system("hostname", intern=TRUE)
> bplapply(1:4, FUN, BPPARAM = param)
[[1]]
[1] "rhino01"

[[2]]
[1] "rhino01"

[[3]]
[1] "rhino02"

[[4]]
[1] "rhino01"
\end{verbatim}

When creating a cluster across Windows machines \Rcode{workers} must
be IP addresses (e.g., "140.107.218.57") instead of machine names.

\subsubsection{MPI}

An MPI cluster across machines is created with \emph{mpirun} or
\emph{mpiexec} from the command line or a script. A list of machine
names provided as the -hostfile argument defines the mpi universe.

The hostfile requests 2 processors on 3 different machines.
\begin{verbatim}
rhino01 slots=2
rhino02 slots=2
rhino03 slots=2
\end{verbatim}

From the command line, start a single interactive \R{} process on the
current machine.
\begin{verbatim}
mpiexec --np 1 --hostfile hostfile R --vanilla
\end{verbatim}

Load \BiocParallel{} and create an MPI Snow cluster. The number of
\Rcode{workers} in \Rclass{SnowParam} should match the number of slots
requested in the hostfile. Using a smaller number of workers uses a
subset of the slots.
\begin{verbatim}
> library(BiocParallel)
> param <- SnowParam(workers = 6, type = "MPI")
\end{verbatim}

Execute FUN 6 times across the workers.
\begin{verbatim}
> FUN <- function(i) system("hostname", intern=TRUE)
> bplapply(1:6, FUN, BPPARAM = param)
bplapply(1:6, FUN, BPPARAM = param)
        6 slaves are spawned successfully. 0 failed.
[[1]]
[1] "rhino01"

[[2]]
[1] "rhino02"

[[3]]
[1] "rhino02"

[[4]]
[1] "rhino03"

[[5]]
[1] "rhino03"

[[6]]
[1] "rhino01"
\end{verbatim}

Batch jobs can be launched with mpiexec and R CMD BATCH. Code to be
executed is in `Rcode.R`.
\begin{verbatim}
mpiexec --hostfile hostfile R CMD BATCH Rcode.R
\end{verbatim}

\subsection{Clusters with schedulers}

Computer clusters are far from standardized, so the following may
require significant adaptation; it is written from experience here at
FHCRC, where we have a large cluster managed via SLURM. Nodes on the
cluster have shared disks and common system images, minimizing
complexity about making data resources available to individual nodes.
There are two simple models for use of the cluster, Cluster-centric
and R-centric.

\subsubsection{Cluster-centric}

The idea is to use cluster management software to allocate resources, and
then arrange for an \R{} script to be evaluated in the context of
allocated resources. NOTE: Depending on your cluster configuration it may
be necessary to add a line to the template file instructing workers to
use the version of R on the master / head node. Otherwise the default R on
the worker nodes will be used.

For SLURM, we might request space for 4 tasks
(with \verb+salloc+ or \verb+sbatch+), arrange to start the MPI
environment (with \verb+orterun+) and on a single node in that universe
run an \R{} script \verb+BiocParallel-MPI.R+. The command is
\begin{verbatim}
$ salloc -N 4 orterun -n 1 R -f BiocParallel-MPI.R
\end{verbatim}

The \R{} script might do the following, using MPI for parallel evaluation. Start
by loading necessary packages and defining \Rcode{FUN} work to be done
<<cluster-MPI-work, eval=FALSE>>=
library(BiocParallel)
library(Rmpi)
FUN <- function(i) system("hostname", intern=TRUE)
@
%%
Create a \Rclass{SnowParam} instance with the number of nodes equal to
the size of the MPI universe minus 1 (let one node dispatch jobs to
workers), and register this instance as the default
<<cluster-MPI, eval=FALSE>>=
param <- SnowParam(mpi.universe.size() - 1, "MPI")
register(param)
@
%%
Evaluate the work in parallel, process the results, clean up, and quit
<<cluster-MPI-do, eval=FALSE>>=
xx <- bplapply(1:100, FUN)
table(unlist(xx))
mpi.quit()
@
%%
The entire session is as follows:
\begin{knitrout}
\definecolor{shadecolor}{rgb}{0.969, 0.969, 0.969}\color{fgcolor}\begin{kframe}
\begin{verbatim}
$ salloc -N 4 orterun -n 1 R --vanilla -f BiocParallel-MPI.R
salloc: Job is in held state, pending scheduler release
salloc: Pending job allocation 6762292
salloc: job 6762292 queued and waiting for resources
salloc: job 6762292 has been allocated resources
salloc: Granted job allocation 6762292
## ...
> FUN <- function(i) system("hostname", intern=TRUE)
>
> library(BiocParallel)
> library(Rmpi)
> param <- SnowParam(mpi.universe.size() - 1, "MPI")
> register(param)
> xx <- bplapply(1:100, FUN)
        4 slaves are spawned successfully. 0 failed.
> table(unlist(xx))

gizmof13 gizmof71 gizmof86 gizmof88
      25       25       25       25
>
> mpi.quit()
salloc: Relinquishing job allocation 6762292
salloc: Job allocation 6762292 has been revoked.
\end{verbatim}
\end{kframe}
\end{knitrout}
One advantage of this approach is that the responsibility for managing
the cluster lies firmly with the cluster management software -- if one
wants more nodes, or needs special resources, then adjust parameters
to \verb+salloc+ (or \verb+sbatch+).

Notice that workers are spawned within the \Rcode{bplapply} function;
it might often make sense to more explicitly manage workers with
\Rfunction{bpstart} and \Rfunction{bpstop}, e.g.,
<<cluster-MPI-bpstart, eval=FALSE>>=
param <- bpstart(SnowParam(mpi.universe.size() - 1, "MPI"))
register(param)
xx <- bplapply(1:100, FUN)
bpstop(param)
mpi.quit()
@

\subsubsection{R-centric}

A more \R-centric approach might start an \R{}
script on the head node, and use \Rpackage{batchtools} to submit jobs
from within the \R{} session. One way of doing this is to create a
file containing a template for the job submission step, e.g., for
SLURM; a starting point might be found at
<<slurm>>=
tmpl <- system.file(package="batchtools", "templates", "slurm-simple.tmpl")
noquote(readLines(tmpl))
@ 
The \R{} script, run interactively or from the command line, might
then look like
<<cluster-batchtools, eval=FALSE>>=
## define work to be done
FUN <- function(i) system("hostname", intern=TRUE)

library(BiocParallel)

## register SLURM cluster instructions from the template file
param <- BatchtoolsParam(workers=5, cluster="slurm", template=tmpl)
register(param)

## do work
xx <- bplapply(1:100, FUN)
table(unlist(xx))
@
%%
The code runs on the head node until \Rcode{bplapply}, where the \R{}
script interacts with the SLURM scheduler to request a SLURM
allocation, run jobs, and retrieve results. The argument \Rcode{4} to
\Rcode{BatchtoolsParam} specifies the number of workers to request from
the scheduler; \Rcode{bplapply} divides the 100 jobs among the 4
workers. If \Rcode{BatchtoolsParam} had been created without specifying
any workers, then 100 jobs implied by the argument to \Rcode{bplapply}
would be associated with 100 tasks submitted to the scheduler.

Because cluster tasks are running in independent \R{} instances, and
often on physically separate machines, a convenient `best practice'
is to write \Rcode{FUN} in a `functional programming' manner, such
that all data required for the function is passed in as arguments or
(for large data) loaded implicitly or explicitly (e.g., via an \R{}
library) from disk.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Analyzing genomic data in \Bioconductor{}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

General strategies exist for handling large genomic data that are well
suited to \R{} programs. A manuscript titled \emph{Scalable Genomics with
\R{} and \Bioconductor{}}
(\url{http://arxiv.org/abs/1409.2864}) by Michael Lawrence
and Martin Morgan, reviews several of these approaches and demonstrate
implementation with \Bioconductor{} packages. Problem areas include
scalable processing, summarization and visualization. The techniques
presented include restricting queries, compressing data, iterating, and
parallel computing.

Ideas are presented in an approachable fashion within a framework of
common use cases. This is a benificial read for anyone anyone tackling
genomics problems in \R{}.


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{For developers}
\label{sec:developers}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Developers wishing to use \BiocParallel{} in their own packages should
include \BiocParallel{} in the \texttt{DESCRIPTION} file
\begin{verbatim}
    Imports: BiocParallel
\end{verbatim}
and import the functions they wish to use in the \texttt{NAMESPACE}
file, e.g.,
\begin{verbatim}
    importFrom(BiocParallel, bplapply)
\end{verbatim}
Then invoke the desired function in the code, e.g.,
<<devel-bplapply>>=
system.time(x <- bplapply(1:3, function(i) { Sys.sleep(i); i }))
unlist(x)
@
%%
This will use the back-end returned by \Rcode{bpparam()}, by default a
\Rcode{MulticoreParam()} instance or the user's preferred back-end if
they have used \Rcode{register()}.  The \Rcode{MulticoreParam} back-end
does not require any special configuration or set-up and is therefore
the safest option for developers. Unfortunately,
\Rcode{MulticoreParam} provides only serial evaluation on Windows.

Developers should document that their function uses \BiocParallel{}
functions on the man page, and should perhaps include in their
function signature an argument \Rcode{BPPARAM=bpparam()}. Developers
should NOT use `register()` in package code -- this sets a preference
that influences use of `bplapply()` and friends in all packages, not
just their package.

Developers wishing to invoke back-ends other than
\Rcode{MulticoreParam}, or to write code that works across Windows,
macOS and Linux, need to take special care to ensure that required
packages, data, and functions are available and loaded on the remote
nodes.

In \Rcode{bplapply()}, the environment of \Rcode{FUN} (other than the
global environment) is serialized to the workers. A consequence is
that, when \Rcode{FUN} is inside a package name space, other functions
available in the name space are available to \Rcode{FUN} on the
workers.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{\Rcode{sessionInfo()}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

<<sessionInfo, results="asis">>=
toLatex(sessionInfo())
@

\end{document}