首先,向给我留言及私信的小伙伴们说一声抱歉,我这段时间太忙了,很久很久没登录简书,所以留言和私信都没回复,抱歉啦;我看着留言和私信的时间比较久远了,也不知道小伙伴们的问题有没有通过其他途径解决,这里我就不再废话了(如果仍然需要我解答的,并且我会的问题,麻烦再次私信我,没有的话就继续疾病本体分析啦。有建议也可以给我提,比如有小伙伴说没有R基础,让我把代码解释的再详细些,那我后面就将解释写的再详细一点;如果需要简洁一点,我就写简洁一点呀;在我力所能及的范围内,尽量满足一路生花的小伙伴们;如果你们的问题我解决不了,我也会去请教可以解决这些问题的朋友,不要怕,我们一起进步;最后祝大家都快乐学习,越来越顺利。)
一、举例回顾
本节所使用GSE1009数据集,已经用limma包进行差异分析,对DEGs做GO、KEGG富集分析;现对DEGs做DO(疾病本体)分析
GSE1009数据集介绍:
样本量:共6个样本,其中后3个为糖尿病肾病(DN)肾小球样本,前3个为正常肾小球样本。
使用芯片:Affymetrix Human Genome U95 Version 2 Array。
平台:GPL8300。
DEGs:共有66个DEGs(diffsig),22个上调(diffup),44个上调(diffDown)(详见13章).
二、需要准备的文件:
包含差异基因名字+logFC值的文本文件,命名为symbol(上一章有介绍详细做法。)
三、具体步骤:
1. ID转换(同14章,我就直接将代码粘贴过来了)
setwd("D:\\Rfile")
rm(list = ls())
options(stringsAsFactors=F)
#老规矩,先设置工作目录。
library("clusterProfiler")
library("org.Hs.eg.db")
library("enrichplot")
library("ggplot2")
##在ID转换前,将要分析的基因名字+logFC保存为symbol.text文件。
##ID转换
rt=read.table("symbol.txt",sep="\t",check.names=F,header=T) #读取文件
genes=as.vector(rt[,1])
entrezIDs <- mget(genes, org.Hs.egSYMBOL2EG, ifnotfound=NA) #找出基因对应的id
entrezIDs <- as.character(entrezIDs)
out=cbind(rt,entrezID=entrezIDs)
write.table(out,file="id.txt",sep="\t",quote=F,row.names=F) #输出结果
##读取ID转换后文件
rt=read.table("id.txt",sep="\t",header=T,check.names=F) #读取id.txt文件
rt=rt[is.na(rt[,"entrezID"])==F,] #去除基因id为NA的基因
gene=rt$entrezID
2. DO疾病本体分析
library(DOSE) #加载需要使用的包
erich.do<-DOSE::enrichDO(gene=gene,
ont = "DO",
pvalueCutoff = 0.05,
qvalueCutoff = 0.05,
readable = T) #DO代码,其中P值和Q值可以自己设置,自己按照需求设置,我这里设为0.05
#########下面是画图
barplot(erich.do) #柱状图如下,颜色越红,疾病越显著
dotplot(erich.do) #点图如下,点越大、颜色越红,相应疾病的基因所占比越显著重要
#将相应疾病以网格的形式进行展示
library(ggnewscale)
cnetplot(erich.do)
DO分析到这里就结束了,建议小伙伴们图形导出格式选择TIFF格式,质量高,并且SCI投稿合并的时候合并图片更方便。下一章是PPI(蛋白互作分析)。