口腔医学小站

口腔文献阅读与资源分享

R做GSEA富集分析

2021年3月22日 发表者 kimi1006 | 暂无评论

R做GSEA富集分析

22020.06.30 16:30:37字数 577阅读 7,023

首先感谢Y叔的clusterprofiler神包，做富集分析优点是在线爬取数据，结果很可信，但是缺点也是网络问题，网络差点就要等很久，不过GSEA有自带GMT文件，因此下载好离线数据，这些就可以摆脱在线的问题，单机就可以操作GSEA了
GSEA也就是基因集富集分析，不需要分析是上调基因还是下调基因，分析的是所有基因，所以结果应该更可靠，根据官方操作说明， 需要有两组数据，第一组是表达矩阵，第二组是分组，然后用软件操作，但是至始至终出的图奇丑无比，而且还只能是png格式，远达不到300dpi的发表级，虽然目前有很多再次作图方法，但依然比较曲折，Y叔的包解决了这个问题，下面分享一下教程

  • 第一步：数据格式
    表格需要两列， 第一列是基因名，第二列是Foldchange或者Log2FC
    以前我一直不明显第二列是如何来的，明明是矩阵和分组，怎么算Fc，但是也怪自己学艺不精，其实有列矩阵和分组就可以算出来了，懂R的可以用各种差异分析包如“DESq”，“limma”和“edgeR”算出来，不懂R的用Excel也可以，可以看我之前的教程，其实Foldchange就是两组数据的均值之比，也就是先算出各组的平均值，然后拿比较组除以对照组就算出来了 (没事还是要多读书啊)
  • 第二步，如基因名是symbol，需要基因ID转换为entrezid格式
gene<-read.csv("你的文件.csv")  #csv可读性较好，带表头，需要有一列是symbol
library(clusterProfiler)
library(org.Hs.eg.db)
library(stringr)
gene<-str_trim(d$symbol,"both") #定义gene
#开始ID转换
gene=bitr(gene,fromType="SYMBOL",toType="ENTREZID",OrgDb="org.Hs.eg.db") #会有部分基因数据丢失，或者ENSEMBL
## 去重
gene <- dplyr::distinct(gene,SYMBOL,.keep_all=TRUE)
gene_df <- data.frame(logFC=gene$logFC, #可以是foldchange
SYMBOL = gene$symbol) #记住你的基因表头名字
gene_df <- merge(gene_df,gene,by="SYMBOL")
24681012141618202224
  • 第三步，定义基因列表，对logFC进行从高到低排序
geneList<-gene_df $logFC #第二列可以是folodchange，也可以是logFC
names(geneList)=gene_df $ENTREZID #使用转换好的ID
geneList=sort(geneList,decreasing = T) #从高到低排序
246

  • 第四步，下载离线GMT文件
    https://www.gsea-msigdb.org/gsea/downloads.jsp （需要先注册，不过一个邮箱地址即可）
    包含GO，KEGG，REACTOME，HALLMARKS等等，需要什么下载什么

    图片.png

  • 第五步，读取GMT格式，直接GSEA分析并出图，以Kegg为例

 kegmt<-read.gmt("c2.cp.kegg.v7.1.entrez.gmt") #读gmt文件
 KEGG<-GSEA(geneList,TERM2GENE = kegmt) #GSEA分析
 library(ggplot2)
 dotplot(KEGG) #出点图 
dotplot(KEGG,color="pvalue")  #按p值出点图 
246810
默认p.ajust

p值
dotplot(KEGG,split=".sign")+facet_grid(~.sign) #出点图，并且分面激活和抑制
2
图片.png
 dotplot(KEGG,split=".sign")+facet_wrap(~.sign,scales = "free") #换个显示方式
2
图片.png
 library(enrichplot)
#特定通路作图
 gseaplot2(KEGG,1,color="red",pvalue_table = T) # 按第一个做二维码图，并显示p值
246
图片.png
 gseaplot2(KEGG,1:10,color="red") #按第一到第十个出图，不显示p值
2
图片.png

23人点赞

更多精彩内容下载简书APP
“小礼物走一走，来简书关注我”

还没有人赞赏，支持一下

欧阳松石河子大学医学院第一附属医院泌尿外科主治医师<br>华中科技大学同济医学院附属同济医院泌尿外科…
总资产7 (约0.54元)共写了7433字获得174个赞共186个粉丝

点击数：0

发表回复

*为必填字段！