Code

    library(BulkSignalR)

## Cluster-based differential analysis

By default, `BulkSignalR` inference of L-R interactions relies on correlation analysis across all the samples. That is, a putative L-R interaction is assessed for coherent expression of the ligand, the receptor, and target genes in a pathway below receptor across all the samples.

Alternatively, it is possible to define sample clusters based on an independent analysis (outside of `BulkSignalR`) to group similar samples. Based on cluster definitions, it is for instance possible to assess whether L-R interactions are significantly stronger in a cluster compared to all the samples, to all the samples but the chosen cluster, or in comparison with another cluster. In practice, the user must simply provide a table of differential gene or protein expression analysis that is relevant for the selected comparison.

In this mode, L-R interactions are inferred based on gene or protein regulation-associated P-values instead of P-values resulting from the correlation analysis in the default mode. As already stated, differential P-value estimations are left to the user to compute. With bulk data, commonly used libraries such as `DESeq2`, `EdgeR`, or `limma` constitute obvious options.

In the next chunk of code, we illustrate the differential mode using Salivary Duct Carcinoma (SDC) samples where we compare two clusters of patients.

We first create a BSRDataModelComp object as follows:

Code

    data(sdc, package = "BulkSignalR")
    normal <- grep("^N", names(sdc))
    bsrdm.comp <- BSRDataModelComp(sdc[, -normal])

To make a very small example we generate random values, but user should provide his own logFC and associated pvalues from DGE ouputs.

Code

    colA <- as.integer(1:3)
    colB <- as.integer(12:15)

    n <- nrow(ncounts(bsrdm.comp))
    stats <- data.frame(pval = runif(n),
    logFC = rnorm(n, 0, 2),
    expr = runif(n, 0, 10))
    rownames(stats) <- rownames(ncounts(bsrdm.comp))

We define the cluster comparison and add it.

Code

    bsrcc <- BSRClusterComp(bsrdm.comp, colA, colB, stats)
    bsrdm.comp <- addClusterComp(bsrdm.comp, bsrcc, "random.example")

Finally we infer ligand-receptor interactions from the comparison. We use a subset of the reference to speed up inference in the context of the vignette.

Code

    subset <- c("REACTOME_BASIGIN_INTERACTIONS",
    "REACTOME_SYNDECAN_INTERACTIONS",
    "REACTOME_ECM_PROTEOGLYCANS",
    "REACTOME_CELL_JUNCTION_ORGANIZATION")

    reactSubset <- BulkSignalR:::.SignalR$BulkSignalR_Reactome[
    BulkSignalR:::.SignalR$BulkSignalR_Reactome$`Reactome name` %in% subset,]

    resetPathways(dataframe = reactSubset,
    resourceName = "Reactome")

    bsrinf.comp <- BSRInferenceComp(bsrdm.comp,
    reference="REACTOME",
    max.pval = 1,
    "random.example")

    head(LRinter(bsrinf.comp))

    ##        L     R         pw.id                             pw.name         pval
    ## 3  ADAM2 ITGA9 R-HSA-3000178          REACTOME_ECM_PROTEOGLYCANS 0.6101759164
    ## 46  CSF2 ITGA9 R-HSA-3000178          REACTOME_ECM_PROTEOGLYCANS 0.1530488692
    ## 50 CYR61  CAV1  R-HSA-210991       REACTOME_BASIGIN_INTERACTIONS 0.0023845886
    ## 69   FGA  CDH5  R-HSA-446728 REACTOME_CELL_JUNCTION_ORGANIZATION 0.0002551668
    ## 70   FGA ITGAX R-HSA-3000178          REACTOME_ECM_PROTEOGLYCANS 0.1424734945
    ## 87   FN1 ITGA9 R-HSA-3000178          REACTOME_ECM_PROTEOGLYCANS 0.5626697942
    ##           qval  L.logFC  R.logFC    LR.pval LR.corr rank len rank.pval
    ## 3  0.639231912 2.364800 1.501243 0.63964106       1    8  16 0.7255144
    ## 46 0.306097738 1.159814 1.501243 0.16043954       1    8  16 0.7255144
    ## 50 0.026230475 1.928128 1.739945 0.01539104       1    4   8 0.3875056
    ## 69 0.005613669 1.051807 2.652095 0.00315793       1   18  33 0.3522131
    ## 70 0.306097738 1.051807 2.517955 0.14582282       1    9  17 0.7255144
    ## 87 0.639231912 1.103627 1.501243 0.58984089       1    8  16 0.7255144
    ##    rank.corr  LR.score   L.expr   R.expr
    ## 3          1 0.2012374 5.600278 8.173384
    ## 46         1 0.2010613 5.573067 8.173384
    ## 50         1 0.2147591 7.505326 8.998575
    ## 69         1 0.2238587 8.932501 9.888175
    ## 70         1 0.1439432 8.932501 1.251274
    ## 87         1 0.1999241 5.400985 8.173384

