tcga中的gsva分析

对样本进行分级后，通路的分析不适合传统的kegg和gsea分析，下面把单细胞的gsva分析带入tcga分析，原理一摸一样。

library(ggplot2)
library(clusterProfiler)
library(org.Hs.eg.db)
library(dplyr)
expr <- tcga_panmRNA_expr[,names]
expr$ID <- rownames(expr)
s2e <- bitr(expr$ID, 
            fromType = "SYMBOL",
            toType = "ENTREZID",
            OrgDb = org.Hs.eg.db)#人类 转换成ENTREZID
#其他物种http://bioconductor.org/packages/release/BiocViews.html#___OrgDb
expr <- inner_join(expr,s2e,by=c("ID"="SYMBOL"))
rownames(expr) <- expr$ENTREZID
expr <- expr[,-453]
expr <- expr[,-452]
h.all.v7.2.entrez.gmt
kegggmt1 <- read.gmt("c2.cp.kegg.v7.2.entrez.gmt")
kegggmt2 <- read.gmt("c2.cp.reactome.v7.2.entrez.gmt")
kegggmt3 <- read.gmt("h.all.v7.2.entrez.gmt")
kegggmt <- rbind(kegggmt1,kegggmt2,kegggmt3)
colnames(kegggmt3)
kegg_list = split(kegggmt$gene, kegggmt$ont)
library(GSVA)
expr=as.matrix(expr)
kegg2 <- gsva(expr, kegg_list, kcdf="Gaussian",method = "gsva",parallel.sz=4)
write.table(kegg2,"kegg2.txt",sep="\t",quote = F)
data <- data.table::fread('kegg2.txt',data.table = F)
rownames(data) <- data$V1
data <- data[,-1]
identical(colnames(data),rownames(dd1))
table(dd1$CIC)
tonglu <- c('KEGG_MISMATCH_REPAIR','KEGG_BASE_EXCISION_REPAIR','KEGG_NUCLEOTIDE_EXCISION_REPAIR','KEGG_CELL_CYCLE',
            'HALLMARK_MITOTIC_SPINDLE','HALLMARK_G2M_CHECKPOINT','HALLMARK_E2F_TARGETS','HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION','HALLMARK_ANGIOGENESIS','KEGG_VASCULAR_SMOOTH_MUSCLE_CONTRACTION','KEGG_ECM_RECEPTOR_INTERACTION','HALLMARK_WNT_BETA_CATENIN_SIGNALING',
            'HALLMARK_TGF_BETA_SIGNALING','HALLMARK_NOTCH_SIGNALING','KEGG_MTOR_SIGNALING_PATHWAY','KEGG_FOCAL_ADHESION','KEGG_ADHERENS_JUNCTION','KEGG_TIGHT_JUNCTION','KEGG_GAP_JUNCTION','HALLMARK_CHOLESTEROL_HOMEOSTASIS','HALLMARK_FATTY_ACID_METABOLISM','HALLMARK_GLYCOLYSIS'
            ,'HALLMARK_ADIPOGENESIS','KEGG_GLYCOLYSIS_GLUCONEOGENESIS','KEGG_FRUCTOSE_AND_MANNOSE_METABOLISM','KEGG_CYTOKINE_CYTOKINE_RECEPTOR_INTERACTION','KEGG_CHEMOKINE_SIGNALING_PATHWAY','HALLMARK_COMPLEMENT','HALLMARK_INTERFERON_ALPHA_RESPONSE','HALLMARK_INTERFERON_GAMMA_RESPONSE',
            'HALLMARK_INFLAMMATORY_RESPONSE','HALLMARK_TNFA_SIGNALING_VIA_NFKB','HALLMARK_IL6_JAK_STAT3_SIGNALING')

data1 <- data[tonglu,1:451]
pheatmap(data1,
         cluster_rows = F,
         #annotation_col=dd1$CIC, #标注样本分类
         cluster_cols = F,
         show_rownames = T,
         show_colnames = T,
         scale = "row",
         color =colorRampPalette(c("blue", "white","red"))(100),
         cellwidth = 10, cellheight = 10,
         fontsize = 10)

library(dplyr)
data$A <- apply(data[,1:127], 1, mean)
data$B <- apply(data[,128:234], 1, mean)
data$C <- apply(data[,235:349], 1, mean)
data$D <- apply(data[,350:451], 1, mean)

table(meta$celltype3)
test <- data[tonglu,452:455]
library(pheatmap)
W3 <- pheatmap(test,
               cluster_rows = F,
               cluster_cols = F,
               show_rownames = T,
               show_colnames = T,
             scale = "row",
               color =colorRampPalette(c("blue", "white","red"))(100),
               cellwidth =30, cellheight = 15,
               fontsize = 10)
pdf(("W3.pdf"),width = 8,height =8)
print(W3)
dev.off()

分别取kegg和代谢的子集

test1 <- test[grep(pattern="METABOLISM",rownames(test)),]
G3 <- pheatmap(test1,
               cluster_rows = F,
               cluster_cols = F,
               show_rownames = T,
               show_colnames = T,
               scale = 'row',
               color =colorRampPalette(c("blue", "white","red"))(100),
               cellwidth = 10, cellheight = 10,
               fontsize = 10)
pdf(("G3.pdf"),width = 20,height =30)
print(G3)
dev.off()


test3 <- test[grep(pattern="KEGG",rownames(test)),]
G4 <- pheatmap(test3,
               cluster_rows = F,
               cluster_cols = F,
               show_rownames = T,
               show_colnames = T,
               scale = 'row',
               color =colorRampPalette(c("blue", "white","red"))(100),
               cellwidth = 10, cellheight = 10,
               fontsize = 10)
pdf(("G4.pdf"),width = 20,height =30)
print(G4)
dev.off()