【后续来了】有了这个包,猪的GSEA和GSVA分析也不在话下(第二集)

书接上文,gsva后续的处理很简单,如果熟悉差异基因分析包limma的话,更是简单。之前我们得到的gsva分数的矩阵就类似于基因表达矩阵,按照这个思路继续往下即可:
从通路的表达矩阵开始,我们进行差异分析,首先将之前保存的文件读入:

T_gsva <- read.csv("T_gsva.csv", header = T, row.names = 1)

分析得到的数据结构是行为GO terms, 列为样品(单细胞中为单个细胞)。

group <- c(rep("control", 50), rep("test", 71)) %>% as.factor()#设置分组,对照在前desigN <- model.matrix(~ 0 + group) #构建比较矩阵colnames(desigN) <- levels(group)fit = lmFit(test_control, desigN)fit2 <- eBayes(fit)diff=topTable(fit2,adjust='fdr',coef=2,number=Inf)#校准按照需求修改 ?topTablewrite.csv(diff, file = "Diff.csv")

差异分析结果有logFC,P-value,t值等等,根据阈值筛选差异的通路即可,最后挑选与自己研究相关的或者感兴趣的通路进行可视化。可视化的类型有热图展示,个人认为不好看,因为细胞数太多,也有用平均值做的,效果一般!
这里看到了一篇文章中的可视化,感觉不错,大多数文章也是这种类型的柱状图!

文章引用:[1] Lambrechts, D. , et al. "Phenotype molding of stromal cells in the lung tumor microenvironment." Nature Medicine 24.8(2018).

up <- c("GOBP_EGG_ACTIVATION", "GOBP_TENDON_DEVELOPMENT", "GOBP_SOMITE_SPECIFICATION", "GOBP_THREONINE_CATABOLIC_PROCESS", "GOBP_REGULATION_OF_GLUTAMATE_RECEPTOR_CLUSTERING", "GOBP_NEGATIVE_CHEMOTAXIS", "GOBP_NEGATIVE_REGULATION_OF_FAT_CELL_PROLIFERATION", "GOBP_REGULATION_OF_T_HELPER_17_CELL_LINEAGE_COMMITMENT", "GOBP_REGULATION_OF_ANTIMICROBIAL_HUMORAL_RESPONSE")down <- c("GOBP_DETERMINATION_OF_PANCREATIC_LEFT_RIGHT_ASYMMETRY", "GOBP_MITOTIC_DNA_REPLICATION", "GOBP_EOSINOPHIL_CHEMOTAXIS", "GOBP_NEUTROPHIL_MEDIATED_CYTOTOXICITY", "GOBP_POTASSIUM_ION_EXPORT_ACROSS_PLASMA_MEMBRANE", "GOBP_REGULATION_OF_LEUKOCYTE_MEDIATED_CYTOTOXICITY", "GOBP_REGULATION_OF_SEQUESTERING_OF_ZINC_ION", "GOBP_ENDOTHELIN_RECEPTOR_SIGNALING_PATHWAY", "GOBP_PRE_REPLICATIVE_COMPLEX_ASSEMBLY_INVOLVED_IN_CELL_CYCLE_DNA_REPLICATION", "GOBP_ESTABLISHMENT_OF_PLANAR_POLARITY_OF_EMBRYONIC_EPITHELIUM")TEST <- c(up,down)diff$ID <- rownames(diff) Q <- diff[TEST,]group1 <- c(rep("treat", 9), rep("control", 10)) df <- data.frame(ID = Q$ID, score = Q$t,group=group1 )# 按照t score排序sortdf <- df[order(df$score),]sortdf$ID <- factor(sortdf$ID, levels = sortdf$ID)#增加通路ID那一列

用ggplot画图(ggplot-YYDS)

ggplot(sortdf, aes(ID, score,fill=group)) + geom_bar(stat = 'identity',alpha = 0.7) + coord_flip() + theme_bw() + #去除背景色 theme(panel.grid =element_blank())+ theme(panel.border = element_rect(size = 0.6))+ labs(x = "", y="t value of GSVA score")+ scale_fill_manual(values = c("#008020","#08519C"))#设置颜色

有闲工夫了可以自己修饰修饰即可!
GSVA分析到此结束!
欢迎分享交流!

(0)

相关推荐