In this document, we will show how to predict master regulators that potential drive AML fitness using TCGA AML gene expression data.
library("ggplot2")
library("viper")
source("/sibcb2/bioinformatics/Script/PathwayEnrichment.R")
load("/sibcb2/bioinformatics/KnowledgeBase/Firehose_Expression/Firehose_Expression_Tumor.RData")
load("/sibcb2/bioinformatics/Placenta/RData/TCGA_regulonDB.RData")
load("/sibcb2/bioinformatics/LSC/RData/Signature_Top500.RData")
pw = list(up = Signature_Top500$sigUP_Hs, dn = Signature_Top500$sigDN_Hs)
regulon = TCGA_regulonDB[["regulonlaml"]]
Ex = Firehose_Expression_Tumor[["LAML"]]
ref = apply(Ex, 1, mean)
dM = Ex - ref
zscore = rep(NA, ncol(dM))
for(j in 1:ncol(dM)){
temp = SimpleRankTest(dM[,j], pw)
zscore[j] = temp[1] - temp[2]
}
names(zscore) = colnames(dM)
Enrichment = rep("others", length(zscore))
Enrichment[zscore >= quantile(zscore, 1 - 50/length(zscore))] = "high"
Enrichment[zscore <= quantile(zscore, 50/length(zscore))] = "low"
temp = data.frame(rnk = rank(zscore), zscore, Enrichment)
p <- ggplot(temp, aes(x = rnk, y = zscore, color = Enrichment, fill = Enrichment)) + geom_bar(stat="identity") + theme_bw()
p <- p + theme(plot.title = element_text(hjust = 0.5))
p <- p + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
p <- p + ylab("Enrichment z-score") + xlab("") + theme(axis.text.x=element_blank(), axis.ticks.x=element_blank())
p <- p + scale_colour_manual(values = c("red", "blue", "grey"))
print(p)
ExH = Ex[, names(sort(zscore, decreasing = TRUE)[1:50])]
ExL = Ex[, names(sort(zscore, decreasing = FALSE )[1:50])]
signature = rowTtest(ExH, ExL)
signature = (qnorm(signature$p.value / 2, lower.tail = FALSE) * sign(signature$statistic))[, 1]
nullmodel = ttestNull(ExH, ExL, per = 500, repos = TRUE, verbose = FALSE)
mrs = msviper(signature, regulon, nullmodel, verbose = FALSE)
knitr::kable(data.frame(summary(mrs, 10)))
Regulon | Size | NES | p.value | FDR | |
---|---|---|---|---|---|
YWHAE | YWHAE | 95 | 4.91 | 9.00e-07 | 0.00273 |
SDHD | SDHD | 89 | 4.65 | 3.30e-06 | 0.00411 |
ZBTB8OS | ZBTB8OS | 94 | 4.60 | 4.20e-06 | 0.00411 |
PRDM15 | PRDM15 | 72 | -4.29 | 1.76e-05 | 0.01040 |
LCK | LCK | 33 | -4.45 | 8.80e-06 | 0.00574 |
SEPT5 | SEPT5 | 108 | -4.46 | 8.10e-06 | 0.00574 |
EPHA4 | EPHA4 | 45 | -4.48 | 7.40e-06 | 0.00574 |
KLF8 | KLF8 | 137 | -4.60 | 4.10e-06 | 0.00411 |
RUNDC3A | RUNDC3A | 70 | -4.79 | 1.60e-06 | 0.00325 |
TREML1 | TREML1 | 67 | -4.92 | 9.00e-07 | 0.00273 |
sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-conda_cos6-linux-gnu (64-bit)
Running under: Red Hat Enterprise Linux Server 7.7 (Maipo)
Matrix products: default
BLAS: /sibcb2/bioinformatics/software/BcbioNG/anaconda/lib/R/lib/libRblas.so
LAPACK: /sibcb2/bioinformatics/software/BcbioNG/anaconda/lib/R/lib/libRlapack.so
locale:
[1] C
attached base packages:
[1] parallel grid stats graphics grDevices utils datasets
[8] methods base
other attached packages:
[1] viper_1.16.0 scales_1.1.0 pheatmap_1.0.12
[4] dplyr_1.0.5 survcomp_1.32.0 prodlim_2019.11.13
[7] RColorBrewer_1.1-2 gridExtra_2.3 survHD_0.99.1
[10] survC1_1.0-2 Hmisc_4.3-0 Formula_1.2-3
[13] lattice_0.20-38 penalized_0.9-51 survival_3.1-8
[16] Biobase_2.42.0 BiocGenerics_0.28.0 VennDiagram_1.6.20
[19] futile.logger_1.4.3 ggrepel_0.8.1 ggplot2_3.2.1
[22] rmarkdown_2.0
loaded via a namespace (and not attached):
[1] mixtools_1.1.0 splines_3.5.1 rmeta_3.0
[4] highr_0.8 latticeExtra_0.6-28 bootstrap_2019.6
[7] yaml_2.2.0 survivalROC_1.0.3 pillar_1.5.1
[10] backports_1.1.5 glue_1.4.2 digest_0.6.23
[13] checkmate_1.9.4 colorspace_1.4-1 htmltools_0.5.1.1
[16] Matrix_1.2-18 pkgconfig_2.0.3 purrr_0.3.3
[19] lava_1.6.7 htmlTable_1.13.3 tibble_3.1.0
[22] generics_0.0.2 farver_2.0.1 ellipsis_0.3.0
[25] withr_2.1.2 nnet_7.3-12 lazyeval_0.2.2
[28] magrittr_1.5 crayon_1.3.4 evaluate_0.14
[31] fansi_0.4.0 segmented_1.1-0 MASS_7.3-51.5
[34] class_7.3-15 foreign_0.8-74 SuppDists_1.1-9.5
[37] tools_3.5.1 data.table_1.12.8 formatR_1.7
[40] lifecycle_1.0.0 stringr_1.4.0 munsell_0.5.0
[43] cluster_2.1.0 lambda.r_1.2.4 e1071_1.7-3
[46] compiler_3.5.1 rlang_0.4.10 rstudioapi_0.10
[49] htmlwidgets_1.5.1 base64enc_0.1-3 labeling_0.3
[52] gtable_0.3.0 DBI_1.1.0 R6_2.4.1
[55] knitr_1.26 utf8_1.1.4 futile.options_1.0.1
[58] KernSmooth_2.23-16 stringi_1.4.3 Rcpp_1.0.3
[61] vctrs_0.3.7 rpart_4.1-15 acepack_1.4.1
[64] tidyselect_1.1.0 xfun_0.18