基因组变异注释-利用annovar软件

崔耳又又

ResBang生信交流网官方公众号的聊天框可以对网站内容快速检索。有位读者朋友对“annovar”进行检索，发现没有相关文章，所以我赶紧更新了一期有关基因组变异结构注释的文章。

通常变异检测后会得SNP和Indel变异数据，通过注释可以了解具体哪些等位变异会改变基因结构，哪些发生在基因启动子上，这些信息对下游基因的筛选非常重要。话不多说，我们直接开始教程。

1 软件安装

需要两个软件：gff3ToGenePred和annovar

（1）gff3ToGenePred安装

gff3ToGenePred可直接通过conda安装

conda install bioconda::ucsc-gff3togenepred

（2）annovar 安装

我这边直接给到大家一个安装包，大家下载后传到服务器解压即可直接使用

百度网盘链接：https://pan.baidu.com/s/1CELFS716wGZQQOw6uwZhgA?pwd=2024

提取码：2024

2 对基因组建索引

这一步我们需要对参考基因组fasta文件和参考基因组gff3注释文件建索引

参考基因组fasta文件：REF.fa

参考基因组gff3注释文件：REF.gff3

（1）建一个文件夹用于存放索引文件

mkdir annolib

（2）参考基因组gff3注释文件建索引

gff3ToGenePred REF.gff3 annolib/reference_refGene.txt 
#只需将REF.gff3更改为自己的参考基因组gff3注释文件即可

（3）参考基因组fasta文件建索引

perl annovar/retrieve_seq_from_fasta.pl --format refGene --seqfile REF.fa annolib/reference_refGene.txt --out annolib/reference_refGeneMrna.fa
#只需将REF.fa更改为自己的参考基因组fasta文件即可

需要注意的是这一步非常容易报错，大家可根据报错信息排除错误，也可在留言区或者ResBang生信交流网问答专区写下你的问题。

3 将vcf文件修改为annovar软件输入格式

annovar软件输入格式如下：

这一步通过awk即可实现，我这边以SNP.vcf为例

grep -v '##' snp.vcf | awk '{print $1"\t"$2"\t"$2"\t"$4"\t"$5"\t"$4$4"\t"$5$5}' > snp.avinput

4 进行变异注释

直接运行以下代码

perl annovar/annotate_variation.pl \
        -buildver reference \
        snp.avinput \
        annolib/ \
        -geneanno  \
        -dbtype refGene\
        -outfile snp.anno
#snp.avinput为输入文件，snp.anno为输出文件前缀

5 结果解读

结果会输出4个文件：

snp.anno.variant_function，snp.anno.log，snp.anno.invalid_input，snp.anno.exonic_variant_function

（1）snp.anno.variant_function文件为所有snp的变异注释

每一列的信息如下：

A：发生变异的类型

B：发生变异的基因

C：发生变异的染色体

D：发生变异的起始位置

E：发生变异的终止位置

F：REF基因型

G：ALT基因型

H：样本1的基因型

I：样本2的基因型

（2）snp.anno.exonic_variant_function文件为所有发生在外显子上的变异注释

每一列的信息如下：

A：输入文件中第几行的变异，注意输入文件的表头也算一行

B：突变类型，如“同义突变”、“非同义突变”、“提前终止”等

C：发生变异的基因，“AetAY61_1Dv1G0000100:AetAY61_1Dv1G0000100.1:exon1:c.G108C:p.L36K”的意思为发生变异的基因为AetAY61_1Dv1G0000100，发生变异的转录本为AetAY61_1Dv1G0000100.1，发生在第一个外显子上，是第108个碱基的G突变位C，导致第36个氨基酸的L突变为K

D-J与上一个文件的C-I一致

微信扫一扫分享文章