install.packages("BiocManager")

BiocManager::install("GEOquery")
BiocManager::install("DO.db")
BiocManager::install("clusterProfiler")

BiocManager::install("KEGG.db")
BiocManager::install("reactome.db")
BiocManager::install("KEGGREST.db")

BiocManager::install("topGO") #画GO图用的
BiocManager::install("Rgraphviz")
BiocManager::install("pathview") #看KEGG pathway的
BiocManager::install("org.Hs.eg.db") #这个包里存有人的注释文件

library(GEOquery)
library(clusterProfiler)
library(org.Hs.eg.db)

library(topGO)
library(Rgraphviz)
library(pathview)
library(KEGG.db)
library(reactome.db)
library(KEGGREST.db)

setwd("/Users/mumu/Documents/Documents/R_DataAnalyzed/TREK/ProcessedData")
getwd()

load("HT294D_Regualtion_TFgenes_Activation&Inhibition_isodice.rda")
ls()

R.mat<-R.mat[,which(colSums(abs(R.mat))>5)]
enrich<-Enrich(R.mat)
enrich<-MergeMatrix(enrich)

Enrich<-function(R.mat){

enrich<-list()

for (i in 1:10){

print (i)

tf<-colnames(R.mat)[i]
deg<-R.mat[,i]
deg<-names(deg)[which(abs(deg)!=0)]

DEG.gene_symbol<-as.character(deg)

DEG.entrez_id = mapIds(x = org.Hs.eg.db,
keys = DEG.gene_symbol,
keytype = "SYMBOL",
column = "ENTREZID")

head(DEG.entrez_id)

DEG.entrez_id = na.omit(DEG.entrez_id)

KEGG analysis

kegg <- enrichKEGG(gene = DEG.entrez_id, organism = 'hsa', keyType = 'kegg', pvalueCutoff = 0.05,pAdjustMethod = 'BH',

minGSSize = 10,maxGSSize = 500,qvalueCutoff = 0.2,use_internal_data = FALSE)

ekegg = enrichKEGG(gene = DEG.entrez_id,
organism = 'hsa',
keyType = 'kegg',
pAdjustMethod = 'BH',
#pvalueCutoff = 0.05,
# qvalueCutoff = 0.2,
use_internal_data=FALSE
)

mi<-as.matrix(cbind(ekeggp.adjust))
colnames(mi)<-c("Description",tf)

rownames(mi)<-ekegg$Description

mi<-list(tf=mi)
names(mi)<-tf
enrich<-append(enrich, mi)

}

return (enrich)

}

merge(enrich[[1]], enrich[[2]], by="Description",all.x=TRUE, all.y=TRUE)

MergeMatrix<-function(enrich){

mat0<-enrich[[1]]

for (i in 1:length(enrich)){

mat0<-merge(mat0, enrich[[i]], by="Description",all.x=TRUE, all.y=TRUE)

}

mat0["NA"]<-1

return (mat0)

}