差異基因富集分析及作圖

前面介紹了怎麼分析差異基因，那麼肯定需要對差異基因進行功能富集分析，今天simplify老師來給大家介紹一下怎麼做差異基因富集分析。

 

這裡利用KOBAS 3.0進行富集分析，這個軟體使用起來很簡單，幾乎看一眼就會，而且它3.0版本2019年進行了更新，裡面收錄的資料庫也比較全比較新。

1、進入網站

http://kobas.cbi.pku.edu.cn/kobas3/genelist/

2、選擇物種（人或小鼠）

3、把差異基因gene symbol複製粘貼進去，這裡可以提供多種ID

4、選擇用來富集的庫，可以直選KEGG一個庫，選多個庫的話富集結果會有很多重複的功能項目。建議pathway和GO分開運行分析，這樣便於區分和後面的作圖。

5、點擊run提交，等待一會會出現下面界面

6、點擊download下載結果文件

7、結果表格如下

8、利用R畫氣泡圖

注意：前面導出的文件開頭和末尾會有一些#注釋的行和空行，記得先去掉再保存。另外Term前面的#取掉。表頭的input number改成gene_number，P-Value改成PValue，Corrected pvalue改成Qvalue。要不然後面讀取文件會報錯。

在這裡我將ptahway和GO富集的修改後的文件分別命名為：pathway_enrichment.txt和GO_enrichment.txt，方便後面腳本讀入。

繪圖代碼如下，直接複製粘貼另存為enrichment_plot.R：

library(ggplot2)
kegg = read.table("pathway_enrichment.txt",header=TRUE,sep="\t",check.names = FALSE)kegg = kegg[1:30,]
g = ggplot(kegg,aes(-1*log10(PValue),Term))+geom_point(aes(size=gene_number,color=-1*log10(Qvalue)))g = g+scale_color_gradient(low="blue",high ="red")g = g+labs(color=expression(-log[10](Qvalue)),size="Genes",x="-log[10]Pvalue",y="Pathway",title="Pathway enrichment")ggsave("pathway_inrichment.pdf", width=6, height=6)
go = read.table("GO_enrichment.txt",header=TRUE,sep="\t",check.names = FALSE)go = go[1:30,]
g = ggplot(go,aes(-1*log10(PValue),Term))+geom_point(aes(size=gene_number,color=-1*log10(Qvalue)))g = g+scale_color_gradient(low="blue",high ="red")g = g+labs(size="Genes",x="-log10(Pvalue)",y="GO term",title="GO enrichment")ggsave("GO_inrichment.pdf", width=7, height=6)
最後直接運行代碼：
Rscript enrichment_plot.R
9、富集結果
Pathway富集結果
GO富集結果
今天就分析到這了，下期繼續。
往期精彩
轉錄組系列分析(一)：一行命令完成RNAseq差異分析+畫圖
轉錄組系列分析(二)：一行命令完成RNAseq數據標記重要基因
轉錄組系列分析(三)：如何利用網絡圖給自己轉錄組數據的文章添彩