Example gene set analysis: The case of HDAC inhibitor SAHA
Mark Ziemann & Kaumadi Wijesooriya
2022-01-25
Source: https://github.com/markziemann/combined_enrichment
Intro
Here we are performing an analysis of some gene expression data to demonstrate the difference between ORA and FCS methods and to highlight the differences caused by improper background gene set use.
The dataset being used is SRP096178 and we are comparing the cells grown in normal condition (control) to those grown with addition of SAHA (case).
Data are obtained from http://dee2.io/
suppressPackageStartupMessages({
library("getDEE2")
library("DESeq2")
library("clusterProfiler")
library("mitch")
library("kableExtra")
library("eulerr")
})Get expression data
I’m using some RNA-seq data looking at the effect of SAHA on HAEC cells.
name="SRP037718"
mdat<-getDEE2Metadata("hsapiens")
samplesheet <- mdat[grep("SRP037718",mdat$SRP_accession),]
samplesheet<-samplesheet[order(samplesheet$SRR_accession),]
samplesheet$trt<-as.factor(c(1,1,1,0,0,0))
s1 <- samplesheet
s1 %>% kbl(caption = "sample sheet") %>% kable_paper("hover", full_width = F)| SRR_accession | QC_summary | SRX_accession | SRS_accession | SRP_accession | Sample_name | GEO_series | Library_name | trt | |
|---|---|---|---|---|---|---|---|---|---|
| 238589 | SRR1168225 | PASS | SRX469930 | SRS557162 | SRP037718 | GSM1326469 | GSE37378 | 1 | |
| 238590 | SRR1168226 | PASS | SRX469931 | SRS557163 | SRP037718 | GSM1326470 | GSE37378 | 1 | |
| 238591 | SRR1168227 | PASS | SRX469932 | SRS557164 | SRP037718 | GSM1326471 | GSE37378 | 1 | |
| 238592 | SRR1168228 | PASS | SRX469933 | SRS557165 | SRP037718 | GSM1326472 | GSE37378 | 0 | |
| 238593 | SRR1168229 | PASS | SRX469934 | SRS557166 | SRP037718 | GSM1326473 | GSE37378 | 0 | |
| 238594 | SRR1168230 | PASS | SRX469935 | SRS557167 | SRP037718 | GSM1326474 | GSE37378 | 0 |
w<-getDEE2("hsapiens",s1$SRR_accession,metadata=mdat,legacy = TRUE)## For more information about DEE2 QC metrics, visit
## https://github.com/markziemann/dee2/blob/master/qc/qc_metrics.mdx<-Tx2Gene(w)
x<-x$Tx2Gene
# save the genetable for later
gt<-w$GeneInfo[,1,drop=FALSE]
gt$accession<-rownames(gt)
# counts
x1<-x[,which(colnames(x) %in% samplesheet$SRR_accession)]Here show the number of genes in the annotation set, and those detected above the detection threshold.
# filter out lowly expressed genes
x1<-x1[which(rowSums(x1)/ncol(x1)>=(10)),]
nrow(x)## [1] 39297nrow(x1)## [1] 15477Now multidimensional scaling (MDS) plot to show the correlation between the datasets. If the control and case datasets are clustered separately, then it is likely that there will be many differentially expressed genes with FDR<0.05.
plot(cmdscale(dist(t(x1))), xlab="Coordinate 1", ylab="Coordinate 2", pch=19, col=s1$trt, main="MDS")
Differential expression
Now run DESeq2 for control vs case.
y <- DESeqDataSetFromMatrix(countData = round(x1), colData = s1, design = ~ trt)## converting counts to integer modey <- DESeq(y)## estimating size factors## estimating dispersions## gene-wise dispersion estimates## mean-dispersion relationship## final dispersion estimates## fitting model and testingde <- results(y)
de<-as.data.frame(de[order(de$pvalue),])
rownames(de)<-sapply(strsplit(rownames(de),"\\."),"[[",1)
head(de) %>% kbl() %>% kable_paper("hover", full_width = F)| baseMean | log2FoldChange | lfcSE | stat | pvalue | padj | |
|---|---|---|---|---|---|---|
| ENSG00000099250 | 6148.277 | -4.163885 | 0.1204942 | -34.55671 | 0 | 0 |
| ENSG00000187193 | 3478.958 | 4.219910 | 0.1248838 | 33.79071 | 0 | 0 |
| ENSG00000126785 | 2075.476 | -3.356182 | 0.1085961 | -30.90517 | 0 | 0 |
| ENSG00000166741 | 1046.262 | -5.101344 | 0.1661858 | -30.69662 | 0 | 0 |
| ENSG00000102010 | 1153.368 | -5.552842 | 0.1817905 | -30.54528 | 0 | 0 |
| ENSG00000159167 | 1898.444 | 7.629754 | 0.2499001 | 30.53122 | 0 | 0 |
Now let’s have a look at some of the charts showing differential expression. In particular, an MA plot and volcano plot.
maplot <- function(de,contrast_name) {
sig <-subset(de, padj < 0.05 )
up <-rownames(subset(de, padj < 0.05 & log2FoldChange > 0))
dn <-rownames(subset(de, padj < 0.05 & log2FoldChange < 0))
GENESUP <- length(up)
GENESDN <- length(dn)
DET=nrow(de)
SUBHEADER = paste(GENESUP, "up, ", GENESDN, "down", DET, "detected")
ns <-subset(de, padj > 0.05 )
plot(log2(de$baseMean),de$log2FoldChange,
xlab="log2 basemean", ylab="log2 foldchange",
pch=19, cex=0.5, col="dark gray",
main=contrast_name, cex.main=0.7)
points(log2(sig$baseMean),sig$log2FoldChange,
pch=19, cex=0.5, col="red")
mtext(SUBHEADER,cex = 0.7)
}
make_volcano <- function(de,name) {
sig <- subset(de,padj<0.05)
N_SIG=nrow(sig)
N_UP=nrow(subset(sig,log2FoldChange>0))
N_DN=nrow(subset(sig,log2FoldChange<0))
DET=nrow(de)
HEADER=paste(N_SIG,"@5%FDR,", N_UP, "up", N_DN, "dn", DET, "detected")
plot(de$log2FoldChange,-log10(de$padj),cex=0.5,pch=19,col="darkgray",
main=name, xlab="log2 FC", ylab="-log10 pval", xlim=c(-6,6))
mtext(HEADER)
grid()
points(sig$log2FoldChange,-log10(sig$padj),cex=0.5,pch=19,col="red")
}
maplot(de,name)
make_volcano(de,name)
Gene sets from Reactome
In order to perform gene set analysis, we need some gene sets.
if (! file.exists("ReactomePathways.gmt")) {
download.file("https://reactome.org/download/current/ReactomePathways.gmt.zip",
destfile="ReactomePathways.gmt.zip")
unzip("ReactomePathways.gmt.zip")
}
genesets<-gmt_import("ReactomePathways.gmt")FCS with Mitch
Mitch uses rank-ANOVA statistics for enrichment detection.
Here I’m using the standard approach
m <- mitch_import(de,DEtype = "DEseq2", geneTable = gt)## The input is a single dataframe; one contrast only. Converting
## it to a list for you.## Note: Mean no. genes in input = 15477## Note: no. genes in output = 14488## Note: estimated proportion of input genes in output = 0.936msep <- mitch_calc(m,genesets = genesets)## Note: When prioritising by significance (ie: small
## p-values), large effect sizes might be missed.ms_up <- subset(msep$enrichment_result,p.adjustANOVA<0.05 & s.dist > 0)[,1]
ms_dn <- subset(msep$enrichment_result,p.adjustANOVA<0.05 & s.dist < 0)[,1]
message(paste("Number of up-regulated pathways:",length(ms_up) ))## Number of up-regulated pathways: 49message(paste("Number of down-regulated pathways:",length(ms_dn) ))## Number of down-regulated pathways: 266head(msep$enrichment_result,10) #%>% kbl() %>% kable_paper("hover", full_width = F)## set
## 624 Metabolism of RNA
## 1283 Translation
## 149 Cell Cycle
## 151 Cell Cycle, Mitotic
## 1365 rRNA processing
## 150 Cell Cycle Checkpoints
## 636 Metabolism of proteins
## 1367 rRNA processing in the nucleus and cytosol
## 728 Nonsense Mediated Decay (NMD) enhanced by the Exon Junction Complex (EJC)
## 730 Nonsense-Mediated Decay (NMD)
## setSize pANOVA s.dist p.adjustANOVA
## 624 661 8.016722e-34 -0.2779789 1.101498e-30
## 1283 273 2.072635e-22 -0.3429991 1.423900e-19
## 149 599 5.288637e-22 -0.2319664 2.422196e-19
## 151 484 2.718405e-21 -0.2524928 9.337720e-19
## 1365 208 3.176411e-18 -0.3505784 8.728776e-16
## 150 251 7.305882e-18 -0.3161425 1.673047e-15
## 636 1562 1.816513e-17 -0.1313830 3.565555e-15
## 1367 187 2.818720e-16 -0.3472750 4.841152e-14
## 728 113 4.289967e-16 -0.4428306 5.894415e-14
## 730 113 4.289967e-16 -0.4428306 5.894415e-14#mitch_report(msep,outfile="mitch_separate.html",overwrite=TRUE)Here I’m using the combined approach.
mcom <- mitch_calc(abs(m),genesets = genesets)## Note: When prioritising by significance (ie: small
## p-values), large effect sizes might be missed.mc_up <- subset(mcom$enrichment_result,p.adjustANOVA<0.05 & s.dist > 0)[,1]
mc_dn <- subset(mcom$enrichment_result,p.adjustANOVA<0.05 & s.dist < 0)[,1]
message(paste("Number of up-regulated pathways:",length(mc_up) ))## Number of up-regulated pathways: 61message(paste("Number of down-regulated pathways:",length(mc_dn) ))## Number of down-regulated pathways: 23head(mcom$enrichment_result,10) #%>% kbl() %>% kable_paper("hover", full_width = F)## set setSize pANOVA
## 1074 Signal Transduction 1887 2.279744e-14
## 487 Hemostasis 433 1.936932e-12
## 276 Developmental Biology 712 5.807837e-08
## 718 Neuronal System 244 5.492600e-07
## 807 Platelet activation, signaling and aggregation 196 1.465562e-06
## 510 Immune System 1399 1.535736e-06
## 1148 Signaling by Receptor Tyrosine Kinases 423 5.199623e-06
## 1114 Signaling by GPCR 352 9.722372e-06
## 92 Axon guidance 462 1.401792e-05
## 902 RNA Polymerase II Pre-transcription Events 79 1.911207e-05
## s.dist p.adjustANOVA
## 1074 0.10869146 3.132369e-11
## 487 0.19813369 1.330672e-09
## 276 0.12031733 2.659989e-05
## 718 0.18662445 1.886708e-04
## 807 0.19989776 3.516836e-04
## 510 0.07803283 3.516836e-04
## 1148 0.12978463 1.020612e-03
## 1114 0.13776600 1.669817e-03
## 92 0.11854947 2.140068e-03
## 902 -0.27838877 2.625998e-03#mitch_report(mcom,outfile="mitch_combined.html",overwrite=TRUE)Let’s look at the significant ones based on the combined analysis. There weren’t many gene sets classe d as significant. Let’s see how many have a direction agnostic enrichment score which is larger in magnitude than the direction informed enrichment score. There are only 11 such sets which would benefit from such combined analysis.
Euler diagram of the significant pathways found with each approach.
l0 <- list("sep up"=ms_up,"sep dn"=ms_dn,"comb up"=mc_up,"comb dn"=mc_dn)
par(cex.main=0.5)
plot(euler(l0),quantities = TRUE, edges = "gray", main="FCS: combined vs separated")
length(ms_up)## [1] 49length(ms_dn)## [1] 266length(ms_up)+length(ms_dn)## [1] 315length(mc_up)## [1] 61length(mc_dn)## [1] 23length(mc_up)+length(mc_dn)## [1] 84( length(ms_up)+length(ms_dn) ) / ( length(mc_up)+length(mc_dn) )## [1] 3.75List gene sets which are specific to each approach.
ms <- c(ms_up,ms_dn)
# in sep but not comb
setdiff(ms,mc_up)## [1] "Dopamine Neurotransmitter Release Cycle"
## [2] "Neurotransmitter release cycle"
## [3] "Ion channel transport"
## [4] "Acetylcholine Neurotransmitter Release Cycle"
## [5] "Sensory Perception"
## [6] "Norepinephrine Neurotransmitter Release Cycle"
## [7] "Cilium Assembly"
## [8] "Metabolism of lipids"
## [9] "Sensory processing of sound by outer hair cells of the cochlea"
## [10] "Neurexins and neuroligins"
## [11] "GABA synthesis, release, reuptake and degradation"
## [12] "Assembly and cell surface presentation of NMDA receptors"
## [13] "Ion homeostasis"
## [14] "Protein-protein interactions at synapses"
## [15] "Glutamate Neurotransmitter Release Cycle"
## [16] "Branched-chain amino acid catabolism"
## [17] "Negative regulation of TCF-dependent signaling by WNT ligand antagonists"
## [18] "Phospholipid metabolism"
## [19] "NGF-stimulated transcription"
## [20] "Inositol phosphate metabolism"
## [21] "Heparan sulfate/heparin (HS-GAG) metabolism"
## [22] "Synthesis of IP2, IP, and Ins in the cytosol"
## [23] "A tetrasaccharide linker sequence is required for GAG synthesis"
## [24] "Glycerophospholipid biosynthesis"
## [25] "Nuclear Events (kinase and transcription factor activation)"
## [26] "HS-GAG biosynthesis"
## [27] "Transport of small molecules"
## [28] "Signaling by NTRKs"
## [29] "Defective B3GAT3 causes JDSSDHD"
## [30] "Regulation of IFNG signaling"
## [31] "Defective EXT1 causes exostoses 1, TRPS2 and CHDS"
## [32] "Defective EXT2 causes exostoses 2"
## [33] "Defective B4GALT7 causes EDS, progeroid type"
## [34] "GPCR ligand binding"
## [35] "Insulin receptor recycling"
## [36] "Plasma lipoprotein assembly, remodeling, and clearance"
## [37] "Antigen activates B Cell Receptor (BCR) leading to generation of second messengers"
## [38] "Signaling by NTRK1 (TRKA)"
## [39] "Metabolism of RNA"
## [40] "Translation"
## [41] "Cell Cycle"
## [42] "rRNA processing"
## [43] "Cell Cycle Checkpoints"
## [44] "Metabolism of proteins"
## [45] "rRNA processing in the nucleus and cytosol"
## [46] "Nonsense Mediated Decay (NMD) enhanced by the Exon Junction Complex (EJC)"
## [47] "Nonsense-Mediated Decay (NMD)"
## [48] "SRP-dependent cotranslational protein targeting to membrane"
## [49] "Mitotic Metaphase and Anaphase"
## [50] "Cap-dependent Translation Initiation"
## [51] "Eukaryotic Translation Initiation"
## [52] "L13a-mediated translational silencing of Ceruloplasmin expression"
## [53] "Mitotic Anaphase"
## [54] "GTP hydrolysis and joining of the 60S ribosomal subunit"
## [55] "Nonsense Mediated Decay (NMD) independent of the Exon Junction Complex (EJC)"
## [56] "Major pathway of rRNA processing in the nucleolus and cytosol"
## [57] "Formation of a pool of free 40S subunits"
## [58] "M Phase"
## [59] "Signaling by ROBO receptors"
## [60] "Regulation of expression of SLITs and ROBOs"
## [61] "Viral mRNA Translation"
## [62] "Response of EIF2AK4 (GCN2) to amino acid deficiency"
## [63] "Peptide chain elongation"
## [64] "Eukaryotic Translation Termination"
## [65] "Eukaryotic Translation Elongation"
## [66] "Cellular responses to stress"
## [67] "Influenza Infection"
## [68] "Selenocysteine synthesis"
## [69] "Selenoamino acid metabolism"
## [70] "Separation of Sister Chromatids"
## [71] "Influenza Viral RNA Transcription and Replication"
## [72] "DNA Replication"
## [73] "Processing of Capped Intron-Containing Pre-mRNA"
## [74] "S Phase"
## [75] "G2/M Checkpoints"
## [76] "Cellular responses to stimuli"
## [77] "Mitotic G1 phase and G1/S transition"
## [78] "Synthesis of DNA"
## [79] "G1/S Transition"
## [80] "mRNA Splicing - Major Pathway"
## [81] "Infectious disease"
## [82] "HATs acetylate histones"
## [83] "DNA Replication Pre-Initiation"
## [84] "mRNA Splicing"
## [85] "Gene expression (Transcription)"
## [86] "Resolution of Sister Chromatid Cohesion"
## [87] "APC/C-mediated degradation of cell cycle proteins"
## [88] "Regulation of mitotic cell cycle"
## [89] "Mitotic Spindle Checkpoint"
## [90] "Regulation of APC/C activators between G1/S and early anaphase"
## [91] "Ribosomal scanning and start codon recognition"
## [92] "Class I MHC mediated antigen processing & presentation"
## [93] "Mitotic Prophase"
## [94] "Switching of origins to a post-replicative state"
## [95] "Cellular response to starvation"
## [96] "Translation initiation complex formation"
## [97] "Activation of the mRNA upon binding of the cap-binding complex and eIFs, and subsequent binding to 43S"
## [98] "Activation of APC/C and APC/C:Cdc20 mediated degradation of mitotic proteins"
## [99] "Assembly of the pre-replicative complex"
## [100] "Host Interactions of HIV factors"
## [101] "Mitochondrial translation initiation"
## [102] "Mitochondrial translation termination"
## [103] "APC/C:Cdc20 mediated degradation of mitotic proteins"
## [104] "Cdc20:Phospho-APC/C mediated degradation of Cyclin A"
## [105] "Formation of the ternary complex, and subsequently, the 43S complex"
## [106] "APC:Cdc20 mediated degradation of cell cycle proteins prior to satisfation of the cell cycle checkpoint"
## [107] "EML4 and NUDC in mitotic spindle formation"
## [108] "Mitochondrial translation elongation"
## [109] "Amplification of signal from unattached kinetochores via a MAD2 inhibitory signal"
## [110] "Amplification of signal from the kinetochores"
## [111] "Orc1 removal from chromatin"
## [112] "Transcriptional regulation by RUNX3"
## [113] "RNA Polymerase II Transcription"
## [114] "APC/C:Cdh1 mediated degradation of Cdc20 and other APC/C:Cdh1 targeted proteins in late mitosis/early G1"
## [115] "RNA Polymerase II Transcription Termination"
## [116] "Signaling by NOTCH4"
## [117] "RHO GTPases Activate Formins"
## [118] "Antigen processing-Cross presentation"
## [119] "Mitochondrial translation"
## [120] "Transport of Mature Transcript to Cytoplasm"
## [121] "Transcriptional regulation by RUNX2"
## [122] "The role of GTSE1 in G2/M progression after G2 checkpoint"
## [123] "ER-Phagosome pathway"
## [124] "RUNX1 regulates transcription of genes involved in differentiation of HSCs"
## [125] "Post-translational protein modification"
## [126] "tRNA processing"
## [127] "DNA Repair"
## [128] "Autodegradation of Cdh1 by Cdh1:APC/C"
## [129] "CDK-mediated phosphorylation and removal of Cdc6"
## [130] "APC/C:Cdc20 mediated degradation of Securin"
## [131] "Transport of Mature mRNA derived from an Intron-Containing Transcript"
## [132] "mRNA 3'-end processing"
## [133] "Antigen processing: Ubiquitination & Proteasome degradation"
## [134] "Adaptive Immune System"
## [135] "Stabilization of p53"
## [136] "Cyclin E associated events during G1/S transition"
## [137] "Deubiquitination"
## [138] "Transcriptional Regulation by TP53"
## [139] "DNA Double-Strand Break Repair"
## [140] "Transcriptional regulation by RUNX1"
## [141] "Generic Transcription Pathway"
## [142] "G1/S DNA Damage Checkpoints"
## [143] "Cyclin A:Cdk2-associated events at S phase entry"
## [144] "Mitotic Prometaphase"
## [145] "Polo-like kinase mediated events"
## [146] "Regulation of RUNX3 expression and activity"
## [147] "AUF1 (hnRNP D0) binds and destabilizes mRNA"
## [148] "SCF(Skp2)-mediated degradation of p27/p21"
## [149] "Negative regulation of NOTCH4 signaling"
## [150] "Ubiquitin Mediated Degradation of Phosphorylated Cdc25A"
## [151] "p53-Independent DNA Damage Response"
## [152] "p53-Independent G1/S DNA damage checkpoint"
## [153] "Activation of the pre-replicative complex"
## [154] "Homology Directed Repair"
## [155] "Ubiquitin-dependent degradation of Cyclin D"
## [156] "Neddylation"
## [157] "Vif-mediated degradation of APOBEC3G"
## [158] "Mitotic G2-G2/M phases"
## [159] "Nuclear Envelope Breakdown"
## [160] "Nuclear Envelope (NE) Reassembly"
## [161] "HIV Infection"
## [162] "Autodegradation of the E3 ubiquitin ligase COP1"
## [163] "Mitochondrial protein import"
## [164] "ROS sensing by NFE2L2"
## [165] "Cross-presentation of soluble exogenous antigens (endosomes)"
## [166] "Downstream TCR signaling"
## [167] "RMTs methylate histone arginines"
## [168] "SCF-beta-TrCP mediated degradation of Emi1"
## [169] "Regulation of activated PAK-2p34 by proteasome mediated degradation"
## [170] "FBXL7 down-regulates AURKA during mitotic entry and in early mitosis"
## [171] "PTEN Regulation"
## [172] "Initiation of Nuclear Envelope (NE) Reformation"
## [173] "MAPK6/MAPK4 signaling"
## [174] "Regulation of mRNA stability by proteins that bind AU-rich elements"
## [175] "Degradation of beta-catenin by the destruction complex"
## [176] "Ub-specific processing proteases"
## [177] "p53-Dependent G1 DNA Damage Response"
## [178] "p53-Dependent G1/S DNA damage checkpoint"
## [179] "G2/M Transition"
## [180] "Defective CFTR causes cystic fibrosis"
## [181] "SUMOylation"
## [182] "Interactions of Rev with host cellular proteins"
## [183] "UCH proteinases"
## [184] "SUMO E3 ligases SUMOylate target proteins"
## [185] "GLI3 is processed to GLI3R by the proteasome"
## [186] "FCERI mediated NF-kB activation"
## [187] "DNA strand elongation"
## [188] "Oxygen-dependent proline hydroxylation of Hypoxia-inducible Factor Alpha"
## [189] "Disorders of transmembrane transporters"
## [190] "ABC transporter disorders"
## [191] "Transcriptional Regulation by VENTX"
## [192] "Degradation of DVL"
## [193] "CDT1 association with the CDC6:ORC:origin complex"
## [194] "Downstream signaling events of B Cell Receptor (BCR)"
## [195] "Activation of ATR in response to replication stress"
## [196] "G2/M DNA damage checkpoint"
## [197] "Regulation of Apoptosis"
## [198] "RHO GTPase Effectors"
## [199] "Degradation of AXIN"
## [200] "HDR through Homologous Recombination (HRR) or Single Strand Annealing (SSA)"
## [201] "Regulation of HMOX1 expression and activity"
## [202] "Activation of NF-kappaB in B cells"
## [203] "Regulation of PTEN stability and activity"
## [204] "Cellular response to hypoxia"
## [205] "Vpu mediated degradation of CD4"
## [206] "Degradation of GLI2 by the proteasome"
## [207] "Regulation of RAS by GAPs"
## [208] "TCR signaling"
## [209] "Degradation of GLI1 by the proteasome"
## [210] "Hedgehog 'on' state"
## [211] "Chk1/Chk2(Cds1) mediated inactivation of Cyclin B:Cdk1 complex"
## [212] "Phosphorylation of the APC/C"
## [213] "PKMTs methylate histone lysines"
## [214] "Golgi Cisternae Pericentriolar Stack Reorganization"
## [215] "Hh mutants are degraded by ERAD"
## [216] "Rev-mediated nuclear export of HIV RNA"
## [217] "Nuclear import of Rev protein"
## [218] "Nuclear Pore Complex (NPC) Disassembly"
## [219] "Signaling by NOTCH"
## [220] "Hh mutants abrogate ligand secretion"
## [221] "Inactivation of APC/C via direct inhibition of the APC/C complex"
## [222] "Inhibition of the proteolytic activity of APC/C required for the onset of anaphase by mitotic spindle checkpoint components"
## [223] "SUMOylation of DNA replication proteins"
## [224] "Regulation of ornithine decarboxylase (ODC)"
## [225] "Interleukin-1 family signaling"
## [226] "Global Genome Nucleotide Excision Repair (GG-NER)"
## [227] "NS1 Mediated Effects on Host Pathways"
## [228] "Cytoprotection by HMOX1"
## [229] "MASTL Facilitates Mitotic Progression"
## [230] "TNFR2 non-canonical NF-kB pathway"
## [231] "Dectin-1 mediated noncanonical NF-kB signaling"
## [232] "Cellular response to chemical stress"
## [233] "NIK-->noncanonical NF-kB signaling"
## [234] "Signaling by Interleukins"
## [235] "HDR through Homologous Recombination (HRR)"
## [236] "HCMV Infection"
## [237] "Regulation of RUNX2 expression and activity"
## [238] "Metabolism of non-coding RNA"
## [239] "snRNP Assembly"
## [240] "Deposition of new CENPA-containing nucleosomes at the centromere"
## [241] "Nucleosome assembly"
## [242] "Interleukin-1 signaling"
## [243] "Pre-NOTCH Transcription and Translation"
## [244] "rRNA modification in the nucleus and cytosol"
## [245] "rRNA processing in the mitochondrion"
## [246] "APC-Cdc20 mediated degradation of Nek2A"
## [247] "Regulation of RUNX1 Expression and Activity"
## [248] "Regulation of TP53 Activity through Phosphorylation"
## [249] "IRAK4 deficiency (TLR2/4)"
## [250] "Cellular response to heat stress"
## [251] "Diseases associated with N-glycosylation of proteins"
## [252] "ESR-mediated signaling"
## [253] "Condensation of Prophase Chromosomes"
## [254] "Asymmetric localization of PCP proteins"
## [255] "Nucleotide Excision Repair"
## [256] "APC/C:Cdc20 mediated degradation of Cyclin B"
## [257] "ABC-family proteins mediated transport"
## [258] "Nicotinamide salvaging"
## [259] "C-type lectin receptors (CLRs)"
## [260] "Processing of DNA double-strand break ends"
## [261] "Postmitotic nuclear pore complex (NPC) reformation"
## [262] "Aberrant regulation of mitotic cell cycle due to RB1 defects"
## [263] "TCF dependent signaling in response to WNT"
## [264] "tRNA processing in the mitochondrion"
## [265] "Hedgehog ligand biogenesis"
## [266] "tRNA modification in the nucleus and cytosol"
## [267] "NEP/NS2 Interacts with the Cellular Export Machinery"
## [268] "Regulation of TP53 Activity"
## [269] "Regulated Necrosis"
## [270] "HCMV Early Events"
## [271] "Respiratory electron transport"
## [272] "Dual Incision in GG-NER"
## [273] "TP53 Regulates Transcription of Cell Cycle Genes"
## [274] "Nicotinate metabolism"
## [275] "Diseases of mitotic cell cycle"
## [276] "Metabolism of polyamines"
## [277] "Transcriptional Regulation by E2F6"
## [278] "Homologous DNA Pairing and Strand Exchange"
## [279] "Gene and protein expression by JAK-STAT signaling after Interleukin-12 stimulation"
## [280] "Estrogen-dependent gene expression"
## [281] "Transport of Ribonucleoproteins into the Host Nucleus"
## [282] "Processing of Capped Intronless Pre-mRNA"
## [283] "Diseases of DNA repair"
## [284] "CLEC7A (Dectin-1) signaling"
## [285] "Export of Viral Ribonucleoproteins from Nucleus"
## [286] "Metabolism of nucleotides"
## [287] "Oncogene Induced Senescence"
## [288] "Complex I biogenesis"
## [289] "Cyclin A/B1/B2 associated events during G2/M transition"
## [290] "Condensation of Prometaphase Chromosomes"# in comb but not sep
setdiff(mc_up,ms)## [1] "Signal Transduction"
## [2] "Hemostasis"
## [3] "Platelet activation, signaling and aggregation"
## [4] "Signaling by Receptor Tyrosine Kinases"
## [5] "Signaling by GPCR"
## [6] "Cell junction organization"
## [7] "GPCR downstream signalling"
## [8] "G alpha (s) signalling events"
## [9] "Response to elevated platelet cytosolic Ca2+"
## [10] "ADORA2B mediated anti-inflammatory cytokines production"
## [11] "Platelet degranulation"
## [12] "Extracellular matrix organization"
## [13] "Regulation of Insulin-like Growth Factor (IGF) transport and uptake by Insulin-like Growth Factor Binding Proteins (IGFBPs)"
## [14] "Leishmania infection"
## [15] "Class B/2 (Secretin family receptors)"
## [16] "L1CAM interactions"
## [17] "Post-translational protein phosphorylation"
## [18] "Chaperonin-mediated protein folding"
## [19] "Diseases of signal transduction by growth factor receptors and second messengers"
## [20] "EPH-Ephrin signaling"
## [21] "Cell-extracellular matrix interactions"
## [22] "Factors involved in megakaryocyte development and platelet production"
## [23] "Dissolution of Fibrin Clot"
## [24] "G alpha (i) signalling events"
## [25] "FCGR3A-mediated phagocytosis"
## [26] "Leishmania phagocytosis"
## [27] "Parasite infection"
## [28] "Netrin-1 signaling"
## [29] "Cell-cell junction organization"
## [30] "G alpha (z) signalling events"
## [31] "Vasopressin regulates renal water homeostasis via Aquaporins"
## [32] "DCC mediated attractive signaling"
## [33] "Aquaporin-mediated transport"
## [34] "Cell-Cell communication"
## [35] "RHO GTPases Activate WASPs and WAVEs"
## [36] "Intrinsic Pathway for Apoptosis"# intersection
intersect(mc_up,ms)## [1] "Developmental Biology"
## [2] "Neuronal System"
## [3] "Immune System"
## [4] "Axon guidance"
## [5] "Nervous system development"
## [6] "Potassium Channels"
## [7] "Muscle contraction"
## [8] "Disease"
## [9] "Transmission across Chemical Synapses"
## [10] "Sensory processing of sound"
## [11] "Programmed Cell Death"
## [12] "Cytokine Signaling in Immune system"
## [13] "Voltage gated Potassium channels"
## [14] "Response to metal ions"
## [15] "Sensory processing of sound by inner hair cells of the cochlea"
## [16] "Apoptosis"
## [17] "Chromatin modifying enzymes"
## [18] "Chromatin organization"
## [19] "Cardiac conduction"
## [20] "Signaling by Rho GTPases"
## [21] "Intraflagellar transport"
## [22] "Signaling by Rho GTPases, Miro GTPases and RHOBTB3"
## [23] "Cell Cycle, Mitotic"
## [24] "Serotonin Neurotransmitter Release Cycle"
## [25] "Unwinding of DNA"If we consider both strategies to be valid, then we can define the significant sets as dysregulated. We can calculate the percent sentitivity of both approaches.
all <- unique(c(ms_up,ms_dn,mc_up))
message("Sensitivity: separate only")## Sensitivity: separate only(length(ms_up)+length(ms_dn))/length(all)## [1] 0.8974359message("Sensitivity: combined only")## Sensitivity: combined onlylength(mc_up)/length(all)## [1] 0.1737892ORA with clusterprofiler
Clusterprofiler uses a hypergeometric test. Firstly I will conduct the analysis separately for up and down regulated genes and with the correct backgound (as intended by the developers).
genesets2 <- read.gmt("ReactomePathways.gmt")
de_up <- rownames(subset(de, padj<0.05 & log2FoldChange > 0))
de_up <- unique(gt[which(rownames(gt) %in% de_up),1])
de_dn <- rownames(subset(de, padj<0.05 & log2FoldChange < 0))
de_dn <- unique(gt[which(rownames(gt) %in% de_dn),1])
de_bg <- rownames(de)
de_bg <- unique(gt[which(rownames(gt) %in% de_bg),1])
o_up <- as.data.frame(enricher(gene = de_up, universe = de_bg, maxGSSize = 5000, TERM2GENE = genesets2, pAdjustMethod="fdr"))
o_up <- rownames(subset(o_up, p.adjust < 0.05))
o_dn <- as.data.frame(enricher(gene = de_dn, universe = de_bg, maxGSSize = 5000, TERM2GENE = genesets2, pAdjustMethod="fdr"))
o_dn <- rownames(subset(o_dn, p.adjust < 0.05))
o_com <- as.data.frame(enricher(gene = union(de_up,de_dn), universe = de_bg, maxGSSize = 5000, TERM2GENE = genesets2, pAdjustMethod="fdr"))
o_com <- rownames(subset(o_com, p.adjust < 0.05))
length(o_up)## [1] 26length(o_dn)## [1] 120length(o_up) + length(o_dn)## [1] 146length(o_com)## [1] 0( length(o_up) + length(o_dn) ) / length(o_com)## [1] Infall <- unique(c(o_up,o_dn,o_com))
message("Sensitivity: separate only")## Sensitivity: separate only(length(o_up)+length(o_dn))/length(all)## [1] 1message("Sensitivity: combined only")## Sensitivity: combined onlylength(o_com)/length(all)## [1] 0Euler diagram of the significant pathways found with each approach.
l2 <- list("sep up"=o_up,"sep dn"=o_dn,"comb"=o_com)
plot(euler(l2),quantities = TRUE, edges = "gray", main="ORA: combined vs separated")## Warning in colSums(id & !empty) == 0 | merged_sets: longer object length is not
## a multiple of shorter object length
List gene sets which are specific to each approach.
o_sep <- c(o_up,o_dn)
# in sep but not comb
setdiff(o_sep,o_com)## [1] "Neuronal System"
## [2] "Intraflagellar transport"
## [3] "Voltage gated Potassium channels"
## [4] "Signaling by GPCR"
## [5] "Potassium Channels"
## [6] "Dopamine Neurotransmitter Release Cycle"
## [7] "Transmission across Chemical Synapses"
## [8] "Ion channel transport"
## [9] "Metabolism of lipids"
## [10] "Serotonin Neurotransmitter Release Cycle"
## [11] "Neurotransmitter release cycle"
## [12] "GPCR ligand binding"
## [13] "Cilium Assembly"
## [14] "Sensory processing of sound"
## [15] "Ion homeostasis"
## [16] "Protein-protein interactions at synapses"
## [17] "GPCR downstream signalling"
## [18] "Cardiac conduction"
## [19] "Assembly and cell surface presentation of NMDA receptors"
## [20] "Branched-chain amino acid catabolism"
## [21] "Acetylcholine Neurotransmitter Release Cycle"
## [22] "Response to metal ions"
## [23] "Phospholipid metabolism"
## [24] "Muscle contraction"
## [25] "Neurexins and neuroligins"
## [26] "Transport of small molecules"
## [27] "Metabolism of RNA"
## [28] "Nonsense Mediated Decay (NMD) enhanced by the Exon Junction Complex (EJC)"
## [29] "Nonsense-Mediated Decay (NMD)"
## [30] "SRP-dependent cotranslational protein targeting to membrane"
## [31] "rRNA processing"
## [32] "Translation"
## [33] "Signaling by ROBO receptors"
## [34] "Regulation of expression of SLITs and ROBOs"
## [35] "Cap-dependent Translation Initiation"
## [36] "Eukaryotic Translation Initiation"
## [37] "Cell Cycle"
## [38] "Response of EIF2AK4 (GCN2) to amino acid deficiency"
## [39] "Nonsense Mediated Decay (NMD) independent of the Exon Junction Complex (EJC)"
## [40] "Peptide chain elongation"
## [41] "Cell Cycle, Mitotic"
## [42] "L13a-mediated translational silencing of Ceruloplasmin expression"
## [43] "Eukaryotic Translation Termination"
## [44] "Viral mRNA Translation"
## [45] "rRNA processing in the nucleus and cytosol"
## [46] "GTP hydrolysis and joining of the 60S ribosomal subunit"
## [47] "Eukaryotic Translation Elongation"
## [48] "Formation of a pool of free 40S subunits"
## [49] "Chromatin modifying enzymes"
## [50] "Chromatin organization"
## [51] "Major pathway of rRNA processing in the nucleolus and cytosol"
## [52] "Selenocysteine synthesis"
## [53] "Mitotic Metaphase and Anaphase"
## [54] "Selenoamino acid metabolism"
## [55] "Mitotic Anaphase"
## [56] "Cell Cycle Checkpoints"
## [57] "Metabolism of proteins"
## [58] "M Phase"
## [59] "Influenza Infection"
## [60] "Influenza Viral RNA Transcription and Replication"
## [61] "HATs acetylate histones"
## [62] "Cellular responses to stress"
## [63] "Mitotic Prophase"
## [64] "Separation of Sister Chromatids"
## [65] "Cellular response to starvation"
## [66] "Mitotic G1 phase and G1/S transition"
## [67] "Cellular responses to stimuli"
## [68] "DNA Replication"
## [69] "Infectious disease"
## [70] "Developmental Biology"
## [71] "S Phase"
## [72] "Synthesis of DNA"
## [73] "RNA Polymerase II Transcription Termination"
## [74] "Processing of Capped Intron-Containing Pre-mRNA"
## [75] "G1/S Transition"
## [76] "Nervous system development"
## [77] "Transport of Mature Transcript to Cytoplasm"
## [78] "Transport of Mature mRNA derived from an Intron-Containing Transcript"
## [79] "G2/M Checkpoints"
## [80] "Axon guidance"
## [81] "DNA Replication Pre-Initiation"
## [82] "Switching of origins to a post-replicative state"
## [83] "Resolution of Sister Chromatid Cohesion"
## [84] "Disease"
## [85] "Programmed Cell Death"
## [86] "APC/C-mediated degradation of cell cycle proteins"
## [87] "Regulation of mitotic cell cycle"
## [88] "Regulation of APC/C activators between G1/S and early anaphase"
## [89] "Mitotic Spindle Checkpoint"
## [90] "Nuclear Envelope (NE) Reassembly"
## [91] "Activation of APC/C and APC/C:Cdc20 mediated degradation of mitotic proteins"
## [92] "Apoptosis"
## [93] "mRNA 3'-end processing"
## [94] "Transcriptional regulation by RUNX3"
## [95] "Cdc20:Phospho-APC/C mediated degradation of Cyclin A"
## [96] "Orc1 removal from chromatin"
## [97] "APC/C:Cdc20 mediated degradation of mitotic proteins"
## [98] "APC:Cdc20 mediated degradation of cell cycle proteins prior to satisfation of the cell cycle checkpoint"
## [99] "RMTs methylate histone arginines"
## [100] "Assembly of the pre-replicative complex"
## [101] "Initiation of Nuclear Envelope (NE) Reformation"
## [102] "Nuclear Envelope Breakdown"
## [103] "Ribosomal scanning and start codon recognition"
## [104] "Translation initiation complex formation"
## [105] "Unwinding of DNA"
## [106] "Cytokine Signaling in Immune system"
## [107] "Transcriptional regulation by RUNX2"
## [108] "Host Interactions of HIV factors"
## [109] "Signaling by NOTCH4"
## [110] "mRNA Splicing - Major Pathway"
## [111] "Golgi Cisternae Pericentriolar Stack Reorganization"
## [112] "Activation of the mRNA upon binding of the cap-binding complex and eIFs, and subsequent binding to 43S"
## [113] "APC/C:Cdh1 mediated degradation of Cdc20 and other APC/C:Cdh1 targeted proteins in late mitosis/early G1"
## [114] "Polo-like kinase mediated events"
## [115] "EML4 and NUDC in mitotic spindle formation"
## [116] "Interactions of Rev with host cellular proteins"
## [117] "mRNA Splicing"
## [118] "Transcriptional regulation by RUNX1"
## [119] "Amplification of signal from unattached kinetochores via a MAD2 inhibitory signal"
## [120] "Amplification of signal from the kinetochores"
## [121] "Autodegradation of Cdh1 by Cdh1:APC/C"
## [122] "RUNX1 regulates transcription of genes involved in differentiation of HSCs"
## [123] "Class I MHC mediated antigen processing & presentation"
## [124] "Regulation of RUNX3 expression and activity"
## [125] "Chk1/Chk2(Cds1) mediated inactivation of Cyclin B:Cdk1 complex"
## [126] "Transcriptional Regulation by VENTX"
## [127] "rRNA processing in the mitochondrion"
## [128] "tRNA processing in the mitochondrion"
## [129] "The role of GTSE1 in G2/M progression after G2 checkpoint"
## [130] "Activation of the pre-replicative complex"
## [131] "APC/C:Cdc20 mediated degradation of Securin"
## [132] "CDK-mediated phosphorylation and removal of Cdc6"
## [133] "tRNA processing"
## [134] "SCF(Skp2)-mediated degradation of p27/p21"
## [135] "Formation of the ternary complex, and subsequently, the 43S complex"
## [136] "RHO GTPases Activate Formins"
## [137] "SUMO E3 ligases SUMOylate target proteins"
## [138] "Antigen processing-Cross presentation"
## [139] "IRAK4 deficiency (TLR2/4)"
## [140] "Rev-mediated nuclear export of HIV RNA"
## [141] "SUMOylation"
## [142] "Stabilization of p53"
## [143] "Phosphorylation of the APC/C"
## [144] "Regulation of mRNA stability by proteins that bind AU-rich elements"
## [145] "Nucleotide biosynthesis"
## [146] "AUF1 (hnRNP D0) binds and destabilizes mRNA"# in comb but not sep
setdiff(o_com,o_sep)## character(0)# intersection
intersect(mc_up,ms)## [1] "Developmental Biology"
## [2] "Neuronal System"
## [3] "Immune System"
## [4] "Axon guidance"
## [5] "Nervous system development"
## [6] "Potassium Channels"
## [7] "Muscle contraction"
## [8] "Disease"
## [9] "Transmission across Chemical Synapses"
## [10] "Sensory processing of sound"
## [11] "Programmed Cell Death"
## [12] "Cytokine Signaling in Immune system"
## [13] "Voltage gated Potassium channels"
## [14] "Response to metal ions"
## [15] "Sensory processing of sound by inner hair cells of the cochlea"
## [16] "Apoptosis"
## [17] "Chromatin modifying enzymes"
## [18] "Chromatin organization"
## [19] "Cardiac conduction"
## [20] "Signaling by Rho GTPases"
## [21] "Intraflagellar transport"
## [22] "Signaling by Rho GTPases, Miro GTPases and RHOBTB3"
## [23] "Cell Cycle, Mitotic"
## [24] "Serotonin Neurotransmitter Release Cycle"
## [25] "Unwinding of DNA"Euler diagrams comparing FCS and ORA methods
par(cex.main=0.5)
par(mar=c(2,2,2,2))
l3 <- list("ORA up"=o_up,"ORA dn"=o_dn,"ORA comb"=o_com,
"FCS up"=ms_up,"FCS dn"=ms_dn,"FCS comb"=mc_up)
plot(euler(l3),quantities = TRUE, edges = "gray", main="FCS compared to ORA")## Warning in colSums(id & !empty) == 0 | merged_sets: longer object length is not
## a multiple of shorter object length
Save data
dat <- list( "FCS_up"=ms_up,
"FCS_dn"=ms_dn,
"FCS_com"=mc_up,
"ORA_up"= o_up,
"ORA_dn"=o_dn,
"ORA_com"=o_com)
str(dat)## List of 6
## $ FCS_up : chr [1:49] "Neuronal System" "Intraflagellar transport" "Dopamine Neurotransmitter Release Cycle" "Response to metal ions" ...
## $ FCS_dn : chr [1:266] "Metabolism of RNA" "Translation" "Cell Cycle" "Cell Cycle, Mitotic" ...
## $ FCS_com: chr [1:61] "Signal Transduction" "Hemostasis" "Developmental Biology" "Neuronal System" ...
## $ ORA_up : chr [1:26] "Neuronal System" "Intraflagellar transport" "Voltage gated Potassium channels" "Signaling by GPCR" ...
## $ ORA_dn : chr [1:120] "Metabolism of RNA" "Nonsense Mediated Decay (NMD) enhanced by the Exon Junction Complex (EJC)" "Nonsense-Mediated Decay (NMD)" "SRP-dependent cotranslational protein targeting to membrane" ...
## $ ORA_com: chr(0)saveRDS(dat,file = "ex3dat.rds")Conclusion
For mitch, it would appear that performing direction informed (DI) analysis clearly yields more differentially regulated pathways (413) as compared to the direction agnostic (DA) method (55).
That being said, there were 18 pathways identified only in the DA method that appeared to be related to the physiology of the model. These gene sets are likely to contain a mix of genes affected by the stimulus in different ways - for example a mix of up and downregulated genes. Are these really real? Not sure.
This pattern was consistent with ORA, where 80 sets were identified with separate analysis and only 23 with the combined analysis.
When comparing ORA to FCS, we found that FCS identified many more sets than ORA. In fact all gene sets that were identified by ORA were also identified by FCS, except for 3 that were specific to the ORA up set.
Let’s look at those now.
myfcs <- c(ms_up, mc_up)
setdiff(o_up,myfcs)## character(0)Session information
sessionInfo()## R version 4.1.2 (2021-11-01)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.3 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
##
## locale:
## [1] LC_CTYPE=en_AU.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_AU.UTF-8 LC_COLLATE=en_AU.UTF-8
## [5] LC_MONETARY=en_AU.UTF-8 LC_MESSAGES=en_AU.UTF-8
## [7] LC_PAPER=en_AU.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_AU.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] eulerr_6.1.1 mitch_1.4.1
## [3] clusterProfiler_4.0.5 DESeq2_1.32.0
## [5] SummarizedExperiment_1.22.0 Biobase_2.52.0
## [7] MatrixGenerics_1.4.3 matrixStats_0.61.0
## [9] GenomicRanges_1.44.0 GenomeInfoDb_1.28.4
## [11] IRanges_2.26.0 S4Vectors_0.30.2
## [13] BiocGenerics_0.38.0 getDEE2_1.2.0
## [15] beeswarm_0.4.0 kableExtra_1.3.4
##
## loaded via a namespace (and not attached):
## [1] shadowtext_0.1.1 fastmatch_1.1-3 systemfonts_1.0.3
## [4] plyr_1.8.6 igraph_1.2.11 lazyeval_0.2.2
## [7] polylabelr_0.2.0 splines_4.1.2 BiocParallel_1.26.2
## [10] ggplot2_3.3.5 digest_0.6.29 yulab.utils_0.0.4
## [13] htmltools_0.5.2 GOSemSim_2.18.1 viridis_0.6.2
## [16] GO.db_3.13.0 fansi_1.0.0 magrittr_2.0.1
## [19] memoise_2.0.1 Biostrings_2.60.2 annotate_1.70.0
## [22] graphlayouts_0.8.0 svglite_2.0.0 enrichplot_1.12.3
## [25] colorspace_2.0-2 blob_1.2.2 rvest_1.0.2
## [28] ggrepel_0.9.1 xfun_0.29 dplyr_1.0.7
## [31] crayon_1.4.2 RCurl_1.98-1.5 jsonlite_1.7.2
## [34] scatterpie_0.1.7 genefilter_1.74.1 survival_3.2-13
## [37] ape_5.6-1 glue_1.6.0 polyclip_1.10-0
## [40] gtable_0.3.0 zlibbioc_1.38.0 XVector_0.32.0
## [43] webshot_0.5.2 htm2txt_2.1.1 DelayedArray_0.18.0
## [46] scales_1.1.1 DOSE_3.18.3 DBI_1.1.2
## [49] GGally_2.1.2 Rcpp_1.0.7 viridisLite_0.4.0
## [52] xtable_1.8-4 gridGraphics_0.5-1 tidytree_0.3.7
## [55] bit_4.0.4 htmlwidgets_1.5.4 httr_1.4.2
## [58] fgsea_1.18.0 gplots_3.1.1 RColorBrewer_1.1-2
## [61] ellipsis_0.3.2 reshape_0.8.8 pkgconfig_2.0.3
## [64] XML_3.99-0.8 farver_2.1.0 sass_0.4.0
## [67] locfit_1.5-9.4 utf8_1.2.2 later_1.3.0
## [70] ggplotify_0.1.0 tidyselect_1.1.1 rlang_0.4.12
## [73] reshape2_1.4.4 AnnotationDbi_1.54.1 munsell_0.5.0
## [76] tools_4.1.2 cachem_1.0.6 downloader_0.4
## [79] generics_0.1.1 RSQLite_2.2.9 evaluate_0.14
## [82] stringr_1.4.0 fastmap_1.1.0 yaml_2.2.1
## [85] ggtree_3.0.4 knitr_1.37 bit64_4.0.5
## [88] tidygraph_1.2.0 caTools_1.18.2 purrr_0.3.4
## [91] KEGGREST_1.32.0 ggraph_2.0.5 nlme_3.1-153
## [94] mime_0.12 aplot_0.1.2 DO.db_2.9
## [97] xml2_1.3.3 compiler_4.1.2 rstudioapi_0.13
## [100] png_0.1-7 treeio_1.16.2 tibble_3.1.6
## [103] tweenr_1.0.2 geneplotter_1.70.0 bslib_0.3.1
## [106] stringi_1.7.6 highr_0.9 lattice_0.20-45
## [109] Matrix_1.4-0 vctrs_0.3.8 pillar_1.6.4
## [112] lifecycle_1.0.1 jquerylib_0.1.4 data.table_1.14.2
## [115] cowplot_1.1.1 bitops_1.0-7 httpuv_1.6.5
## [118] patchwork_1.1.1 qvalue_2.24.0 R6_2.5.1
## [121] promises_1.2.0.1 KernSmooth_2.23-20 echarts4r_0.4.3
## [124] gridExtra_2.3 gtools_3.9.2 MASS_7.3-54
## [127] assertthat_0.2.1 GenomeInfoDbData_1.2.6 grid_4.1.2
## [130] ggfun_0.0.4 tidyr_1.1.4 rmarkdown_2.11
## [133] ggforce_0.3.3 shiny_1.7.1