Creating a PipelineDefinition

The PipelineDefinition object is, minimally, a set of functions consecutively executed on the output of the previous one, and optionally accompanied by evaluation and aggregation functions. A simple pipeline can be constructed as follows:

my_pip <- PipelineDefinition( list( 
  step1=function(x, param1){
    # do something with x and param1
    x
  },
  step2=function(x, method1, param2){
    get(method1)(x, param2) # apply method1 to x, with param2
  },
  step3=function(x, param3){
    x <- some_fancy_function(x, param3)
    # the functions can also output evaluation through the `intermediate_return` slot:
    e <- my_evaluation_function(x)
    list( x=x, intermediate_return=e)
  }
))
my_pip

## A PipelineDefinition object with the following steps:
##   - step1(x, param1)
##   - step2(x, method1, param2)
##   - step3(x, param3)

The PipelineDefinition can also include descriptions of each step or evaluation and aggregation functions. For example:

my_pip <- PipelineDefinition( 
  list( step1=function(x, meth1){ get(meth1)(x) },
        step2=function(x, meth2){ get(meth2)(x) } ),
  evaluation=list( step2=function(x) c(mean=mean(x), max=max(x)) ),
  description=list( step1="This steps applies meth1 to x.",
                    step2="This steps applies meth2 to x.")
)
my_pip

## A PipelineDefinition object with the following steps:
##   - step1(x, meth1)
## This steps applies meth1 to x.
##   - step2(x, meth2) * 
## This steps applies meth2 to x.

Running it with dummy data and functions:

datasets <- list( ds1=1:3, ds2=c(5,10,15))
alternatives <- list(meth1=c("log","sqrt"), meth2="cumsum")
tmpdir1 <- paste0(tempdir(),"/")
res <- runPipeline(datasets, alternatives, my_pip, output.prefix=tmpdir1)

## Running 2 pipeline settings on 2 datasets using 1 threads

## 
##                 Finished running on all datasets

## 
##                 Aggregating results...

## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE

## Aggregating evaluation results for: step2

## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE

res$evaluation$step2

##   dataset meth1  meth2      mean      max
## 1     ds1   log cumsum 0.8283022 1.791759
## 2     ds1  sqrt cumsum 2.5201593 4.146264
## 3     ds2   log cumsum 4.0471780 6.620073
## 4     ds2  sqrt cumsum 5.6352475 9.271329

Computing times can be accessed through res$elapsed; they can either be accessed as the pipeline total for each combination, or in a step-wise fashion. They can also be plotted using plotElapsed (see below). Evaluation results can be accessed through res$evaluation and plotted with evalHeatmap (see below).

Manipulating a PipelineDefinition

A number of generic methods are implemented for PipelineDefinition objects, including show, names, length, [, as.list. This means that, for instance, a step can be removed from a pipeline in the following way:

my_pip[-1]

## A PipelineDefinition object with the following steps:
##   - step2(x, meth2) * 
## This steps applies meth2 to x.

Steps can be added using the addPipelineStep function:

pip2 <- addPipelineStep(my_pip, name="newstep", after="step1")
pip2

## A PipelineDefinition object with the following steps:
##   - step1(x, meth1)
## This steps applies meth1 to x.
##   - newstep(x)
##   - step2(x, meth2) * 
## This steps applies meth2 to x.

Functions for the new step can be specified through the slots argument of addPipelineStep or afterwards through stepFn:

stepFn(pip2, "newstep", type="function") <- function(x, newparam){
  do_something(x, newparam)
}
pip2

## A PipelineDefinition object with the following steps:
##   - step1(x, meth1)
## This steps applies meth1 to x.
##   - newstep(x, newparam)
##   - step2(x, meth2) * 
## This steps applies meth2 to x.

Finally, the arguments() method can be used to extract the arguments for each step, and the defaultArguments methods can be used to get or set the default arguments.

Debugging

Errors typically come from the method wrappers, and debugging them from within runPipeline can be complicated by multithreading. For this reason runPipeline wraps additional information around error messages. This would be an example of the error report:

Error in dataset `mydataset` with parameters:
filt=lenient, norm=norm.test, k=20
in step `normalization`, evaluating command:
`pipDef[["normalization"]](x = x, norm = "norm.test")`
                Error:
Error in get(norm)(x): not enough cells!

The offending dataset, pipeline step and exact set of parameters is reported, enabling you to try to repeat the error outside, where it is easier to debug.

Skipping errors and fixing them afterwards

A useful behavior is for runPipeline to continue when encountering errors, so that successful results are not lost. This behavior can be enabled with the skipErrors=TRUE argument. We can illustrate this with our above dummy pipeline, by adding an alternative that is bound to give rise to an error:

fragile.fun <- function(x){
  if(any(x>10)) stop("Too big!")
  2^x
}
alternatives$meth1 <- c(alternatives$meth1, "fragile.fun")
res <- runPipeline(datasets, alternatives, my_pip, output.prefix=tmpdir1, skipErrors=TRUE)

## Running 3 pipeline settings on 2 datasets using 1 threads

## Error in dataset `ds2` with parameters:
## meth1=fragile.fun, meth2=cumsum
## in step `step1`, evaluating command:
## `pipDef[["step1"]](x = x, meth1 = "fragile.fun")`
##                 Error:
## Error in get(meth1)(x): Too big!

## 
##                 Finished running on all datasets

## 
## Some errors were encountered during the run:

##   dataset              step
## 1     ds2 meth1=fragile.fun

## All results were saved, but only those of analyses successful across all datasets will be retained for aggregation.

## 
##                 Aggregating results...

## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE

## Aggregating evaluation results for: step2

## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE
## Warning in type.convert.default(y[[i]]): 'as.is' should be specified by the
## caller; using TRUE

As indicated in the error message, only the analyses successful across all datasets were aggregated:

res$evaluation$step2

##   dataset meth1  meth2     mean      max
## 1     ds1  log2 cumsum 1.194988 2.584963
## 2     ds1  sqrt cumsum 2.520159 4.146264
## 3     ds2  log2 cumsum 5.838844 9.550747
## 4     ds2  sqrt cumsum 5.635248 9.271329

However, all successful results (e.g. including ‘fragile.fun’ on ‘ds1’) were saved in the non-aggregated files. This means that one can then correct the problematic function, re-run only the analyses that failed (i.e. only on ‘ds2’ in this case), and merge the results with the first batch (see ?mergePipelineResults).

Running times

Running times can be plotted with the plotElapsed function:

plotElapsed(res, agg.by=FALSE)

In this case we show all tested combinations, but in most cases the agg.by argument will be used to aggregate the results by parameters of interest.

Evaluation results can be plotted through the evalHeatmap function:

evalHeatmap(res, what=c("mean", "max"))

By default, the heatmap prints the actual metric values, while the colors is mapped to scaled data. Specifically, pipeComp map colors to the signed square-root of the number of (matrix-wise) Median Absolute Deviations from the (column-wise) median. We found this mapping to be particularly useful because differences in color have comparable meaning across datasets, and the color-scale is not primarily capturing outliers or baseline differences between datasets. This behavior can be changed with the scale parameter of the function. For more complex examples of using evalHeatmap, see its usage in the pipeComp_scRNA vignette, or refer to the function’s help.

For a more complex, ‘real-life’ example of the creation of a PipelineDefinition, see the pipeComp_dea vignette. For a complex example of evaluation outputs, see the pipeComp_scRNA vignette.