## Technical notes

The three basic `BulkSignalR` S4 classes BSRDataModel, BSRInference, and BSRSignature were derived into new classes to add the functionalities required by the differential mode:

- BSRDataModelComp, denoted here as `bsrdm.comp` is a daughter class of BSRDataModel, previously denoted `bsrdm`
- BSRInferenceComp, denoted here as `bsrinf.comp` is a daughter class of BSRInference, previously denoted `bsrinf`
- BSRSignatureComp, denoted here as `bsrsig.comp` is a daughter class of BSRSignature, previously denoted `bsrsig`

The new S4 class BSRClusterComp is meant to represent the comparison between two clusters of samples. It primarily contains the results of the user-provided differential analysis in a table containing gene/protein names, P-values, and log-fold-changes (logFC).

Usually, we recommend instantiating a `BSRDataModel` object to contain the expression matrix underlying the differential analysis. Then, it is promoted into a `BSRDataModelComp` object such that cluster definitions and comparisons can be added to it. This can be for instance done with `as(bsrdm, "BSRDataModelComp")`. The `addClusterComp` method of a `BSRDataModelComp` object allows adding one or several comparisons between clusters, each added comparison being defined in a `BSRClusterComp` object.

Thank you for reading this guide and for using `BulkSignalR`.

## Session Information

