Contents

1 Introduction

1.1 About Cell-Cell Interaction (CCI)

Due to the rapid development of single-cell RNA-Seq (scRNA-Seq) technologies, wide variety of cell types such as multiple organs of healty person, stem cell niche and cancer stem cell have been found. Such complex systems are composed by communication between cells (cell-cell interaction or CCI).

Many CCI studies are based on the ligand-receptor (L-R)-pair list of FANTOM5 project1 Jordan A. Ramilowski, A draft network of ligand-receptor-mediated multicellular signalling in human, Nature Communications, 2015 as the evidence of CCI (http://fantom.gsc.riken.jp/5/suppl/Ramilowski_et_al_2015/data/PairsLigRec.txt). The project proposed the L-R-candidate genes by following two reasons.

  1. Subcellular Localization
    1. Known Annotation (UniProtKB and HPRD) : The term “Secreted” for candidate ligand genes and “Plasma Membrane” for candidate receptor genes
    2. Computational Prediction (LocTree3 and PolyPhobius)
  2. Physical Binding of Proteins : Experimentally validated PPI (protein-protein interaction) information of HPRD and STRING

The project also merged the data with previous L-R database such as IUPHAR/DLRP/HPMR and filter out the list without PMIDs.

Here, we extend the similar approarch to the construction of the L-R-pair lists of 12 organsisms and implemented as multiple R/Bioconductor annotation packages for sustainable maintainance (LRBaseDbi and LRBase.XXX.eg.db-type packages (Figure 1). XXX is the abbreviation of scientific name of organisms such as LRBase.Hsa.eg.db for L-R database of Homo sapiens). Besides, we also developed (scTensor), whichis a method to detect CCI and the CCI-releated L-R pairs simutaneously. This document provides the way to use LRBaseDbi, LRBase.XXX.eg.db-type packages, and scTensor package.

Figure 1 : Workflow of L-R-related packages

Figure 1 : Workflow of L-R-related packages

2 Usage

2.1 LRBase.XXX.eg.db (ligand-receptor database of 12 organisms)

To create the L-R-list of 12 organisms, we used the information about the subcelluar localization of proteins from SWISSPROT (knowledge databaase) and TREMBL (computational prediction). We also used the PPI information from STRING database (Figure 1). The proteins which is assigned to the term “Secreted” and “Cellular Membrane” are retrived as candidate ligand and receptor, respectively. Finally, the L-R-pairs which is registed in STRING are extracted as candidate L-R-lists.

Following 12 organisms are implmented as LRBase.XXX.eg.db-type packages.

Table 1: The summary of 12 organisms of
LRBase\(_\cdot\)XXX\(_\cdot\)eg\(_\cdot\)db packages.
Organisms Package Name # Swiss-Prot (Secreted / Membrane) # TrEMBL (Secreted / Membrane) # STRING (PPI) # Pair (Swiss-Prot × STRING) # Pair (TrEMBL × STRING)
\(\textit{Homo sapiens}\) LRBase.Hsa.eg.db 1592 / 2269 176 / 334 18838 21882 472
\(\textit{Mus musculus}\) LRBase.Mmu.eg.db 1309 / 1806 325 / 1555 19715 16386 476
\(\textit{Arabidopsis thaliana}\) LRBase.Ath.eg.db 1260 / 1001 244 / 80 24174 8697 94
\(\textit{Rattus norvegicus}\) LRBase.Rno.eg.db 643 / 983 232 / 1229 19963 5270 65
\(\textit{Bos taurus}\) LRBase.Bta.eg.db 517 / 448 192 / 390 18349 2220 237
\(\textit{Caenorhabditis elegans}\) LRBase.Cel.eg.db 198 / 247 28 / 60 13545 106 1
\(\textit{Drosophila melanogaster}\) LRBase.Dme.eg.db 249 / 333 89 / 148 11903 384 9
\(\textit{Danio rerio}\) LRBase.Dre.eg.db 119 / 169 318 / 376 21746 99 432
\(\textit{Gallus gallus}\) LRBase.Gga.eg.db 175 / 173 185 / 154 13084 140 105
\(\textit{Pongo abelii}\) LRBase.Pab.eg.db 80 / 134 212 / 211 16691 34 184
\(\textit{Xenopus Silurana tropicalis}\) LRBase.Xtr.eg.db 57 / 83 141 / 114 15338 19 107
\(\textit{Sus scrofa}\) LRBase.Ssc.eg.db 223 / 153 202 / 445 18683 277 130

2.1.1 columns, keytypes, keys, and select

Some data access functions are available for LRBase.XXX.eg.db-type packages. The all data table are retrieved by 4 functions defined by AnnotationDbi; columns, keytypes, keys, and select and commonly implemented by LRBaseDbi package. columns returns the rows which we can retrieve in LRBase.XXX.eg.db-type packages. keytypes returns the rows which can be used as the optional parameter in keys and select functions against LRBase.XXX.eg.db-type packages. keys function returns the value of keytype. select function returns the rows in particular columns, which are having user-specified keys. This function returns the result as a dataframe. See the vignette of AnnotationDbi for more details.

## Loading required package: LRBaseDbi
columns(LRBase.Hsa.eg.db)
## [1] "GENEID_L" "GENEID_R" "SOURCEDB" "SOURCEID"
keytypes(LRBase.Hsa.eg.db)
## [1] "GENEID_L" "GENEID_R" "SOURCEDB" "SOURCEID"
key_HSA <- keys(LRBase.Hsa.eg.db, keytype="GENEID_L")
head(select(LRBase.Hsa.eg.db, keys=key_HSA[1:2],
            columns=c("GENEID_L", "GENEID_R"), keytype="GENEID_L"))
##   GENEID_L GENEID_R
## 1     5652       40
## 2      656       41
## 3      656       91
## 4      656       93
## 5     5652      334
## 6      656      658

2.1.2 Other functions

Other additional functions like species, nomenclature, and listDatabases are available. In each LRBase.XXX.eg.db-type package, species function returns the common name and nomenclature returns the scientific name. listDatabases function returns the source of data. dbInfo returns the information of the package. dbfile returns the directory where sqlite file is stored. dbschema returns the schema of database. dbconn returns the connection to the sqlite database.

lrPackageName(LRBase.Hsa.eg.db)
## [1] "LRBase.Hsa.eg.db"
lrNomenclature(LRBase.Hsa.eg.db)
## [1] "Homo sapiens"
species(LRBase.Hsa.eg.db)
## [1] "Human"
lrListDatabases(LRBase.Hsa.eg.db)
##           SOURCEDB
## 1 SWISSPROT_STRING
## 2    TREMBL_STRING
## 3         SOURCEDB
## 4           IUPHAR
## 5             DLRP
lrVersion(LRBase.Hsa.eg.db)
##        NAME VALUE
## 1 LRVERSION  2018
dbInfo(LRBase.Hsa.eg.db)
##               NAME                                              VALUE
## 1       SOURCEDATE                                        1-June-2018
## 2      SOURCENAME1                                          SWISSPROT
## 3      SOURCENAME2                                             TREMBL
## 4      SOURCENAME3                                             STRING
## 5       SOURCEURL1 http://www.uniprot.org/uniprot/?query=reviewed:yes
## 6       SOURCEURL2  http://www.uniprot.org/uniprot/?query=reviewed:no
## 7       SOURCEURL3              https://string-db.org/cgi/download.pl
## 8         DBSCHEMA                                   LRBase.Hsa.eg.db
## 9  DBSCHEMAVERSION                                                1.0
## 10        ORGANISM                                       Homo sapiens
## 11         SPECIES                                              Human
## 12         package                                      AnnotationDbi
## 13         Db type                                           LRBaseDb
## 14       LRVERSION                                               2018
dbfile(LRBase.Hsa.eg.db)
## [1] "/home/biocbuild/bbs-3.9-bioc/R/library/LRBase.Hsa.eg.db/extdata/LRBase.Hsa.eg.db.sqlite"
dbschema(LRBase.Hsa.eg.db)
## [1] "CREATE TABLE `METADATA` (\n  `NAME` TEXT,\n  `VALUE` TEXT\n)"                                           
## [2] "CREATE TABLE `DATA` (\n  `GENEID_L` TEXT,\n  `GENEID_R` TEXT,\n  `SOURCEID` TEXT,\n  `SOURCEDB` TEXT\n)"
dbconn(LRBase.Hsa.eg.db)
## <SQLiteConnection>
##   Path: /home/biocbuild/bbs-3.9-bioc/R/library/LRBase.Hsa.eg.db/extdata/LRBase.Hsa.eg.db.sqlite
##   Extensions: TRUE

Combined with dbGetQuery function of RSQLite package, more complicated queries also can be submitted.

library("RSQLite")
dbGetQuery(dbconn(LRBase.Hsa.eg.db),
  "SELECT * FROM DATA WHERE GENEID_L = '9068' AND GENEID_R = '14' LIMIT 10")
## [1] GENEID_L GENEID_R SOURCEID SOURCEDB
## <0 rows> (or 0-length row.names)

2.2 LRBaseDbi (Class definition and meta-packaging)

LRBaseDbi regulates the class definition of LRBaseDb object instantiated from LRBaseDb-class. Besides, LRBaseDbi package generates user’s original LRBase.XXX.eg.db-type packages by makeLRBasePackage function. This function is inspired by our previous package MeSHDbi, which constructs user’s original MeSH.XXX.eg.db-type packages. Here we call this function “meta”-packaging. The 12 LRBase.XXX.eg.db-type packages described above are also generated by this “meta”-packaging. In this case, only user have to specify are 1. a L-R-list containing the columns “GENEID_L” (ligand NCBI Gene IDs) and “GENEID_R” (receptor NCBI Gene IDs) and 2. a meta information table describing the L-R-list. makeLRBasePackage function generates LRBase.XXX.eg.db like below. The gene identifier is limited as NCBI Gene ID for now.

example("makeLRBasePackage")
## 
## mkLRBP> if(interactive()){
## mkLRBP+     ## makeLRBasePackage enable users to construct
## mkLRBP+     ## user's own custom LRBase package
## mkLRBP+     data(FANTOM5)
## mkLRBP+     head(FANTOM5)
## mkLRBP+ 
## mkLRBP+     # We are also needed to prepare meta data as follows.
## mkLRBP+     data(metaFANTOM5)
## mkLRBP+     metaFANTOM5
## mkLRBP+ 
## mkLRBP+     ## sets up a temporary directory for this example
## mkLRBP+     ## (users won't need to do this step)
## mkLRBP+     tmp <- tempfile()
## mkLRBP+     dir.create(tmp)
## mkLRBP+ 
## mkLRBP+     ## makes an Organism package for human called Homo.sapiens
## mkLRBP+     makeLRBasePackage(pkgname = "FANTOM5.Hsa.eg.db",
## mkLRBP+         data = FANTOM5,
## mkLRBP+         metadata = metaFANTOM5,
## mkLRBP+         organism = "Homo sapiens",
## mkLRBP+         version = "0.99.0",
## mkLRBP+         maintainer = "Koki Tsuyuzaki <k.t.the-answer@hotmail.co.jp>",
## mkLRBP+         author = "Koki Tsuyuzaki",
## mkLRBP+         destDir = tmp,
## mkLRBP+         license="Artistic-2.0")
## mkLRBP+ }

Although any package name is acceptable, note that if the organism that the user summarized L-R-list is also described above (Table 1), same XXX-character is recommended. This is because the HTML report function descibed later identifies the XXX-character and if the XXX is corresponding to the 12 organisms, the gene annotation of generated HTML report will become rich.

2.3 scTensor (CCI-tensor construction, decomposition, and HTML reporting)

Combined with LRBase.XXX.eg.db-type package and user’s gene expression matrix of scRNA-Seq, scTensor detects CCIs and generates HTML reports for exploratory data inspection. The algorithm of scTensor is as follows.

Firstly, scTensor calculates the celltype-level mean vectors, searches the corresponding pair of genes in the rownames of the matrix, and extracted as vector.

Next, the celltype-level mean vectors of ligand expression and that of receptor expression are multiplied as outer product and converted to celltype \(\times\) celltype matrix. Here, the multiple matrices can be represented as a three-order “tensor” (Ligand-Cell * Receptor-Cell * L-R-Pair). scTensor decomposes the tensor into a small tensor (core tensor) and three factor matrices. Tensor decomposition is very similar to the matrix decomposition like PCA (principal component analysis). The core tensor is similar to eigenvalue of PCA; this means that how much the pattern is outstanding. Likewise, three matrices is similar to the PC scores/loadings of PCA; These represent which ligand-cell/receptor-cell/L-R-pair are informative. When the matrices have negative values, interpreting which direction (+/-) is important and which is not, is difficult and laboring task. That’s why, scTensor performs non-negative Tucker decomposition (NTD), which is non-negative version of tensor decomposition (c.f. nnTensor).

Finaly, the result of NTD is summarized as HTML report. Because most of plots are visualized by plotly package, the presise information of the plot can be interactively confirmed by user’s on-site web browser. The three factor matrices can be interactively viewed and which celltypes and which L-R-pairs are likely to be interacted each other. The mode-3 (LR-pair direction) sum of the core tensor is calculated and visualized as Ligand-Receptor Patterns. Detail of (Ligand-Cell, Receptor-Cell, L-R-pair) Patterns are also visualized.

2.3.1 Creating a SingleCellExperiment object

Here, we use the scRNA-Seq dataset of male gemline cells and somaticcells\(^{3}\) GSE86146 as demo data. For saving the package size, the number of genes are strictlly reduced by the standard of highlly variable genes with threshold of p-value is 1E-150 (c.f. Identifying highly variable genes). That’s why we won’t argue about the scientific discussion of the data here.

We assume that user have a scRNA-Seq data matrix containing expression count data summarised at the level of gene. First, we create a SingleCellExperiment object containing the data. The rows of the object correspond to features, and the columns correspond to cells. The gene identifier is limited as NCBI Gene ID for now.

To improve the interpretability of following HTML report, we highly recommend that user specifies the two-dimensional data of input data (e.g. PCA, t-SNE, or UMAP). Such information is easily specified by reducedDims function of SingleCellExperiment package and is saved to reducedDims slot of SingleCellExperiment object (1).

## Registered S3 methods overwritten by 'ggplot2':
##   method         from 
##   [.quosures     rlang
##   c.quosures     rlang
##   print.quosures rlang
## Registered S3 method overwritten by 'seriation':
##   method         from 
##   reorder.hclust gclus
## 
## 
## Registered S3 method overwritten by 'enrichplot':
##   method               from
##   fortify.enrichResult DOSE
## Loading required package: SummarizedExperiment
## Loading required package: GenomicRanges
## Loading required package: stats4
## Loading required package: BiocGenerics
## Loading required package: parallel
## 
## Attaching package: 'BiocGenerics'
## The following objects are masked from 'package:parallel':
## 
##     clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
##     clusterExport, clusterMap, parApply, parCapply, parLapply,
##     parLapplyLB, parRapply, parSapply, parSapplyLB
## The following objects are masked from 'package:stats':
## 
##     IQR, mad, sd, var, xtabs
## The following objects are masked from 'package:base':
## 
##     Filter, Find, Map, Position, Reduce, anyDuplicated, append,
##     as.data.frame, basename, cbind, colnames, dirname, do.call,
##     duplicated, eval, evalq, get, grep, grepl, intersect,
##     is.unsorted, lapply, mapply, match, mget, order, paste, pmax,
##     pmax.int, pmin, pmin.int, rank, rbind, rownames, sapply,
##     setdiff, sort, table, tapply, union, unique, unsplit, which,
##     which.max, which.min
## Loading required package: S4Vectors
## 
## Attaching package: 'S4Vectors'
## The following object is masked from 'package:base':
## 
##     expand.grid
## Loading required package: IRanges
## Loading required package: GenomeInfoDb
## Loading required package: Biobase
## Welcome to Bioconductor
## 
##     Vignettes contain introductory material; view with
##     'browseVignettes()'. To cite Bioconductor, see
##     'citation("Biobase")', and for packages 'citation("pkgname")'.
## Loading required package: DelayedArray
## Loading required package: matrixStats
## 
## Attaching package: 'matrixStats'
## The following objects are masked from 'package:Biobase':
## 
##     anyMissing, rowMedians
## Loading required package: BiocParallel
## 
## Attaching package: 'DelayedArray'
## The following objects are masked from 'package:matrixStats':
## 
##     colMaxs, colMins, colRanges, rowMaxs, rowMins, rowRanges
## The following objects are masked from 'package:base':
## 
##     aperm, apply, rowsum
## Loading required package: MeSHDbi
## 
## Attaching package: 'MeSHDbi'
## The following object is masked from 'package:utils':
## 
##     packageName
data(GermMale)
data(labelGermMale)
data(tsneGermMale)

sce <- SingleCellExperiment(assays=list(counts = GermMale))
reducedDims(sce) <- SimpleList(TSNE=tsneGermMale$Y)
plot(reducedDims(sce)[[1]], col=labelGermMale, pch=16, cex=2,
  xlab="Dim1", ylab="Dim2", main="Germline, Male, GSE86146")
legend("topleft", legend=c(paste0("FGC_", 1:3), paste0("Soma_", 1:4)),
  col=c("#9E0142", "#D53E4F", "#F46D43", "#ABDDA4", "#66C2A5", "#3288BD", "#5E4FA2"),
  pch=16)
Germline, Male, GSE86146

Figure 1: Germline, Male, GSE86146

Note that if you want to use scTensor framework against other species such as mouse or rat, load corresponding LRBase.XXX.eg.db and MeSH.XXX.eg.db packages.

For example, if your scData-Seq dataset is sampled from Mouse, load LRBase.Mmu.eg.db and MeSH.Mmu.eg.db instead of LRBase.Hsa.eg.db and MeSH.Hsa.eg.db.

2.3.2 Parameter setting : cellCellSetting

To perform the tensor decomposition and HTML report, the user is supposed to specify 1. LRBase.XXX.eg.db, 2. cell type vector of each cell, 3. color vector of each cell, to SingleCellExperiment object. The corresponding information is registed to metadata slot of SingleCellExperiment object by cellCellSetting function.

cellCellSetting(sce, LRBase.Hsa.eg.db, labelGermMale, names(labelGermMale))

2.3.3 CCI-tensor construction and decomposition : cellCellDecomp

After cellCellSetting, we can perform tensor decomposition by cellCellDecomp. Here the parameter ranks is specificed as dimension of core tensor. For example, c(2, 3, 4) means The data tensor is decomposed to 2 ligand-patterns, 3 receptor-patterns, and 4 L-R-pairs-patterns.

cellCellDecomp(sce, ranks=c(2,3,4))
## Input data matrix may contains 7 gene symbols because the name contains some alphabets.
## scTensor uses only NCBI Gene IDs for now.
## Here, the gene symbols are removed and remaining 235 NCBI Gene IDs are used for scTensor next step.
## 7 * 7 * 56 Tensor is created

Although the user have to specify the rank to perform cellCellDecomp, we implemented a simple rank estimation function based on the eigenvalues destribution of PCA in the matricised tensor in each mode.

cellCellRanks(sce)
## $selected
## [1]  3  2 10
## 
## $mode1
## [1] 15006.5633  8634.8288  6181.0632  2803.1701  1538.8738   671.7450
## [7]   264.4335
## 
## $mode2
## [1] 16939.51870  6620.09742  4111.89934   872.28806   433.18982   272.60479
## [7]    79.13515
## 
## $mode3
##  [1] 1.369633e+04 8.807127e+03 5.981817e+03 4.172650e+03 3.807329e+03
##  [6] 2.524074e+03 1.582185e+03 1.394027e+03 1.280911e+03 1.007136e+03
## [11] 7.986600e+02 6.772709e+02 5.798093e+02 5.150652e+02 4.996971e+02
## [16] 3.667198e+02 2.948783e+02 2.792895e+02 2.287873e+02 1.885954e+02
## [21] 1.258305e+02 1.069260e+02 8.208425e+01 7.288187e+01 6.309528e+01
## [26] 6.001535e+01 5.128678e+01 4.619312e+01 4.252876e+01 3.640923e+01
## [31] 2.844623e+01 2.638881e+01 2.412549e+01 2.255260e+01 1.820522e+01
## [36] 1.718524e+01 1.111941e+01 8.148261e+00 6.440623e+00 4.963200e+00
## [41] 4.313959e+00 2.911273e+00 2.567474e+00 1.565723e+00 1.254188e+00
## [46] 9.133304e-01 5.996757e-01 2.240669e-01 1.079619e-01

2.3.4 HTML Report : cellCellReport

If cellCellDecomp is properly finished, we can perform cellCellReport function to output the HTML report like belows. Please type example(cellCellReport) and the report will be generated in the temporary directory (it costs 5 to 10 minutes). After cellCellReport, multiple R markdown files, compiled HTML files, figures, and R binary file containing the result of analysis are saved to out.dir (Figure 2). For more details, open the index.html by your web browser. Combined with cloud strage servise such as Amazon Simple Storage Service (S3), it can be simple web application and multiple person like collaborators can confirm the same report simutaneously.

Figure2 : cellCellReport function of scTensor

Figure2 : cellCellReport function of scTensor

Session information

## R version 3.6.0 (2019-04-26)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.2 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.9-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.9-bioc/R/lib/libRlapack.so
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        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       
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] MeSH.Mmu.eg.db_1.12.0       LRBase.Mmu.eg.db_1.0.0     
##  [3] MeSH.Hsa.eg.db_1.12.0       MeSHDbi_1.20.0             
##  [5] SingleCellExperiment_1.6.0  SummarizedExperiment_1.14.0
##  [7] DelayedArray_0.10.0         BiocParallel_1.18.0        
##  [9] matrixStats_0.54.0          Biobase_2.44.0             
## [11] GenomicRanges_1.36.0        GenomeInfoDb_1.20.0        
## [13] IRanges_2.18.0              S4Vectors_0.22.0           
## [15] BiocGenerics_0.30.0         scTensor_1.0.0             
## [17] RSQLite_2.1.1               LRBase.Hsa.eg.db_1.0.0     
## [19] LRBaseDbi_1.2.0             BiocStyle_2.12.0           
## 
## loaded via a namespace (and not attached):
##   [1] rappdirs_0.3.1            rtracklayer_1.44.0       
##   [3] AnnotationForge_1.26.0    prabclus_2.2-7           
##   [5] tidyr_0.8.3               ggplot2_3.1.1            
##   [7] acepack_1.4.1             bit64_0.9-7              
##   [9] knitr_1.22                data.table_1.12.2        
##  [11] rpart_4.1-15              RCurl_1.95-4.12          
##  [13] AnnotationFilter_1.8.0    GenomicFeatures_1.36.0   
##  [15] cowplot_0.9.4             europepmc_0.3            
##  [17] bit_1.1-14                enrichplot_1.4.0         
##  [19] webshot_0.5.1             xml2_1.2.0               
##  [21] assertthat_0.2.1          viridis_0.5.1            
##  [23] xfun_0.6                  hms_0.4.2                
##  [25] evaluate_0.13             TSP_1.1-6                
##  [27] DEoptimR_1.0-8            progress_1.2.0           
##  [29] caTools_1.17.1.2          dendextend_1.10.0        
##  [31] Rgraphviz_2.28.0          igraph_1.2.4.1           
##  [33] DBI_1.0.0                 htmlwidgets_1.3          
##  [35] MeSH.db_1.12.0            purrr_0.3.2              
##  [37] dplyr_0.8.0.1             backports_1.1.4          
##  [39] bookdown_0.9              trimcluster_0.1-2.1      
##  [41] annotate_1.62.0           biomaRt_2.40.0           
##  [43] ensembldb_2.8.0           abind_1.4-5              
##  [45] ggforce_0.2.2             Gviz_1.28.0              
##  [47] triebeard_0.3.0           BSgenome_1.52.0          
##  [49] robustbase_0.93-4         checkmate_1.9.1          
##  [51] GenomicAlignments_1.20.0  gclus_1.3.2              
##  [53] fdrtool_1.2.15            prettyunits_1.0.2        
##  [55] mclust_5.4.3              cluster_2.0.9            
##  [57] DOSE_3.10.0               dotCall64_1.0-0          
##  [59] lazyeval_0.2.2            crayon_1.3.4             
##  [61] genefilter_1.66.0         pkgconfig_2.0.2          
##  [63] tweenr_1.0.1              ProtGenerics_1.16.0      
##  [65] seriation_1.2-3           nnet_7.3-12              
##  [67] rlang_0.3.4               diptest_0.75-7           
##  [69] meshr_1.20.0              registry_0.5-1           
##  [71] MeSH.PCR.db_1.12.0        rTensor_1.4              
##  [73] GOstats_2.50.0            dichromat_2.0-0          
##  [75] polyclip_1.10-0           graph_1.62.0             
##  [77] Matrix_1.2-17             urltools_1.7.3           
##  [79] base64enc_0.1-3           whisker_0.3-2            
##  [81] ggridges_0.5.1            viridisLite_0.3.0        
##  [83] MeSH.AOR.db_1.12.0        bitops_1.0-6             
##  [85] KernSmooth_2.23-15        spam_2.2-2               
##  [87] MeSH.Bsu.168.eg.db_1.12.0 Biostrings_2.52.0        
##  [89] blob_1.1.1                stringr_1.4.0            
##  [91] qvalue_2.16.0             nnTensor_0.99.4          
##  [93] gridGraphics_0.3-0        reactome.db_1.68.0       
##  [95] scales_1.0.0              memoise_1.1.0            
##  [97] graphite_1.30.0           GSEABase_1.46.0          
##  [99] magrittr_1.5              plyr_1.8.4               
## [101] gplots_3.0.1.1            gdata_2.18.0             
## [103] zlibbioc_1.30.0           compiler_3.6.0           
## [105] RColorBrewer_1.1-2        plotrix_3.7-5            
## [107] Rsamtools_2.0.0           XVector_0.24.0           
## [109] Category_2.50.0           MeSH.Aca.eg.db_1.12.0    
## [111] htmlTable_1.13.1          Formula_1.2-3            
## [113] MASS_7.3-51.4             tidyselect_0.2.5         
## [115] stringi_1.4.3             highr_0.8                
## [117] MeSH.Syn.eg.db_1.12.0     yaml_2.2.0               
## [119] GOSemSim_2.10.0           latticeExtra_0.6-28      
## [121] ggrepel_0.8.0             grid_3.6.0               
## [123] VariantAnnotation_1.30.0  fastmatch_1.1-0          
## [125] tools_3.6.0               rstudioapi_0.10          
## [127] foreach_1.4.4             foreign_0.8-71           
## [129] tagcloud_0.6              outliers_0.14            
## [131] gridExtra_2.3             farver_1.1.0             
## [133] ggraph_1.0.2              digest_0.6.18            
## [135] rvcheck_0.1.3             BiocManager_1.30.4       
## [137] fpc_2.1-11.2              Rcpp_1.0.1               
## [139] org.Hs.eg.db_3.8.2        httr_1.4.0               
## [141] cummeRbund_2.26.0         AnnotationDbi_1.46.0     
## [143] biovizBase_1.32.0         kernlab_0.9-27           
## [145] colorspace_1.4-1          XML_3.98-1.19            
## [147] splines_3.6.0             fields_9.7               
## [149] RBGL_1.60.0               ggplotify_0.0.3          
## [151] flexmix_2.3-15            plotly_4.9.0             
## [153] xtable_1.8-4              jsonlite_1.6             
## [155] heatmaply_0.15.2          UpSetR_1.3.3             
## [157] modeltools_0.2-22         R6_2.4.0                 
## [159] Hmisc_4.2-0               pillar_1.3.1             
## [161] htmltools_0.3.6           glue_1.3.1               
## [163] class_7.3-15              codetools_0.2-16         
## [165] maps_3.3.0                fgsea_1.10.0             
## [167] mvtnorm_1.0-10            lattice_0.20-38          
## [169] tibble_2.1.1              curl_3.3                 
## [171] gtools_3.8.1              ReactomePA_1.28.0        
## [173] misc3d_0.8-4              GO.db_3.8.2              
## [175] survival_2.44-1.1         rmarkdown_1.12           
## [177] munsell_0.5.0             DO.db_2.9                
## [179] GenomeInfoDbData_1.2.1    plot3D_1.1.1             
## [181] iterators_1.0.10          reshape2_1.4.3           
## [183] gtable_0.3.0