Status

The integration between NicheNet and OmniPath is a fresh development, as of May 2021 we have already tested the basic functionality. To use the features described here you will need the latest git version of OmnipathR, which we update often as we test and improve the elements of the pipeline. Bug reports, questions and suggestions are welcome, please open an issue on github.

Introduction

The OmnipathR package contains methods to build NichNet[1] format networks from OmniPath[2] data and from the resources used in the original NicheNet study. With a few exceptions, almost all resources used originally by NicheNet are available by OmnipathR both in their raw format or after minimal preprocessing, and in the format suitable for NicheNet. However, the networks built by OmnipathR are not completely identical with the ones used by NicheNet due to the complexity of processing such a large amount of data. The advantage of building these networks with OmnipathR is the transparent and reproducible process, as the origin of any data record can be easily traced back to its original source. The data used in the NicheNet publication is deposited at Zenodo. The other part of the prior knowledge used by NicheNet is a collection of expression data from ligand perturbation experiments. In OmnipathR this is taken directly from the Zenodo repository.

Apart from building the prior knowledge, OmnipathR provides a number of glue methods which connect elements of the NicheNet workflow, ultimately making it possible to run everything with a single function call, starting from the processed transcriptomics data and going until the predicted ligand activities. To implement such a workflow we followed the analysis from the case study in the OmniPath paper and the vignettes from NicheNet. The workflow consists of thin wrappers around some of the NicheNet methods. Parameters for these methods can be provided from the top level call.

The transcriptomics data is the main input for the pipeline which is specific for the study. The transcriptomics data can be processed by standard tools, such as DESeq[3] for bulk and Seurat[4] for single-cell transcriptomics. For the latter more guidance is available in the NicheNet vignettes. The pipeline presented here only requires the list of genes expressed in the transmitter and receiver cells (in case of autocrine signaling these are the same); and a list with genes of interest. Genes of interest might come from the investigation subject or gene set enrichment analysis (see an example here).

Run the workflow by a single call

The nichenet_main function executes all the workflow elements, hence it can be used to run the complete workflow. Here we show it as an example, and in the next sections we take a closer look at the individual steps. The first four arguments comes from the study data and objective, as discussed above. All parameters of the workflow can be overridden here, for example, for the mlrbo_optimization we override the ncores argument, and exclude CPDB from the signaling network and set the confidence threshold for EVEX.

library(OmnipathR)
library(nichenetr)

nichenet_results <- nichenet_main(
    expressed_genes_transmitter = expressed_genes_transmitter,
    expressed_genes_receiver = expressed_genes_receiver,
    genes_of_interest = genes_of_interest,
    background_genes = background_genes,
    signaling_network = list(
        cpdb = NULL,
        evex = list(min_confidence = 1.0)
    ),
    gr_network = list(only_omnipath = TRUE),
    n_top_ligands = 20,
    mlrmbo_optimization_param = list(ncores = 4)
)

Only model building

Without transcriptomics data, it is possible to build a NicheNet model which later can be used for the ligand activity prediction. Hence if the first four arguments above are NULL the pipeline will run only the model building. Still, in this tutorial we wrap all blocks into eval = FALSE because the model optimization is computationally intensive, it might take hours to run. To create the NicheNet model it’s enough to call nichenet_main without arguments, all further arguments are optional override of package defaults.

nichenet_model <- nichenet_main()

Installing packages and loading nichenetr

Before running the pipeline it’s necessary to install nichenetr and its dependencies. The nichenetr package is available only from github, while the dependencies are distributed in standard repositories like CRAN.

require(devtools)
install_github('saeyslab/nichenetr')

This procedure should install the dependencies as well. Some more specific dependencies are mlrMBO, parallelMap, ParamHelpers and smoof.

The nichenetr package is suggested by OmnipathR but not loaded and attached by default. This might cause issues, for example NicheNet requires access to some of its lazy loaded external data. To address this, you can either attach nichenetr:

library(nichenetr)

Or call a function which loads the datasets into the global environment, so NicheNet can access them:

nichenet_workarounds()

Testing the pipeline

