Signatures
source("code/MR.R")
load("Data/LungMR_GRN.RData")
load("Data/GSE123066_DGE.RData")
load("Data/GSE121634_DGE.RData")
set1 = rownames(GSE123066_DGE$uTable)
set2 = rownames(GSE121634_DGE$uTable)
groupColor = c("#08519c", "#ce1256")
venn.plot <- draw.pairwise.venn(
area1 = length(set1),
area2 = length(set2),
cross.area = length(intersect(set1, set2)),
category = c("GSE123066", "GSE121634"),
fill = c("white", "white"),
lty = "solid",
lwd = 3,
cex = 2,
cat.cex = 1,
col = groupColor,
cat.col = groupColor,
euler = TRUE)
load("/sibcb2/bioinformatics/KnowledgeBase/MSigDB/MSigDB.RData")
signature = intersect(set1, set2)
EMT = MSigDB$MSigDB_h_all_v5$HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION
venn.plot <- draw.pairwise.venn(
area1 = length(signature),
area2 = length(EMT),
cross.area = length(intersect(signature, EMT)),
category = c("signature", "EMT"),
fill = c("white", "white"),
lty = "solid",
lwd = 3,
cex = 2,
cat.cex = 1,
col = groupColor,
cat.col = groupColor,
euler = TRUE)
Master regulators
T1 = mrMulti(GEX, GRN, gs = list(EMT = EMT), nsample = 50)
write.table(T1, file = "out/LungGRN_EMT.txt", row.names = FALSE, col.names = TRUE, sep = "\t")
T2 = mrMulti(GEX, GRN, gs = list(TKI = signature), nsample = 50)
write.table(T2, file = "out/LungGRN_TKI.txt", row.names = FALSE, col.names = TRUE, sep = "\t")
T1 = read.table(file = "out/LungGRN_EMT.txt", row.names = NULL, header = TRUE, sep = "\t")
T1$signature = "EMT"
T2 = read.table(file = "out/LungGRN_TKI.txt", row.names = NULL, header = TRUE, sep = "\t")
T2$signature = "signature"
mergedT = rbind(T1, T2)
uTable = filter(mergedT, p.value < 0.01, NES > 0)
dTable = filter(mergedT, p.value < 0.01, NES < 0)
x = table(uTable$Regulon, uTable$signature)
y = as.data.frame.matrix(x) %>% filter(EMT + signature > 0)
z = melt(updateTable(table(y[,1], y[,2]), space = 0:7, N00 = TRUE)) %>% filter(Var1 + Var2 > 0)
p = ggplot(z, aes(Var1, Var2)) + theme_bw() +
labs(x = "# of data sets\n(EMT pathway)", y = "# of data sets\n(EGFR TKI signatue)", title = "") +
geom_tile(aes(fill = value)) +
geom_text(aes(label = round(value, 1))) +
scale_fill_gradient(low = "white", high = "red")
print(p)
filter(y, EMT >=6, signature >=6)
EMT signature
ADGRA2 7 6
AXL 7 7
CAVIN1 7 7
GAS1 7 6
GAS7 6 7
LMCD1 6 6
LRP1 6 7
PDGFRB 7 6
TGFB1I1 7 6
ZCCHC24 6 6
ZFPM2 7 6
ZNF423 6 7
write.table(y, file = "out/TKI_MR_Activated.txt", row.names = TRUE, col.names = TRUE, sep = "\t")
Pathway enrichment
source("/sibcb2/bioinformatics/Script/PathwayEnrichment.R")
PG = as.character(read.table(file = "Data/VIPER_Regulators.txt", row.names = NULL, col.names=TRUE,sep = "\t")[,1])
pwList = MSigDB$MSigDB_h_all_v5
res = as.data.frame.matrix(x)
TopSignature = rownames(filter(res, signature >= 3))
TopEMT = rownames(filter(res, EMT >= 3))
EnrichSignature = pHyperList(TopSignature, pwList, PG) %>% arrange(pvalue)
EnrichEMT = pHyperList(TopEMT, pwList, PG) %>% arrange(pvalue)
filter(EnrichSignature, qvalue < 0.05)[,1:5]
sg sr n pvalue qvalue
HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION 221 64 33 6.412226e-33 3.206113e-31
HALLMARK_COMPLEMENT 221 93 16 2.883499e-08 7.208746e-07
HALLMARK_MYOGENESIS 221 84 14 2.678771e-07 4.464618e-06
HALLMARK_TNFA_SIGNALING_VIA_NFKB 221 128 15 1.422406e-05 1.778008e-04
HALLMARK_HYPOXIA 221 79 11 2.443860e-05 2.443860e-04
HALLMARK_KRAS_SIGNALING_UP 221 97 12 4.414370e-05 3.678642e-04
HALLMARK_UV_RESPONSE_DN 221 84 11 4.588782e-05 3.277702e-04
HALLMARK_INFLAMMATORY_RESPONSE 221 141 15 4.811660e-05 3.007288e-04
HALLMARK_APICAL_JUNCTION 221 83 10 1.903429e-04 1.057461e-03
HALLMARK_COAGULATION 221 44 6 9.504916e-04 4.752458e-03
HALLMARK_ANGIOGENESIS 221 13 3 9.527980e-04 4.330900e-03
HALLMARK_TGF_BETA_SIGNALING 221 47 6 1.421659e-03 5.923578e-03
HALLMARK_PANCREAS_BETA_CELLS 221 21 3 6.339816e-03 2.438391e-02
HALLMARK_WNT_BETA_CATENIN_SIGNALING 221 36 4 9.273279e-03 3.311885e-02
HALLMARK_APOPTOSIS 221 103 8 1.251605e-02 4.172018e-02
HALLMARK_NOTCH_SIGNALING 221 27 3 1.567437e-02 4.749808e-02
HALLMARK_HEDGEHOG_SIGNALING 221 27 3 1.567437e-02 4.749808e-02
write.table(EnrichSignature, file = "out/TKI_signature_pathway_Activated.txt",
row.names = TRUE, col.names = TRUE, sep = "\t")
filter(EnrichEMT, qvalue < 0.05)[,1:5]
sg sr n pvalue qvalue
HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION 235 64 37 2.961978e-38 1.480989e-36
HALLMARK_TNFA_SIGNALING_VIA_NFKB 235 128 23 7.209090e-11 1.802273e-09
HALLMARK_HYPOXIA 235 79 18 8.362838e-11 1.393806e-09
HALLMARK_KRAS_SIGNALING_UP 235 97 19 5.132346e-10 6.415432e-09
HALLMARK_INFLAMMATORY_RESPONSE 235 141 21 1.703149e-08 1.703149e-07
HALLMARK_MYOGENESIS 235 84 14 5.941812e-07 4.951510e-06
HALLMARK_ANGIOGENESIS 235 13 5 4.382058e-06 3.130042e-05
HALLMARK_COMPLEMENT 235 93 13 1.152288e-05 7.201802e-05
HALLMARK_UV_RESPONSE_DN 235 84 12 1.753377e-05 9.740984e-05
HALLMARK_ALLOGRAFT_REJECTION 235 121 13 2.231515e-04 1.115758e-03
HALLMARK_COAGULATION 235 44 7 2.414010e-04 1.097277e-03
HALLMARK_APOPTOSIS 235 103 11 5.879662e-04 2.449859e-03
HALLMARK_INTERFERON_GAMMA_RESPONSE 235 93 10 8.724607e-04 3.355618e-03
HALLMARK_TGF_BETA_SIGNALING 235 47 6 2.025044e-03 7.232301e-03
HALLMARK_IL2_STAT5_SIGNALING 235 107 10 2.752383e-03 9.174611e-03
HALLMARK_P53_PATHWAY 235 108 10 2.962071e-03 9.256470e-03
HALLMARK_APICAL_JUNCTION 235 83 8 4.611105e-03 1.356207e-02
HALLMARK_APICAL_SURFACE 235 21 3 7.865004e-03 2.125677e-02
HALLMARK_PANCREAS_BETA_CELLS 235 21 3 7.865004e-03 2.125677e-02
HALLMARK_ANDROGEN_RESPONSE 235 46 5 8.180400e-03 2.045100e-02
HALLMARK_CHOLESTEROL_HOMEOSTASIS 235 22 3 9.325421e-03 2.220338e-02
HALLMARK_IL6_JAK_STAT3_SIGNALING 235 66 6 1.346515e-02 3.060261e-02
HALLMARK_HEDGEHOG_SIGNALING 235 27 3 1.923565e-02 4.181663e-02
HALLMARK_MTORC1_SIGNALING 235 57 5 2.229334e-02 4.644446e-02
HALLMARK_ESTROGEN_RESPONSE_EARLY 235 91 7 2.434416e-02 4.868831e-02
write.table(EnrichEMT, file = "out/TKI_EMT_pathway_Activated.txt",
row.names = TRUE, col.names = TRUE, sep = "\t")
shownPathway = c("HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION",
"HALLMARK_COMPLEMENT",
"HALLMARK_TNFA_SIGNALING_VIA_NFKB",
"HALLMARK_HYPOXIA",
"HALLMARK_INFLAMMATORY_RESPONSE",
"HALLMARK_ANGIOGENESIS",
"HALLMARK_TGF_BETA_SIGNALING",
"HALLMARK_WNT_BETA_CATENIN_SIGNALING",
"HALLMARK_APOPTOSIS",
"HALLMARK_NOTCH_SIGNALING")
EnrichMAP = data.frame(pathway = gsub("HALLMARK_", "", c(shownPathway, shownPathway)),
qvalue = c(EnrichSignature[shownPathway,"qvalue"], EnrichEMT[shownPathway, "qvalue"]),
group = c(rep("signature", length(shownPathway)), rep("EMT", length(shownPathway))))
EnrichMAP$pathway = factor(EnrichMAP$pathway, levels = gsub("HALLMARK_", "", shownPathway))
p <- ggplot(EnrichMAP, aes(x = pathway, y = -log10(qvalue), fill = group)) +
geom_bar(stat="identity", position = "dodge") + theme_bw() +
theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank()) +
ylab("-log10(FDR)") + xlab("") + theme(axis.text.x=element_text(angle = 60, hjust = 1, size = 8)) +
geom_hline(yintercept = -log10(0.05), linetype="longdash", color = "black", size=0.5) +
scale_fill_manual(values = c("signature" = "red", "EMT" = "dodgerblue4"))
print(p)
