用Annovar注释非人类基因组,如小鼠mm9
annovar一般只包含人类基因组注释数据库,其他的物种如小鼠需要自己进行建立注释信息。
第一步:下载annovar软件
上Annovar官网下载(http://annovar.openbioinformatics.org/en/latest/user-guide/download/),现在要邮件注册后才能下载。邮件注册后会给你最新版软件下载地址,下载后文件为annovar.latest.tar.gz。
第二步:安装Annovar
linux系统下用该命令解压
tar zxvf annovar.latest.tar.gz
解压后生成annovar文件夹,里面有6个perl脚本程序和两个文件夹,其中一个是example文件夹,另一个是已经建立好的hg19或者GRCh37的humandb的数据库文件夹,可用于人的注释。
第三步:使用Annovar
人的注释方法,官网介绍的很详细,但仅仅有人的数据库肯定是满足不了大家的需求。
一般如果你想看是否有某种物种,如小鼠mm9的注释库时,命令行运行
perl annotate_variation.pl -builder mm9 -downdb avdblist -webfrom annovar ./
会生成一个mm9开头的文件,里面包含小鼠mm9有多少注释数据库,然后自己可以构建一个mousedb数据库
先在annovar文件夹里面创建mousedb文件夹(名字可自取),命令mkdir mousedb
然后使用annovar文件夹下的perl程序annotate_variation.pl
perl annotate_variation.pl -downdb -buildver mm9 -webfrom annovar refGene mousedb/
这个命令能实现的是帮忙下载mm9的refGene的文件,保存在mousedb文件下,自动解压后文件名为mm9_refGene.txt。
然后程序会提示使用以下两个命令继续建库
annotate_variation.pl --buildver mm9 --downdb seq mousedb/mm9_seq retrieve_seq_from_fasta.pl mousedb/mm9_refGene.txt -seqdir mousedb/mm9_seq -format refGene -outfile mousedb/mm9_refGeneMrna.fa
同样在annovar文件下运行这两个perl程序
perl annotate_variation.pl --buildver mm9 --downdb seq mousedb/mm9_seq
通过这个命令,会在mousedb下创建文件夹mm9_seq,并且在里面下载mm9的基因组文件chromFa.tar.gz,perl程序帮忙解压后是按染色体分开的fasta格式文件。
然后继续运行perl程序
perl retrieve_seq_from_fasta.pl mousedb/mm9_refGene.txt -seqdir mousedb/mm9_seq -format refGene -outfile mousedb/mm9_refGeneMrna.fa
该程序会会在mousedb下创建mm9_refGeneMrna.fa文件,是根据mm9_refGene.txt的信息,重新构建成的老鼠转录表达基因fasta格式文件
这样老鼠mm9 annovar gene based注释库就弄好了
以文本文件test.input为案例进行测试
生成test.input的txt格式文件,根据annovar官网介绍,只要这最基本的五列信息就可以进行注释,五列分别染色体名称,染色体上的位置,染色体上的位置,参考基因组碱基,变异碱基。
1 19215217 19215217 T C 1 33803084 33803084 A G 1 33803198 33803198 A G 1 37499237 37499237 T C 1 37499238 37499238 T C 1 37500003 37500003 T C 1 43826936 43826936 T C 1 58853960 58853960 A G 1 58854487 58854487 A G 1 60436865 60436865 T C
然后使用perl程序进行gene based的注释
perl annotate_variation.pl -out test -build mm9 test.input mousedb
注释后会生成test.variant_function,test.exonic_variant_function和test.log文件,前两个即为所需要的文件。用这个例子输出test.exonic_variant_function文件输出为空文件,因为这些位点没有在exonic区域的,所以没有结果。如果有位点在exonic中,则在test.exonic_variant_function中会更具体的描述为同义突变还是非同义突变
intronic Tfap2b 1 19215217 19215217 T C UTR3 Bag2 1 33803084 33803084 A G UTR3 Bag2 1 33803198 33803198 A G UTR3 Mgat4a 1 37499237 37499237 T C UTR3 Mgat4a 1 37499238 37499238 T C UTR3 Mgat4a 1 37500003 37500003 T C intronic Uxs1 1 43826936 43826936 T C intronic Casp8 1 58853960 58853960 A G intronic Casp8 1 58854487 58854487 A G intronic Cyp20a1 1 60436865 60436865 T C