Code

    sessionInfo()

    ## R version 4.5.1 (2025-06-13 ucrt)
    ## Platform: x86_64-w64-mingw32/x64
    ## Running under: Windows Server 2022 x64 (build 20348)
    ##
    ## Matrix products: default
    ##   LAPACK version 3.12.1
    ##
    ## locale:
    ## [1] LC_COLLATE=C
    ## [2] LC_CTYPE=English_United States.utf8
    ## [3] LC_MONETARY=English_United States.utf8
    ## [4] LC_NUMERIC=C
    ## [5] LC_TIME=English_United States.utf8
    ##
    ## time zone: America/New_York
    ## tzcode source: internal
    ##
    ## attached base packages:
    ## [1] stats     graphics  grDevices utils     datasets  methods   base
    ##
    ## other attached packages:
    ## [1] BulkSignalR_1.1.4
    ##
    ## loaded via a namespace (and not attached):
    ##   [1] DBI_1.2.3                   bitops_1.0-9
    ##   [3] grr_0.9.5                   gridExtra_2.3
    ##   [5] httr2_1.2.1                 rlang_1.1.6
    ##   [7] magrittr_2.0.3              clue_0.3-66
    ##   [9] GetoptLong_1.0.5            matrixStats_1.5.0
    ##  [11] compiler_4.5.1              RSQLite_2.4.2
    ##  [13] gprofiler2_0.2.3            png_0.1-8
    ##  [15] vctrs_0.6.5                 ggalluvial_0.12.5
    ##  [17] pkgconfig_2.0.3             SpatialExperiment_1.19.1
    ##  [19] shape_1.4.6.1               crayon_1.5.3
    ##  [21] fastmap_1.2.0               backports_1.5.0
    ##  [23] dbplyr_2.5.0                magick_2.8.7
    ##  [25] XVector_0.49.0              rmarkdown_2.29
    ##  [27] purrr_1.1.0                 bit_4.6.0
    ##  [29] xfun_0.52                   glmnet_4.1-10
    ##  [31] cachem_1.1.0                aplot_0.2.8
    ##  [33] jsonlite_2.0.0              blob_1.2.4
    ##  [35] DelayedArray_0.35.2         broom_1.0.9
    ##  [37] parallel_4.5.1              cluster_2.1.8.1
    ##  [39] R6_2.6.1                    bslib_0.9.0
    ##  [41] RColorBrewer_1.1-3          car_3.1-3
    ##  [43] GenomicRanges_1.61.1        jquerylib_0.1.4
    ##  [45] Rcpp_1.1.0                  Seqinfo_0.99.2
    ##  [47] SummarizedExperiment_1.39.1 iterators_1.0.14
    ##  [49] knitr_1.50                  IRanges_2.43.0
    ##  [51] Matrix_1.7-3                splines_4.5.1
    ##  [53] igraph_2.1.4                tidyselect_1.2.1
    ##  [55] dichromat_2.0-0.1           abind_1.4-8
    ##  [57] yaml_2.3.10                 doParallel_1.0.17
    ##  [59] codetools_0.2-20            curl_6.4.0
    ##  [61] lattice_0.22-7              tibble_3.3.0
    ##  [63] withr_3.0.2                 treeio_1.33.0
    ##  [65] Biobase_2.69.0              evaluate_1.0.4
    ##  [67] Rtsne_0.17                  gridGraphics_0.5-1
    ##  [69] survival_3.8-3              BiocFileCache_2.99.5
    ##  [71] circlize_0.4.16             ggtree_3.17.1
    ##  [73] pillar_1.11.0               ggpubr_0.6.1
    ##  [75] filelock_1.0.3              carData_3.0-5
    ##  [77] MatrixGenerics_1.21.0       foreach_1.5.2
    ##  [79] stats4_4.5.1                plotly_4.11.0
    ##  [81] ggfun_0.2.0                 generics_0.1.4
    ##  [83] RCurl_1.98-1.17             S4Vectors_0.47.0
    ##  [85] ggplot2_3.5.2               tidytree_0.4.6
    ##  [87] scales_1.4.0                glue_1.8.0
    ##  [89] orthogene_1.15.0            lazyeval_0.2.2
    ##  [91] tools_4.5.1                 data.table_1.17.8
    ##  [93] ggsignif_0.6.4              babelgene_22.9
    ##  [95] RANN_2.6.2                  fs_1.6.6
    ##  [97] grid_4.5.1                  ape_5.8-1
    ##  [99] tidyr_1.3.1                 colorspace_2.1-1
    ## [101] SingleCellExperiment_1.31.1 patchwork_1.3.1
    ## [103] nlme_3.1-168                homologene_1.4.68.19.3.27
    ## [105] Formula_1.2-5               cli_3.6.5
    ## [107] rappdirs_0.3.3              viridisLite_0.4.2
    ## [109] S4Arrays_1.9.1              ComplexHeatmap_2.25.2
    ## [111] dplyr_1.1.4                 gtable_0.3.6
    ## [113] stabledist_0.7-2            yulab.utils_0.2.0
    ## [115] rstatix_0.7.2               sass_0.4.10
    ## [117] digest_0.6.37               BiocGenerics_0.55.1
    ## [119] ggplotify_0.1.2             SparseArray_1.9.1
    ## [121] ggrepel_0.9.6               htmlwidgets_1.6.4
    ## [123] rjson_0.2.23                farver_2.1.2
    ## [125] memoise_2.0.1               htmltools_0.5.8.1
    ## [127] multtest_2.65.0             lifecycle_1.0.4
    ## [129] httr_1.4.7                  GlobalOptions_0.1.2
    ## [131] bit64_4.6.0-1               MASS_7.3-65