1 CBNplot: Bayesian network plot for clusterProfiler results

1.1 Introduction

The R package to infer and plot Bayesian networks. The network are inferred from expression data based on clusterProfiler or ReactomePA results. It makes use of libraries including clusterProfiler, ReactomePA, bnlearn, graphite and depmap. In this vignette, the description of functions and several use cases are depicted using GSE133624, which contains RNA-Seq data of bladder cancer. The more detail can be found on the book (https://noriakis.github.io/CBNplot/).

1.2 Installation

BiocManager::install("CBNplot")

1.3 Usage

1.3.1 The preprocessing and DEG identification of GSE133624

library(CBNplot)
library(bnlearn)
library(DESeq2)
library(org.Hs.eg.db)
library(GEOquery)
## Load dataset and make metadata
filePaths <- getGEOSuppFiles("GSE133624")
counts = read.table(rownames(filePaths)[1], header=1, row.names=1)
meta = sapply(colnames(counts), function (x) substring(x,1,1))
meta = data.frame(meta)
colnames(meta) = c("Condition")

dds <- DESeqDataSetFromMatrix(countData = counts,
                              colData = meta,
                              design= ~ Condition)
## Prefiltering
filt <- rowSums(counts(dds) < 10) > dim(meta)[1]*0.9
dds <- dds[!filt,]

## Perform DESeq2()
dds = DESeq(dds)
res = results(dds, pAdjustMethod = "bonferroni")

## apply variance stabilizing transformation
v = vst(dds, blind=FALSE)
vsted = assay(v)

## Define the input genes, and use clusterProfiler::bitr to convert the ID.
sig = subset(res, padj<0.05)
cand.entrez = clusterProfiler::bitr(rownames(sig),
  fromType="ENSEMBL", toType="ENTREZID", OrgDb=org.Hs.eg.db)$ENTREZID

## Perform enrichment analysis
pway = ReactomePA::enrichPathway(gene = cand.entrez)
pway = clusterProfiler::setReadable(pway, org.Hs.eg.db)

## Define including samples
incSample = rownames(subset(meta, Condition=="T"))

1.3.2 The use of CBNplot

1.4 bngeneplot

Then use CBNplot. Basically, you need to supply the enrichment analysis result, normalized expression value and samples to be included. For bngeneplot, the pathway number in the result slot of enrichment analysis results must be given.

bngeneplot(results = pway,exp = vsted,
  expSample = incSample, pathNum = 15)

Data frame of raw values used in the inference, data frame containing strength and direction, averaged network, and plot can be obtained by specifying returnNet=TRUE

ret <- bngeneplot(results = pway,exp = vsted,
  expSample = incSample, pathNum = 15, returnNet=TRUE)
ret$str |> head()
FALSE    from     to strength direction
FALSE 1 CENPQ  SPDL1     0.40 0.7500000
FALSE 2 CENPQ NUP133     0.30 0.5000000
FALSE 3 CENPQ  NUP37     0.45 0.6666667
FALSE 4 CENPQ  NDC80     0.25 0.4000000
FALSE 5 CENPQ   XPO1     0.40 0.5000000
FALSE 6 CENPQ  BIRC5     0.30 0.5000000

The resulting network can be converted to igraph object using bnlearn::as.igraph().

g <- bnlearn::as.igraph(ret$av)
igraph::evcent(g)$vector
##     CENPQ     SPDL1    NUP133     NUP37     NDC80      XPO1     BIRC5     PSMB5 
## 0.4043407 0.2454753 0.5489709 0.2138337 0.3534869 0.5340994 0.3925472 0.1279739 
##     PSMA6     PSMA7     CENPI     PSMA2     UBE2S    NUP107     CENPA    PSMD14 
## 0.4173347 0.3032090 0.2188213 0.3115492 0.3782319 0.2363872 0.7303752 0.3345827 
##     CDC20     CENPF     CENPL    KIF18A     PDS5A     ZWINT     CENPK     NUP85 
## 0.5345431 0.4975080 0.1678548 0.6729121 0.5689742 0.3450721 0.3997586 0.4273428 
##      SGO1     CDCA8     ESPL1      KNL1     CENPO     CENPE     KIF2C      NUF2 
## 0.4666727 0.8794683 0.6096135 0.2543809 0.6283373 0.4181057 0.6744517 0.3627395 
##     CDCA5    INCENP      DSN1     CENPU     SPC25     CENPH    ANAPC1    AHCTF1 
## 0.8284538 0.6608254 0.3850940 0.9349226 0.5191835 0.3064180 0.1884831 0.1738264 
##      BUB3      SKA1     BUB1B     PSMD4     SPC24      SGO2    MAD2L1     PTTG1 
## 0.3595813 0.8919014 0.5674659 0.2278502 0.5604701 0.2294218 1.0000000 0.5845318 
##     RAD21      PLK1    TUBA1C      BUB1    ZWILCH     UBE2C     CKAP5     AURKB 
## 0.2761026 0.1164351 0.3741177 0.7620561 0.3926130 0.6734771 0.4457428 0.5683693 
##      RCC2     KNTC1    ERCC6L 
## 0.4743328 0.5470337 0.5785199

1.5 bnpathplot

The relationship between pathways can be drawn by bnpathplot. The number to be included in the inference can be specified by nCategory.

bnpathplot(results = pway,exp = vsted,
  expSample = incSample, nCategory=10, shadowText = TRUE)

1.6 bngeneplotCustom and bnpathplotCustom

bngeneplotCustom and bnpathplotCustom can be used to customize visualization with more flexibility, like highlighting the nodes and edges of interest by glowEdgeNum and hub.

bnpathplotCustom(results = pway, exp = vsted, expSample = incSample,
  fontFamily="serif", glowEdgeNum=3, hub=3)

bngeneplotCustom(results = pway, exp = vsted, expSample = incSample,
  pathNum=15, fontFamily="sans", glowEdgeNum=NULL, hub=3)

The detailed usage for the package, like including covariates to the plot and probabilistic reasoning is available in the package documentation (https://noriakis.github.io/CBNplot/).

sessionInfo()
## R version 4.3.2 Patched (2023-11-13 r85521)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.04.3 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.18-bioc/R/lib/libRblas.so 
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_GB              LC_COLLATE=C              
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: America/New_York
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] GEOquery_2.70.0             org.Hs.eg.db_3.18.0        
##  [3] AnnotationDbi_1.64.1        DESeq2_1.42.0              
##  [5] SummarizedExperiment_1.32.0 Biobase_2.62.0             
##  [7] MatrixGenerics_1.14.0       matrixStats_1.2.0          
##  [9] GenomicRanges_1.54.1        GenomeInfoDb_1.38.3        
## [11] IRanges_2.36.0              S4Vectors_0.40.2           
## [13] BiocGenerics_0.48.1         bnlearn_4.9.1              
## [15] CBNplot_1.2.1               BiocStyle_2.30.0           
## 
## loaded via a namespace (and not attached):
##   [1] splines_4.3.2                 later_1.3.2                  
##   [3] bitops_1.0-7                  ggplotify_0.1.2              
##   [5] filelock_1.0.3                tibble_3.2.1                 
##   [7] oaqc_1.0                      polyclip_1.10-6              
##   [9] graph_1.80.0                  lifecycle_1.0.4              
##  [11] lattice_0.22-5                MASS_7.3-60                  
##  [13] ggdist_3.3.1                  magrittr_2.0.3               
##  [15] limma_3.58.1                  sass_0.4.8                   
##  [17] rmarkdown_2.25                jquerylib_0.1.4              
##  [19] yaml_2.3.8                    httpuv_1.6.13                
##  [21] cowplot_1.1.2                 DBI_1.2.0                    
##  [23] RColorBrewer_1.1-3            abind_1.4-5                  
##  [25] zlibbioc_1.48.0               purrr_1.0.2                  
##  [27] ggraph_2.1.0                  RCurl_1.98-1.13              
##  [29] yulab.utils_0.1.1             tweenr_2.0.2                 
##  [31] rappdirs_0.3.3                pvclust_2.2-0                
##  [33] GenomeInfoDbData_1.2.11       enrichplot_1.22.0            
##  [35] ggrepel_0.9.4                 tidytree_0.4.6               
##  [37] reactome.db_1.86.2            codetools_0.2-19             
##  [39] DelayedArray_0.28.0           xml2_1.3.6                   
##  [41] DOSE_3.28.2                   ggforce_0.4.1                
##  [43] tidyselect_1.2.0              aplot_0.2.2                  
##  [45] farver_2.1.1                  gmp_0.7-3                    
##  [47] viridis_0.6.4                 BiocFileCache_2.10.1         
##  [49] jsonlite_1.8.8                ellipsis_0.3.2               
##  [51] tidygraph_1.3.0               tools_4.3.2                  
##  [53] treeio_1.26.0                 Rcpp_1.0.11                  
##  [55] glue_1.6.2                    gridExtra_2.3                
##  [57] SparseArray_1.2.2             xfun_0.41                    
##  [59] qvalue_2.34.0                 distributional_0.3.2         
##  [61] dplyr_1.1.4                   withr_2.5.2                  
##  [63] BiocManager_1.30.22           fastmap_1.1.1                
##  [65] fansi_1.0.6                   digest_0.6.33                
##  [67] R6_2.5.1                      mime_0.12                    
##  [69] gridGraphics_0.5-1            colorspace_2.1-0             
##  [71] GO.db_3.18.0                  RSQLite_2.3.4                
##  [73] utf8_1.2.4                    tidyr_1.3.0                  
##  [75] generics_0.1.3                data.table_1.14.10           
##  [77] graphlayouts_1.0.2            httr_1.4.7                   
##  [79] S4Arrays_1.2.0                scatterpie_0.2.1             
##  [81] graphite_1.48.0               pkgconfig_2.0.3              
##  [83] gtable_0.3.4                  Rmpfr_0.9-4                  
##  [85] blob_1.2.4                    XVector_0.42.0               
##  [87] clusterProfiler_4.10.0        shadowtext_0.1.2             
##  [89] htmltools_0.5.7               bookdown_0.37                
##  [91] fgsea_1.28.0                  scales_1.3.0                 
##  [93] png_0.1-8                     ggfun_0.1.3                  
##  [95] knitr_1.45                    tzdb_0.4.0                   
##  [97] reshape2_1.4.4                nlme_3.1-164                 
##  [99] curl_5.2.0                    cachem_1.0.8                 
## [101] stringr_1.5.1                 BiocVersion_3.18.1           
## [103] parallel_4.3.2                HDO.db_0.99.1                
## [105] ReactomePA_1.46.0             pillar_1.9.0                 
## [107] grid_4.3.2                    vctrs_0.6.5                  
## [109] promises_1.2.1                dbplyr_2.4.0                 
## [111] xtable_1.8-4                  evaluate_0.23                
## [113] magick_2.8.2                  readr_2.1.4                  
## [115] cli_3.6.2                     locfit_1.5-9.8               
## [117] compiler_4.3.2                rlang_1.1.2                  
## [119] crayon_1.5.2                  labeling_0.4.3               
## [121] plyr_1.8.9                    fs_1.6.3                     
## [123] stringi_1.8.3                 viridisLite_0.4.2            
## [125] BiocParallel_1.36.0           munsell_0.5.0                
## [127] Biostrings_2.70.1             lazyeval_0.2.2               
## [129] GOSemSim_2.28.0               Matrix_1.6-4                 
## [131] ExperimentHub_2.10.0          hms_1.1.3                    
## [133] patchwork_1.1.3               bit64_4.0.5                  
## [135] ggplot2_3.4.4                 statmod_1.5.0                
## [137] KEGGREST_1.42.0               shiny_1.8.0                  
## [139] highr_0.10                    interactiveDisplayBase_1.40.0
## [141] AnnotationHub_3.10.0          igraph_1.6.0                 
## [143] memoise_2.0.1                 bslib_0.6.1                  
## [145] ggtree_3.10.0                 fastmatch_1.1-4              
## [147] bit_4.0.5                     ape_5.7-1                    
## [149] gson_0.1.0