metaunicox.R

m<-read.table("metaunicox")
mbb<-cbind(rownames(m),m$risk)
mbt<-t(mbb)
lncrg<-read.table("chose_lncp_ex_rg")
lnrgt<-t(lncrg)
colnames(lnrgt)<-rownames(lncrg)
rownames(lnrgt)<-colnames(lncrg)
lnrgt<-as.data.frame(lnrgt)
xl=subset(lnrgt,risk_rg=="0")
xh=subset(lnrgt,risk_rg=="1")
xlmean=apply(xl,2,mean)
xhmean=apply(xh,2,mean)
xx<-rbind(xlmean,xhmean)   
xxt<-as.data.frame(t(xx))
xxt$fc=log(xxt$xhmean/xxt$xlmean,2)
write.table(xxt,"lnc_ex_rg_fc",sep="\t") 

xxt<-read.table("lnc_ex_rg_fc",sep="\t")
library(enrichplot)
library(ChIPseeker)
library(GSEABase)
library(clusterProfiler)
library(org.Hs.eg.db)

xxt$entrez=?{probe...annot}

gmt<-"https://wikipathways-data.wmcloud.org/current/gmt/wikipathways-20211110-gmt-Homo_sapiens.gmt"
wp <- read.gmt.wp(gmt)

id= annot[match(rownames(xxt),annot$probe),2]
fc= xxt[match(rownames(xxt),annot$probe),3]
id=id[!is.na(id)]
fc=fc[!is.na(fc)]

a<-data.frame(id=id,fc=fc)
a=t(a[!is.na(a$id),])
a<-as.data.frame(a)
genelist<-as.numeric(a[2,])
names(genelist) = as.character(a[1,])
genelist = sort(genelist, decreasing = TRUE)
ewp <- GSEA(genelist,TERM2GENE = wp[,c("wpid","gene")],TERM2NAME = wp[,c("wpid","name")])



probe<-rownames(xxt)
ensg2entrez=toTable(org.Hs.egENSEMBL)
id= ensg2entrez[match(probe,ensg2entrez$ensembl_id),1]
annot=data.frame(probe=probe,"EntrezGene.ID"=id)   
annot=annot[!is.na(annot$EntrezGene.ID),]