With a fresh setup, as a first step it is recommended to test the pipeline, to make sure it is able to run in your particular environment, all packages are installed and loaded properly. The optimization process takes hours to run, while the test included in OmnipathR takes only around 10 minutes. This way you can avoid some errors to come up after running the analysis already for many hours. Although it’s never guaranteed that you won’t have errors in the pipeline as the outcome depends on the input data and certain conditions are not handled in nichenetr, mlr and mlrMBO. If you get an error during the optimization process, start it over a few times to see if it’s consistent. Sometimes it helps also to restart your R session. I would recommend to try it around 10 times before suspecting a bug. See more details in the docs of nichenet_test.

The nichenet_test function builds a tiny network and runs the pipeline with a highly reduced number of parameters and iteration steps in the optimization. Also it inputs dummy gene sets as genes of interest and background genes.

require(OmnipathR)
nichenet_workarounds()
nichenet_test()

Workflow steps

networks <- nichenet_networks()

networks <- nichenet_networks(only_omnipath = TRUE)

networks <- nichenet_networks(
    only_omnipath = TRUE,
    signaling_network = list(
        omnipath = list(
            resources = c('SIGNOR', 'PhosphoSite')
        )
    ),
    lr_network = list(
        omnipath = list(
            transmitter_param = list(parent = c('ligand', 'secreted_enzyme')),
            receiver_param = list(parent = c('receptor', 'transporter'))
        )
    ),
    gr_network = list(
        omnipath = list(
            resources = 'DoRothEA',
            dorothea_levels = 'A'
        )
    )
)

signaling_network <- nichenet_signaling_network(
    cpdb = list(
        complex_max_size = 1,
        min_score = .98
    ),
    evex = list(
        min_confidence = 2
    ),
    pathwaycommons = NULL
)

evex_signaling <- nichenet_signaling_network_evex(top_confidence = .9)

evex <- evex_download(remove_negatives = FALSE)

nichenet_results <- nichenet_main(
    quality_filter_param = list(
        min_curation_effort = 1,
        consensus_percentile = 30
    )
)

nichenet_results <- nichenet_main(small = TRUE)

expression <- nichenet_expression_data()

lr_network <- nichenet_lr_network()
expression <- nichenet_remove_orphan_ligands(
    expression = expression,
    lr_network = lr_network
)

optimization_results <- nichenet_optimization(
    networks = networks,
    expression = expression,
    mlrmbo_optimization_param = list(ncores = 4)
)

options(omnipath.nichenet_results_dir = 'my/nichenet/dir')
nichenet_results_dir()
# [1] "my/nichenet/dir"

nichenet_model <- nichenet_build_model(
    optimization_results = optimization_results,
    networks = networks,
)

lt_matrix <- nichenet_ligand_target_matrix(
    nichenet_model$weighted_networks,
    networks$lr_network,
    nichenet_model$optimized_parameters
)

genes_of_interest <- sample(rownames(ligand_target_matrix), 50)
background_genes <- setdiff(
    rownames(ligand_target_matrix),
    genes_of_interest
)
expressed_genes_transmitter <- union(
    unlist(purrr::map(networks, 'from')),
    unlist(purrr::map(networks, 'to'))
)
expressed_genes_receiver <- expressed_genes_transmitter

ligand_activities <- nichenet_ligand_activities(
    ligand_target_matrix = lt_matrix,
    lr_network = networks$lr_network,
    expressed_genes_transmitter = expressed_genes_transmitter,
    expressed_genes_receiver = expressed_genes_receiver,
    genes_of_interest = genes_of_interest
)

lt_links <- nichenet_ligand_target_links(
    ligand_activities = ligand_activities,
    ligand_target_matrix = lt_matrix,
    genes_of_interest = genes_of_interest,
    n_top_ligands = 20,
    n_top_targets = 100
)

Further steps

The NicheNet workflow in OmnipathR is implemented until this point. The results can be visualized by the methods included in nichenetr, can be analysed for enrichment of functional properties, or the signaling network connecting the ligands to targets can be reconstructed and analysed, as it’s shown in the OmniPath case study.