今天看到一个作业,直接找到教程开始学习:
为R包写一本书(像Y叔致敬)
拿到的基因列表是symbol,但是这个包里面用超几何分布检验的模型弄成的一个对象enrichKK所需要的基因名是Expected input gene ID: 8310,1375,230,5634,2539,8790
这种类型的。
因此去菜鸟团上搜索了一下,曾老师16年的时候有个教程。
用R的shiny包写一个基因的ID转换小程序
我太懒了,直接用一个网页工具来获得gene ID
id转换网页工具
殊途同归,拿到ID就好。
2019-12-29 修改
其实这里可以通过代码来做转换的。
我昨天的内容就不做转换了,因为有小伙伴已经转了跟我用网页作出来的结果差不多,那么代码怎么实现呢,我今天用了另外一个基因列表来做转换因为是自己的课题就不写的后面的结果,前面ID转换需要的包及代码如下:
rm(list = ls())
options(stringsAsFactors = F)
library(ggplot2)
library(clusterProfiler)
library(org.Hs.eg.db)
deg=read.csv(file = "IP_MS_TF_list.csv")
head(deg)
colnames(deg)=c("SYMBOL","ENTRIZID")
deg2=deg[,-3]
ls("package:org.Hs.eg.db")
df <- bitr(unique(deg2$SYMBOL), fromType = "SYMBOL",
toType = c( "ENTREZID"),
OrgDb = org.Hs.eg.db)
df2 <- bitr(unique(deg2$ENTRIZID), fromType = "ENSEMBL",
toType = c( "ENTREZID"),
OrgDb = org.Hs.eg.db)
我的这个列表如果单纯以SYMBOL来转换呢会丢失一些基因信息,而我用了对应的ENSEMBL基因名来转换就不会丢失。这个是由于SYMBOL和ENTREZID的数据库不完全对应。这里需要注意!
拿到转换后的ID列表保存开始按照教程完成任务:
rm(list = ls())
options(stringsAsFactors = F)
library(ggplot2)
library(clusterProfiler)
library(org.Hs.eg.db)
a=read.csv("down-gen.csv")
b=read.csv("up-gene.csv")
gene_up=b$converted_alias
gene_down=a$converted_alias
enrichKK <- enrichKEGG(gene = gene_up,
organism = 'hsa',
#universe = gene_all,
pvalueCutoff = 0.1,
qvalueCutoff =0.1)
head(enrichKK)[,1:6]
browseKEGG(enrichKK, 'hsa04512')
dotplot(enrichKK)
enrichKK=DOSE::setReadable(enrichKK, OrgDb='org.Hs.eg.db',keyType='ENTREZID')
enrichKK
这里得到2个图,一个是KEGG通路https://www.kegg.jp/kegg-bin/show_pathway?hsa04512/NA
然后是气泡图:
2,富集分析
这2个图出现了说明找到的前面的方法是对的,那就开始Y书的表演
最基础的条形图和点图
#条带图
barplot(enrichKK,showCategory=20)
#气泡图
dotplot(enrichKK)
图1
图2
通路与基因之间的关系可视化
cnetplot(enrichKK, categorySize="pvalue", foldChange="qvalue",colorEdge = TRUE)
cnetplot(enrichKK, foldChange="qvalue", circular = TRUE, colorEdge = TRUE)
这里我把原来那个genelist改成了qvalue,我也不知道对不对,但是如果不改会报错的。
图3-1
图3-2
通路与通路之间的连接展示
emapplot(enrichKK)
图4
热图展现通路与基因之间的关系
heatplot(enrichKK)
图5
好了,这个包的功能真强大,主要是画的图好看。
marker一下,第一次完成作业。
网友评论