- 前言
我们在做转录组数据分析时,大多数都会进行功能富集分析,预测目的基因所具有的的功能。富集工具常用到的R语言中clusterProfiler
包,里面包含了上千个功能富集的背景数据文件,功能非常强大,目前已经更新到V4.0版本。
但是,我们使用里面clusterProfiler包中的背景数据文件,你的基因名必须是要与背景数据文件中的一致。在人类等非模式生物中基本可以满足。但对于一些模式植物,并不是那么的友好,有些背景基因数据版本并非是你想要的。以及,如果你是进行无参比对,那基因的功能富集的背景数据就需要你自己制作哦!
在agriGO数据库中下载,
网址:http://systemsbiology.cau.edu.cn/agriGOv2/index.php
1 数据前期准备
我们制作番茄基因组的背景基因数据文件,番茄注释文件ITAG4.1
,下载网址:https://solgenomics.net/ftp//tomato_genome/annotation/ITAG4.1_release/,文件大小104M (.gff),下载番茄基因组数据S_lycopersicum_chromosomes.4.00.fa
,文件大小759M (.fa)。
提取出所有基因的fa
文件
#!/bin/bash
## download the tomato reference geneome to 4.1
wget https://solgenomics.net/ftp//tomato_genome/annotation/ITAG4.1_release/ITAG4.1_gene_models.gff
wget https://solgenomics.net/ftp//tomato_genome/annotation/ITAG4.1_release/ITAG4.1_CDS.fasta
wget https://solgenomics.net/ftp//tomato_genome/assembly/build_4.00/S_lycopersicum_chromosomes.4.00.fa
gffread ITAG4.1_gene_models.gff -T -o ITAG4.1_gene_models.gtf
##
gffread ITAG4.1_gene_models.gff -T -o ITAG4.1_gene_models.gtf
##
gffread -w tomato_4.1.fa -g S_lycopersicum_chromosomes.4.00.fa ITAG4.1_gene_models.gtf\
2 Swiss-Prot蛋白数据库下载及建库
下载网址:https://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/
使用diamond makedb
建库
diamond makedb --in swissprot.gz -p 60 -d swissprot
3 diamond blastx比对
diamond blastx -d ../swissprot.dmnd -q ../Tomato_4.1/tomato_4.1.fa -k 1 -e 0.00001 -o tomato.gene.m8
我这里使用共享服务器,CUP是192个线程,运行时使用内存14.5G左右,运行时间4-5分钟左右。得到结果文件
3.1 筛选出blastx最佳的结果文件
使用脚本进行筛选
4 ID TO GO好的转换
下载idmapping数据库
https://ftp.proteininformationresource.org/databases/idmapping/idmapping.tb.gz
5 提取相关ID
5.1 提取Swissprot blast 结果ID文件
提出Blastx输出Best结果文件
使用python脚本提取tomato-gene.4.1.swissprot.m8.best
第一列和第二列结果,并更剔除第二列后面的.1
,.2
最终得到模式植物GO背景文件,结果如下所示:
功能富集分析
视频教程:
模式植物或无参基因比对 | GO背景基因文件制作 | 详细教程(必须收藏) - 知乎 (zhihu.com)
模式植物或无参基因比对 | GO背景基因文件制作 | 详细教程(公众号)
“小杜的生信筆記” 公众号、知乎、简书平台,主要发表或收录生物信息学的教程,以及基于R的分析和可视化(包括数据分析,图形绘制等);分享感兴趣的文献和学习资料!