###go富集
对populus_simonii_genome_data.emapper.annotations处理
vi chang_go1.py
#!/usr/bin/python
# -*- coding: utf-8 -*-
import sys, os
x = sys.argv[1]
file = open(x, "r")
lines = file.readlines()
for line in lines:
line=line.strip()
if line.startswith("#"):
continue
else:
tmp=line.split("\t")
if tmp[12] == "-":
continue
else:
if(tmp[12].find(',')):
tmp2=tmp[12].split(',')
for i in range(0,len(tmp2)):
print(tmp[0]+"\t"+tmp2[i])
else:
print(tmp[0]+"\t"+tmp[12])
python chang_go1.py populus_simonii_genome_data.emapper.annotations >go.all
然后用R语言画图
library(clusterProfiler)
go_anno <- read.delim("go.all", header=FALSE, stringsAsFactors =FALSE)
names(go_anno) <- c('gene_id','ID')
go_class <- read.delim('go_term.txt', header=FALSE, stringsAsFactors =FALSE)
names(go_class) <- c('ID','Description','Ontology')
go_anno <-merge(go_anno, go_class, by = 'ID', all.x = TRUE)
gene_list <- read.delim("genego.list",header=FALSE,stringsAsFactors = FALSE)
names(gene_list) <- c('gene_id')
gene_select <- gene_list$gene_id
go_rich <- enricher(gene = gene_select,
TERM2GENE = go_anno[c('ID','gene_id')],
TERM2NAME = go_anno[c('ID','Description')],
pvalueCutoff = 0.05,
qvalueCutoff = 0.1,
pAdjustMethod = 'BH',
maxGSSize = 200)
df<-as.data.frame(go_rich)
pdf(file="enrichment.pdf",width=9,height=5)
dotplot(go_rich)
dev.off()
利用eggnog做go富集
最后编辑于 ：
©著作权归作者所有,转载或内容合作请联系作者
- 文/潘晓璐 我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
- 文/花漫 我一把揭开白布。 她就那样静静地躺着，像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上，一...
- 文/苍兰香墨 我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
推荐阅读更多精彩内容
- 首先需要加载的包 需要准备的文件 包含gene列表的txt或csv文件,如“genelist.txt” 格式为EN...
- 自备的GO库，使用R clusterProfiler做GO富集非模式物种和数据库中没有的物种要做GO分析怎么办呢...
- 葡萄的参考基因组下载自NCBI，下载链接是https://ftp.ncbi.nlm.nih.gov/genomes...
- 今天新学会了clusterprofiler，把学习过程记录下来。Clusterprofiler是Y叔开发的一个R包...
- 找到了一个在线工具，直接上传 geneid 就可以 ,链接是 http://wheat.cau.edu.cn/TG...
