结果文件的解读

输出文件1：*.variant_function

$head  ZH13_vs_JYH.refGene.variant_function
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=50203)   GWHAAEV00000001 18052   18052   C       T       1       30      .       GWHAAEV0000000118052    .       C       T       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=47215)   GWHAAEV00000001 21040   21040   A       C       1       30      .       GWHAAEV0000000121040    .       A       C       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=47096)   GWHAAEV00000001 21159   21159   G       A       1       30      .       GWHAAEV0000000121159    .       G       A       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=47056)   GWHAAEV00000001 21199   21199   A       G       1       30      .       GWHAAEV0000000121199    .       A       G       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=47002)   GWHAAEV00000001 21253   21253   A       T       1       30      .       GWHAAEV0000000121253    .       A       T       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=46432)   GWHAAEV00000001 21823   21823   T       C       1       30      .       GWHAAEV0000000121823    .       T       C       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=46279)   GWHAAEV00000001 21976   21976   -       AG      1       30      .       GWHAAEV0000000121976    .       C       CAG     30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=46113)   GWHAAEV00000001 22141   22142   TA      -       1       30      .       GWHAAEV0000000122140    .       ATA     A       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=46016)   GWHAAEV00000001 22239   22239   G       T       1       30      .       GWHAAEV0000000122239    .       G       T       30      PASS    .       GT      1/1
intergenic      NONE(dist=NONE),SoyZH13_01G000100(dist=45816)   GWHAAEV00000001 22439   22439   G       A       1       30      .       GWHAAEV0000000122439    .       G       A       30      PASS    .       GT      1/1

第一个文件包含所有变异的注释，方法是在每个输入行的开头添加两列（例如，第一行下面的“intergenic NONE(dist=NONE),SoyZH13_01G000100(dist=50203)”）：

第一列告诉该变异是否命中外显子、基因间区域、内含子或、非编码RNA基因。如果变异是外显子/内含子/ ncRNA，第二列给出了基因名称（如果命中了多个基因，则在基因名称之间添加逗号）。如果不是，第二列将给出两个相邻基因以及到这些相邻基因的距离。

第一列的可能值汇总如下：

Value	默认优先级	说明	Sequence Ontology
exonic	1	变异与编码区重叠	exon_variant（SO：0001791）
splicing	1	变异在剪接位点的2 bp之内（使用-splicing_threshold进行更改）	splicing_variant（SO：0001568）
ncRNA	2	变异与转录本重叠，而在基因定义中没有编码注释（更多说明请参见下面的注释）	non_coding_transcript_variant（SO：0001619）
UTR5	3	变异与5'非翻译区重叠	5_prime_UTR_variant（SO：0001623）
UTR3	3	变异与3'非翻译区重叠	3_prime_UTR_variant（SO：0001624）
intronic	4	变异与内含子重叠	intron_variant（SO：0001627）
upstream	5	变异在转录起始位点上游重叠1-kb区	upstream_gene_variant（SO：0001631）
downstream	5	变异与转录终止位点下游1-kb区（使用-neargene进行更改）	downstream_gene_variant（SO：0001632）
intergenic	6	变异位于基因间区域	intergenic_variant（SO：0001628）

第一列的值具有以下优先级（截至2010年12月及更高版本的ANNOVAR）：
$\color{red}{exonic = splicing > ncRNA> UTR5/UTR3 > intron > upstream/downstream > intergenic}$
。上面定义的优先级用于确定当该变异适合多个功能类别时要打印什么功能。注意：

exonic仅指编码外显子部分，而不是UTR部分，因为有两个专门为UTR注释保留的关键字（UTR5，UTR3）。
默认情况下，ANNOVAR中的splicing定义为距离外显子/内含子边界2 bp以内的变异，但是可以通过--splicing_threshold参数更改阈值。在2013年2月之前，如果显示“ exonic，splicing”，则表示这是外显子内的一个变异，但接近外显子/内含子边界；此行为是由于历史原因引起的，当时用户要求在splicing位点附近的外显子变异也要加上splicing注释。从2013年2月开始，“splicing”仅指内含子中的2bp，它在外显子附近，如果要具有与以前相同的行为，请添加-exonicsplicing参数。
如果一个变体位于5'UTR和3'UTR区域中（可能是两个不同的基因），则将输出“ UTR5，UTR3”。
考虑到mRNA的链，术语upstream和downstream分别定义为距转录起始位点或转录终止位点1kb。--neargene阈值可用于调整此阈值。如果一个变异同时位于下游和上游区域（可能是两个不同的基因），那么upstream,downstream将被输出。在2011年6月版的ANNOVAR中，对拼接注释进行了改进。如果剪接位点位于内含子，则将打印所有 isoforms和相应的碱基变化。例如，

splicing SMS(NM_004595:c.447+2T>G) X 21895357 21895357 T G hetero 8 15
splicing DMD(NM_004011:c.48+1A>C) X 31803228 31803228 T G homo 117 30
splicing BAGE(NM_001187:c.14+1A>G),BAGE4(NM_181704:c.14+1A>G),BAGE5(NM_182484:c.14+1A>G) 21 10120594 10120594 T C hetero 66 53

下面讨论了几个技术说明。

技术说明：上面的ncRNA是指没有编码注释的RNA。这并不意味着这是一个永远不会被翻译的RNA。这仅意味着用户选择的基因注释系统无法提供编码序列注释。它仍然可以编码蛋白质，并且在将来的基因注释版本或其他基因注释系统中可能具有此类注释。例如，当使用UCSC已知基因注释时，BC039000被ANNOVAR视为ncRNA，但是当使用ENSEMBL注释时，其被ANNOVAR视为蛋白质编码基因。如果用户的目标是找到已知的（经过良好注释的）microRNA或其他已知的（经过良好注释的）非编码RNA，则应使用基于区域的注释，并应选择wgRNA轨道。在此处阅读说明。

技术说明：如果删除了转录本的第一个密码子，则由于该基因无法翻译，因此将被ANNOVAR报告为全基因缺失(wholegene deletion)。

如果用户希望打印出所有功能性结果（而不仅仅是上面优先级定义的最重要的结果），则应使用参数--separate。在这种情况下，每个变异可能会出现几条输出线，代表几种可能的功能影响。

要进一步说明variant_function注释，请检查下图：

在上图中，SNP1是一个基因间变异，因为它与任何基因相距> 1kb； SNP2是下游变异，因为它与NADK基因的3'端相距1kb；SNP3是UTR3变异；SNP4是内含子。SNP5是外显子变异。

类似地，del1是一种基因间变异。del2是下游变异；del3是UTR3变异；del4与UTR3和内含子都重叠，并且基于优先规则，它是UTR3变异。del5是内含子变异；del6与外显子和内含子都重叠，并且基于优先规则，它是外显子变体。

技术说明：有时用户可能想更改默认优先级规则。-precedence参数可用于微调变量功能的优先级。不同的变异功能应根据所需的优先级在命令行中用逗号分隔。变异功能的允许关键字是exonic, intronic, splicing, utr5, utr3, upstream, downstream, splicing, ncrna。

例如，当-precedence intronic,utr5,utr3指定时，内含子变异将优先于UTR，而del4将成为上面的内含子变异。这是唯一的更改，所有其他默认优先级规则仍然适用于此。

技术说明：默认情况下，基因名称打印在variant_function文件的第二栏中。有时，用户可能希望改为查看transcript 名称。--transcript_function参数可用于指定此行为。请注意，很有可能在输出中打印多个转录名称，以逗号分隔，因为每个基因名称通常对应于多个转录名称。

技术说明：基因定义的当前逻辑：基因组很复杂，因此基因定义也很复杂。ANNOVAR完全依靠用户提供的基因定义（例如RefSeq，UCSC基因和Ensembl基因）将转录本映射到基因组并将转录本与基因相关联，并使用以下逻辑处理复杂的情况：

1.如果将一个基因注释为编码和非编码（multiple transcripts, some coding, some non-coding），该基因都将被视为编码基因（非编码转录本将被忽略）。

2.如果一个基因或转录本具有一个或几个非编码但没有编码，则在注释输出中将其视为ncRNA。

3.如果一个转录本作为coding transcripts映射到多个位置，但是有些具有完整的ORF，有些则没有完整的ORF（premature stop codon），则没有完整ORF的密码将被忽略。

4.如果一个转录本映射到多个位置，全部都作为“coding transcripts”，但没有一个完整的ORF，则此转录本将不会在exonic_variant_function注释中使用，并且相应的注释将标记为“ UNKNOWN”。

2014年7月新增功能：如果转录本映射到多个基因组位置，则所有映射都将在注释过程中使用。以前，注释中仅使用“最可能的”映射。

以上规则确实有意义。规则3和4是在2011年11月的ANNOVAR版本中制定的（这样，用户不再向我发送抱怨外显子注解错误的电子邮件，这并不是ANNOVAR本身的错）。但是结果是，您可能在输出文件中看到更多的“ UNKNOWN”外显子注释。

输出文件2（refSeq gene annotation）

第二个输出文件ex1.exonic_variant_function包含由于外显子变异导致的氨基酸变异。在不同版本的ANNOVAR之间，以下输出的确切格式可能会略有不同。

[kaiwang@biocluster ~/]$ cat ex1.exonic_variant_function 
line145 synonymous SNV  SoyZH13_01G000100:SoyZH13_01G000100.m1:exon1:c.A81G:p.L27L,     GWHAAEV00000001 68335   68335   A       G       1       30     .GWHAAEV00000001 68335   .       A       G       30      PASS    .       GT      1/1
line148 synonymous SNV  SoyZH13_01G000100:SoyZH13_01G000100.m1:exon3:c.T216C:p.T72T,    GWHAAEV00000001 68804   68804   T       C       1       30     .GWHAAEV00000001 68804   .       T       C       30      PASS    .       GT      1/1
line149 synonymous SNV  SoyZH13_01G000100:SoyZH13_01G000100.m1:exon3:c.C258T:p.G86G,    GWHAAEV00000001 68846   68846   C       T       1       30     .GWHAAEV00000001 68846   .       C       T       30      PASS    .       GT      1/1
line150 synonymous SNV  SoyZH13_01G000100:SoyZH13_01G000100.m1:exon3:c.G288A:p.K96K,    GWHAAEV00000001 68876   68876   G       A       1       30     .GWHAAEV00000001 68876   .       G       A       30      PASS    .       GT      1/1
line163 frameshift deletion     SoyZH13_01G000200:SoyZH13_01G000200.m1:exon1:c.29_36del:p.R11Wfs*17,    GWHAAEV00000001 71985   71992   AGAGAGAG       -1       30      .       GWHAAEV00000001 71984   .       GAGAGAGAG       G       30      PASS    .       GT      1/1
line164 nonsynonymous SNV       SoyZH13_01G000200:SoyZH13_01G000200.m1:exon1:c.A274G:p.I92V,    GWHAAEV00000001 72230   72230   A       G       1      30       .       GWHAAEV00000001 72230   .       A       G       30      PASS    .       GT      1/1
line165 nonsynonymous SNV       SoyZH13_01G000200:SoyZH13_01G000200.m1:exon1:c.T275C:p.I92T,    GWHAAEV00000001 72231   72231   T       C       1      30       .       GWHAAEV00000001 72231   .       T       C       30      PASS    .       GT      1/1
line166 synonymous SNV  SoyZH13_01G000200:SoyZH13_01G000200.m1:exon1:c.C372T:p.D124D,   GWHAAEV00000001 72328   72328   C       T       1       30     .GWHAAEV00000001 72328   .       C       T       30      PASS    .       GT      1/1
line167 nonsynonymous SNV       SoyZH13_01G000200:SoyZH13_01G000200.m1:exon2:c.G437A:p.R146K,   GWHAAEV00000001 72477   72477   G       A       1      30       .       GWHAAEV00000001 72477   .       G       A       30      PASS    .       GT      1/1
line168 synonymous SNV  SoyZH13_01G000200:SoyZH13_01G000200.m1:exon2:c.A534G:p.G178G,   GWHAAEV00000001 72574   72574   A       G       1       30     .GWHAAEV00000001 72574   .       A       G       30      PASS    .       GT      1/1

请注意，此文件中仅注释了外显子变异，因此：

第一列给出了原始输入文件中的行号。
第二个字段说明了该变异的功能后果（此字段中的可能值包括：非同义SNV，同义SNV，移码插入，移码删除，非移码插入，非移码删除，移码块替换，非移码块替换）。
第三列包含基因名称，转录本标识符和相应转录本中的序列变化。在指定序列更改时使用了标准命名法（您可能需要添加-hgvs参数，以便cDNA级别注释与HGVS命名法兼容）。

这些exonic_variant_functoin annotation的详细说明如下。请注意，stopgain 和stoploss优先于其他注释；例如，每当非同义突变将野生型氨基酸更改为终止密码子时，它将被标注为stopgain，而不是非同义SNV。

注释	优先顺序	说明	Sequence Ontology
frameshift insertion	1	插入一个或多个核苷酸导致蛋白质编码序列移码	frameshift_elongation（SO：0001909）
frameshift deletion	2	一个或多个核苷酸的缺失导致蛋白质编码序列移码	frameshift_truncation（SO：0001910）
frameshift block substitution	3	一个或多个核苷酸取代导致蛋白质编码序列移码	frameshift_variant（SO：0001589）
stopgain	4	非同义SNV，移码插入/删除，非移码插入/删除或块替换，导致在变异位点产生终止密码子。对于移码突变，变异下游的终止密码子的产生将不算作“stopgain”！	stop_gained（SO：0001587）
stoploss	5	非同义SNV，移码插入/删除，非移码插入/删除或块替换，可导致在变异位点消除终止密码子	stop_lost（SO：0001578）
nonframeshift insertion	6	插入3个或3个核苷酸的倍数，不会引起蛋白质编码序列移码变化	inframe_insertion（SO：0001821）
nonframeshift deletion	7	删除3个核苷酸或3倍数个核苷酸，不会引起蛋白质编码序列的移码变化	inframe_deletion（SO：0001822）
nonframeshift block substitution	8	一个或多个不会引起蛋白质编码序列移码变化的核苷酸的替换	inframe_variant（SO：0001650）
nonsynonymous SNV	9	导致氨基酸变化的单个核苷酸变化	missense_variant（SO：0001583）
synonymous SNV	10	不引起氨基酸变化的单个核苷酸变化	synonymous_variant（SO：0001819）
unknown	11	功能未知（由于数据库文件中基因结构定义的各种错误）	sequence_variant（SO：0001060）

在2012年10月版的ANNOVAR中，--aamatrixfile添加了自变量，以便用户可以打印出GRANTHAM分数（或任何其他氨基酸替代矩阵）以用于基于基因的非同义变异注释。参见下文，将AAMatrix=43符号添加到输出中，指示R-> Q更改的Grantham得分为43。

[kaiwang@biocluster ]$ annotate_variation.pl example/ex1.avinput humandb/ -aamatrixfile grantham.matrix -out ex1 -buildver hg19 
[kaiwang@biocluster ]$ head -n 1 ex1.exonic_variant_function
line9 nonsynonymous SNV IL23R:NM_144701:exon9:c.G1142A:p.R381Q:AAMatrix=43, 1 67705958 67705958 G A comments: rs11209026 (R381Q), a SNP in IL23R associated with Crohn's disease

氨基酸变化与转录物的位置有关（而非“基因”）。例如，R702W突变是指称为NM_022162的转录本（对应于NOD2基因）在外显子4位置702处的氨基酸变化。由于只有一个转录本带有NOD2基因注释，因此这里没有歧义。但是，由于可变剪接，如果一个基因有两个或两个以上的转录本全部注释，则氨基酸变化的位置将有所不同，除基因名称外，始终列出这些转录本很重要。

如果用户有兴趣将HGVS命名法用于cDNA，请在基因注释中添加参数-hgvs：

[kaiwang@biocluster ~/]$ annotate_variation.pl -out ex1 -build hg19 -hgvs example/ex1.avinput humandb/ 

[kaiwang@biocluster ~/]$ cat ex1.exonic_variant_function
line9 nonsynonymous SNV IL23R:NM_144701:exon9:c.1142G>A:p.R381Q, 1 67705958 67705958 G A comments: rs11209026 (R381Q), a SNP in IL23R associated with Crohn's disease
line10 nonsynonymous SNV ATG16L1:NM_030803:exon9:c.898A>G:p.T300A,ATG16L1:NM_017974:exon8:c.841A>G:p.T281A,ATG16L1:NM_001190267:exon9:c.550A>G:p.T184A,ATG16L1:NM_001190266:exon9:c.646A>G:p.T216A,ATG16L1:NM_198890:exon5:c.409A>G:p.T137A, 2 234183368 234183368 A G comments: rs2241880 (T300A), a SNP in the ATG16L1 associated with Crohn's disease
line11 nonsynonymous SNV NOD2:NM_001293557:exon3:c.2023C>T:p.R675W,NOD2:NM_022162:exon4:c.2104C>T:p.R702W, 16 50745926 50745926 comments: rs2066844 (R702W), a non-synonymous SNP in NOD2
line12 nonsynonymous SNV NOD2:NM_001293557:exon7:c.2641G>C:p.G881R,NOD2:NM_022162:exon8:c.2722G>C:p.G908R, 16 50756540 50756540 comments: rs2066845 (G908R), a non-synonymous SNP in NOD2
line13 frameshift insertion NOD2:NM_001293557:exon10:c.2936dupC:p.A979fs,NOD2:NM_022162:exon11:c.3017dupC:p.A1006fs, 16 50763778 5076377comments: rs2066847 (c.3016_3017insC), a frameshift SNP in NOD2
line14 frameshift deletion GJB2:NM_004004:exon2:c.35delG:p.G12fs, 13 20763686 20763686 G - comments: rs1801002 (del35G), a frameshift mutation in GJB2, associated with hearing loss
line15 frameshift deletion GJB6:NM_001110220:wholegene,GJB6:NM_001110221:wholegene,GJB6:NM_001110219:wholegene,CRYL1:NM_015974:wholegene,GJB6:NM_006783:wholegene, 13 20797176 21105944 0 - comments: a 342kb deletion encompassing GJB6, associated with hearing loss

技术说明：与variant_function文件类似，该exonic_variant_function文件也遵循优先级规则，但是用户无法更改此规则（无论如何，没有太多生物学上的理由来更改此规则）。例如，突变“ chr7 140453136 140453136 A T”将由ANNOVAR使用ENSEMBL定义标注为stoploss突变X208R，因为stoploss优先于非同义突变（对于ENST00000288602为V600E，对于ENST00000479537为V28E）。如果用户希望获得全面的exonic_variant_function输出，请使用-separate参数！

技术说明：某些基因具有多个转录本，ANNOVAR可能会随机排列输出文件中转录本的顺序。例如，同一突变可以从一个输入文件中注释为“ BRAF：ENST00000288602：exon15：c.T1799A：p.V600E，BRAF：ENST00000479537：exon2：c.T83A：p.V28E”，但另一个文件注释为“ BRAF：ENST00000479537” ：exon2：c.T83A：p.V28E，BRAF：ENST00000288602：exon15：c.T1799A：p.V600E”。一些用户希望注释中的绝对一致性。在这种情况下，您可以-exonsort在命令行中添加参数，以便在输出文件中，exon2始终位于exon15之前。

技术说明：与简单的“移码突变”注释相反，许多用户要求在知道是插入缺失后才知道确切的“新蛋白质序列”。我无法在ANNOVAR中直接解决此问题。为了处理这种情况，我实现了一个新脚本，该脚本从基因注释中获取输出，然后重新计算野生型和突变的蛋白质序列，并推断插入或缺失取代是否引起了蛋白质序列stopgain，stoploss或nonsynonymous 。该脚本coding_change.pl在ANNOVAR软件包中。试试看，看看它是如何工作的！

技术说明：在ANNOVAR的早期版本中，所有外显子注释均基于用户指定的基因定义和用户指定的FASTA序列。但是，这可能会带来一些问题：某些基因定义可能导致终止密码子不正确/不完整的ORF，有时与基因定义相比，FASTA序列已过时。尽管原则上用户可以通过encoding_change.pl脚本轻松地识别这些问题，但是某些用户并不想经历额外的麻烦。因此，在2011年11月的ANNOVAR版本中，我决定标识具有终止密码子的转录本，并且不再注释这些转录本的任何外显子突变（换句话说，外显子注释将标记为“未知”）。

包括RefSeq基因版本号

自2017年6月起，ANNOVAR软件包现在包含hg19_refGeneWithVer.txt文件，以举例说明如何使用带有版本的refGene注释varians。用户可以使用-dbtype refGeneWithVer而不是使用-dbtype refGene，这样结果将包含带有版本的transcript identifiers 。对于所有其他基因组构建，用户需要自己生成这些文件。

注释线粒体变异

自2013年2月起，ANNOVAR即可注释线粒体变体（只要您的染色体标识符为M或MT或chrM或chrMT，线粒体特异性密码子表将用于推断氨基酸变化）。但是，有几个重要警告：

RefSeq没有线粒体基因定义。因此，ANNOVAR用户需要使用UCSC已知基因或Ensembl基因。
UCSC的hg19程序集使用旧版本的线粒体基因组（NC_001807），但1000个基因组联合体已用最新的剑桥参考序列版本（NC_012920）取代了chrM。因此，如果您将序列数据比对并针对NC_012920调用变体，那么您将无法真正使用UCSC的基因定义来注释变体。必须坚持使用相同的坐标。为了更好地说明这一点，当您将原始序列数据作为FASTQ文件获取时，如果将数据与UCSC编译的参考基因组进行比对（通常文件名类似于hg19.fa，并且有一些名为chrx_random的染色体，那么您可以直接使用ANNOVAR一起注释所有变体。线粒体变体将是那些与chrM对齐的变体，ANNOVAR的Feb2013版本可以正确地对其进行注释。来自德国癌症研究中心的Konrad Herbst编写了一个脚本，用于在两个参考序列之间进行位置转换（由Y. Guo等人在MitoSeek软件包中由一些类似的脚本提供），并使用它将GRCh37文件映射到hg19文件。在retrieve_seq_from_fasta.plhg19参考序列（AF347015.1）上使用时，将为基于hg19的线粒体注释生成以下文件。

许多用户不知道的一种复杂情况是Ensemble的线粒体基因存在注释错误（通常有几个碱基对），因此不应使用Ensembl的基因注释。举一个简单的例子，您可以在UCSC基因组浏览器中搜索ENST00000389680：当Gencode将位置列出为chrM：650-1603时，Ensembl注释显示为chrM：646-1599，关闭了4bp。由于这些原因，在hg19 coordiante上调用变体时，应将ANNOVAR提供的文件用于任何线粒体注释。

为了使用户更容易使用，我现在在此处提供两个文件：hg19_MT_ensGene.txt和hg19_MT_ensGeneMrna.fa在ANNOVAR软件包humandb/目录中。该-buildver是hg19_MT和-dbtype是ensGene。
但是，如果将原始FASTQ文件与具有NC_012920的参考基因组进行比对（例如1000 Genomes Project提供的文件，通常文件名类似human_g1k_v37.fasta），那么您需要使用具有正确线粒体的自定义基因定义文件NC_012920的基因定义。美国国家老龄研究所的丁俊博士根据Ensemble定义慷慨地提供了这些文件。请注意，染色体名称通常应为MT（在2013年6月之前，我在文件中使用了chrM，这引起了一些ANNOVAR用户的困惑，因此我决定改用MT并坚持使用GRCh37的标准）。在这种情况下，应使用以下命令注释线粒体变体：annotate_variation.pl -buildver GRCh37_MT -dbtype ensGene mt.avinput humandb/。

为了使用户更轻松，我现在在ANNOVAR软件包的humandb /目录中提供两个文件GRCh37_MT_ensGene.txt.gz和GRCh37_MT_ensGeneMrna.fa.gz。该-buildver是GRCh37_MT和-dbtype是ensGene。

切换到UCSC / Ensembl基因注释

ANNOVAR可以选择处理UCSC已知基因注释或Ensembl基因注释，这两者都比RefSeq更为全面，因为它包含了许多注释不佳或计算预测困难的基因。下面显示了一个使用UCSC已知基因注释变异的示例：

[kaiwang@biocluster ~/]$ annotate_variation.pl -out ex1 -build hg19 example/ex1.avinput humandb/ -dbtype knownGene
NOTICE: The --geneanno operation is set to ON by default
NOTICE: Reading gene annotation from humandb/hg19_knownGene.txt ... Done with 78963 transcripts (including 18502 without coding sequence annotation) for 28495 unique genes
NOTICE: Reading FASTA sequences from humandb/hg19_knownGeneMrna.fa ... Done with 45 sequences
WARNING: A total of 43 sequences will be ignored due to lack of correct ORF annotation
NOTICE: Finished gene-based annotation on 15 genetic variants in example/ex1.avinput
NOTICE: Output files were written to ex1.variant_function, ex1.exonic_variant_function

转录本名称（在ex1.exonic_variant_function文件中）看起来像uc002eg1.1等，这是UCSC基因的标识符。

要使用Ensembl基因注释变异，请使用以下命令。输出格式类似于上述格式。“ ENSG”和“ ENST”是带注释的基因和转录本的Ensembl标识符。

[kaiwang@biocluster ~/]$ annotate_variation.pl -out ex1 -build hg19 ex1.hg19.avinput humandb/ -dbtype ensGene
NOTICE: The --geneanno operation is set to ON by default
NOTICE: Reading gene annotation from humandb/hg19_ensGene.txt ... Done with 196501 transcripts (including 101155 without coding sequence annotation) for 57905 unique genes
NOTICE: Reading FASTA sequences from humandb/hg19_ensGeneMrna.fa ... Done with 20 sequences
WARNING: A total of 6780 sequences will be ignored due to lack of correct ORF annotation
NOTICE: Finished gene-based annotation on 15 genetic variants in example/ex1.avinput
NOTICE: Output files were written to ex1.variant_function, ex1.exonic_variant_function

由于输出仅包含Ensembl标识符，因此，如果您要将其翻译为 gene synonym，则可以下载此文件以获取hg19，并使用两列文件自行翻译。

比较来自三个不同基因定义系统的程序消息（“为yyy独特基因完成xxx转录”），我们可以看到Ensembl注释了人类基因组中最多的基因，而RefSeq注释了最少的基因。

技术说明：从技术上讲，RefSeq基因和UCSC基因是transcript-based gene definitions。他们基于转录本数据建立了基因模型，然后将该基因模型映射回人类基因组。相比之下，Ensemble基因和Gencode基因是assembly-based gene definitions，试图直接从参考人类基因组构建基因模型。他们来自不同的角度，试图做同样的事情：在人类基因组中定义基因。

但是，这有其自身的后果。例如，RefSeq通过组装来自种群的转录本数据来构建基因模型，但是参考人类基因组可能具有等位基因作为种群中的次要等位基因。在这种情况下，转录本可能不会与基因组100％比对，从而导致转录本的FASTA文件与从全基因组序列生成的FASTA文件之间存在差异（通过将外显子连接在一起）。

由于这些原因，外显子变异的准确注释不能依赖于公共数据库中的cDNA序列，而只能基于基因组自身实际的chr：start-end位置。因此，我为一些特定的基因组构建了FASTA序列，用户可以直接从ANNOVAR网站下载序列。我还提供了程序（retrieve_seq_from_fasta.pl）为我不提供预构建文件的任何其他基因组构建FASTA序列。

由于这些原因，我提供的文件中的FASTA序列可能与您从RefSeq获得的FASTA序列不同。ANNOVAR使用的序列是基于特定基因组构建和组装的“理论”序列，但是RefSeq编译的FASTA序列是大型数据库中的“观察到”序列，与特定组装版本没有任何关系。它们可能具有相同的标识符，但它们是不同的东西。

切换到hg38 Ensembl基因注释

请注意，UCSC Genome Browser无法提供hg38的ensGene。一位用户指出，UCSC已从22版开始使用GENCODE替换了ensGene.txt（请参阅https://groups.google.com/a/soe.ucsc.edu/forum/#!topic/genome/uOROZuefx_Y）。因此，如果要基于hg38注释Ensemble基因，则应改用Gencode文件。

2017年9月，应用户要求，我现在使用26 GENCODE Basic版本直接在ANNOVAR中为hg38准备了ensGene。下面显示了我用于构建此文件的命令供您参考：

annotate_variation.pl -downdb wgEncodeGencodeBasicV26 tempdir -build hg38
retrieve_seq_from_fasta.pl -format genericGene -seqfile humandb/hg38_seq/hg38.fa -outfile tempdir/hg38_wgEncodeGencodeBasicV26Mrna.fa tempdir/hg38_wgEncodeGencodeBasicV26.txt
mv hg38_wgEncodeGencodeBasicV26.txt hg38_ensGene.txt
mv hg38_wgEncodeGencodeBasicV26Mrna.fa hg38_ensGeneMrna.fa

现在，我只是将这些文件上传到ANNOVAR数据库存储库，以便用户可以直接-build hg38 -downdb ensGene dir/下载ensGene的hg38版本并使用Ensemble关键字执行基于基因的注释。

切换到GENCODE / CCDS基因注释

当前版本的ANNOVAR没有为GENCODE提供特定的关键字，但是ANNOVAR具有足够的通用性，可以处理GENCODE或其他任何基因定义就可以了。GENCODE系统在过去几年中发生了很大变化，最新版本为2014年2月的V19。在V19中，BASIC轨道现已可用，该轨道包含基于以下描述的高质量基因定义：“ GENCODE Basic集合旨在提供对大多数用户有用的GENCODE成绩单注释的简化子集，目标是拥有一个涵盖所有基因座的高质量基本集合，选择GENCODE注释以包含在基本集合中在每个基因位点分别确定编码和非编码转录本。”

以下是我可以用于GENCODE基因定义以进行变体注释的命令（请注意，如果，全基因组FASTA文件不可用humandb/hg19_seq，则应首先执行annotate_variation.pl -downdb -build hg19 seq humandb/hg19_seq/）。

[kaiwang@biocluster ~/]$ annotate_variation.pl -downdb wgEncodeGencodeBasicV19 humandb/ -build hg19
NOTICE: Web-based checking to see whether ANNOVAR new version is available ... Done
NOTICE: Downloading annotation database http://hgdownload.cse.ucsc.edu/goldenPath/hg19/database/wgEncodeGencodeBasicV19.txt.gz ... OK
NOTICE: Uncompressing downloaded files
NOTICE: Finished downloading annotation files for hg19 build version, with files saved at the 'humandb' directory

[kaiwang@biocluster ~/]$ retrieve_seq_from_fasta.pl -format genericGene -seqdir humandb/hg19_seq/ -outfile humandb/hg19_wgEncodeGencodeBasicV19Mrna.fa humandb/hg19_wgEncodeGencodeBasicV19.txt 

[kaiwang@biocluster ~/]$ annotate_variation.pl -build hg19 -out ex1 -dbtype wgEncodeGencodeBasicV19 example/ex1.avinput humandb/ 
NOTICE: The --geneanno operation is set to ON by default
NOTICE: Reading gene annotation from humandb/hg19_wgEncodeGencodeBasicV19.txt ... Done with 95929 transcripts (including 38291 without coding sequence annotation) for 42594 unique genes
NOTICE: Reading FASTA sequences from humandb/hg19_wgEncodeGencodeBasicV19Mrna.fa ... Done with 16 sequences
WARNING: A total of 303 sequences will be ignored due to lack of correct ORF annotation
NOTICE: Finished gene-based annotation on 15 genetic variants in example/ex1.avinput
NOTICE: Output files were written to ex1.variant_function, ex1.exonic_variant_function
Then check the results below:

[kaiwang@biocluster ~/]$ cat ex1.variant_function
UTR5 ISG15(ENST00000379389.4:c.-33T>C) 1 948921 948921 T C comments: rs15842, a SNP in 5' UTR of ISG15
UTR3 ATAD3C(ENST00000378785.2:c.*91G>T) 1 1404001 1404001 G T comments: rs149123833, a SNP in 3' UTR of ATAD3C
splicing NPHP4(ENST00000378156.4:exon22:c.2818-2T>A) 1 5935162 5935162 A T comments: rs1287637, a splice site variant in NPHP4
intronic DDR2 1 162736463 162736463 C T comments: rs1000050, a SNP in Illumina SNP arrays
intronic DNASE2B 1 84875173 84875173 C T comments: rs6576700 or SNP_A-1780419, a SNP in Affymetrix SNP arrays
intergenic RP13-221M14.3(dist=46825),PRAMEF26(dist=5062) 1 13211293 13211294 TC - comments: rs59770105, a 2-bp deletion
intergenic UBIAD1(dist=48375),PTCHD2(dist=135699) 1 11403596 11403596 - AT comments: rs35561142, a 2-bp insertion
intergenic RP11-364B6.1(dist=872522),RP11-251P6.1(dist=652543) 1 105492231 105492231 A ATAAA comments: rs10552169, a block substitution
exonic IL23R 1 67705958 67705958 G A comments: rs11209026 (R381Q), a SNP in IL23R associated with Crohn's disease
exonic ATG16L1 2 234183368 234183368 A G comments: rs2241880 (T300A), a SNP in the ATG16L1 associated with Crohn's disease
exonic NOD2 16 50745926 50745926 C T comments: rs2066844 (R702W), a non-synonymous SNP in NOD2
exonic NOD2 16 50756540 50756540 G C comments: rs2066845 (G908R), a non-synonymous SNP in NOD2
exonic NOD2 16 50763778 50763778 - C comments: rs2066847 (c.3016_3017insC), a frameshift SNP in NOD2
exonic GJB2 13 20763686 20763686 G - comments: rs1801002 (del35G), a frameshift mutation in GJB2, associated with hearing loss
exonic CRYL1,GJB6 13 20797176 21105944 0 - comments: a 342kb deletion encompassing GJB6, associated with hearing loss

[kaiwang@biocluster ~/]$ cat ex1.exonic_variant_function
line9 nonsynonymous SNV IL23R:ENST00000395227.1:exon4:c.G377A:p.R126Q,IL23R:ENST00000347310.5:exon9:c.G1142A:p.R381Q, 1 67705958 6770595comments: rs11209026 (R381Q), a SNP in IL23R associated with Crohn's disease
line10 nonsynonymous SNV ATG16L1:ENST00000373525.5:exon6:c.A466G:p.T156A,ATG16L1:ENST00000347464.5:exon5:c.A409G:p.T137A,ATG16L1:ENST00000392018.1:exon10:c.A949G:p.T317A,ATG16L1:ENST00000392017.4:exon9:c.A898G:p.T300A,ATG16L1:ENST00000392020.4:exon8:c.A841G:p.T281A, 2 234183368 234183368 comments: rs2241880 (T300A), a SNP in the ATG16L1 associated with Crohn's disease
line11 nonsynonymous SNV NOD2:ENST00000300589.2:exon4:c.C2104T:p.R702W, 16 50745926 50745926 C T comments: rs2066844 (R702W), a non-synonymous SNP in NOD2
line12 nonsynonymous SNV NOD2:ENST00000300589.2:exon8:c.G2722C:p.G908R, 16 50756540 50756540 G C comments: rs2066845 (G908R), a non-synonymous SNP in NOD2
line13 frameshift insertion NOD2:ENST00000300589.2:exon11:c.3017dupC:p.A1006fs, 16 50763778 50763778 - C comments: rs2066847 (c.3016_3017insC), a frameshift SNP in NOD2
line14 frameshift deletion GJB2:ENST00000382848.4:exon2:c.35delG:p.G12fs,GJB2:ENST00000382844.1:exon1:c.35delG:p.G12fs, 13 20763686 2076368comments: rs1801002 (del35G), a frameshift mutation in GJB2, associated with hearing loss
line15 frameshift deletion GJB6:ENST00000241124.6:wholegene,CRYL1:ENST00000298248.7:wholegene,GJB6:ENST00000400065.3:wholegene,CRYL1:ENST00000382812.1:wholegene,GJB6:ENST00000356192.6:wholegene,GJB6:ENST00000400066.3:wholegene, 13 20797176 21105944 0 - comments: a 342kb deletion encompassing GJB6, associated with hearing loss
Similarly, you can switch to Comprehensive annotations and PolyA annotations of GENCODE, by using a different -dbtype argument (wgEncodeGencodeCompV19 and wgEncodeGencodePolyaV19).

使用这些注释时，务必始终注意要使用的正确文件名（如果不确定表名，请使用UCSC Genome Browser中的表浏览器），这一点很重要。

还支持许多其他基因定义系统。例如，用户可以选择CCDS基因定义。（请注意，如果，全基因组FASTA文件不可用humandb/hg19_seq，则应首先执行annotate_variation.pl -downdb -build hg19 seq humandb/hg19_seq/）。

[kaiwang@biocluster ~/]$ annotate_variation.pl -downdb -build hg19 ccdsGene humandb

[kaiwang@biocluster ~/]$ retrieve_seq_from_fasta.pl humandb/hg19_ccdsGene.txt -seqdir humandb/hg19_seq -format refGene -outfile humandb/hg19_ccdsGeneMrna.fa

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -out ex1 -dbtype ccdsGene example/ex1.avinput humandb/
NOTICE: The --geneanno operation is set to ON by default
NOTICE: Reading gene annotation from humandb/hg19_ccdsGene.txt ... Done with 29045 transcripts (including 0 without coding sequence annotation) for 29014 unique genes
NOTICE: Reading FASTA sequences from humandb/hg19_ccdsGeneMrna.fa ... Done with 8 sequences
WARNING: A total of 38 sequences will be ignored due to lack of correct ORF annotation
NOTICE: Finished gene-based annotation on 15 genetic variants in example/ex1.avinput
NOTICE: Output files were written to ex1.variant_function, ex1.exonic_variant_function

技术说明：对于CCDS基因，输出将不包含基因名称，而仅包含CCDS标识符。为了获得基因名称，用户必须编写自己的程序来处理ANNOVAR输出文件。有两个表可用于将CCDS ID转换为其他ID：ccdsNotes将CCDS转换为UCSC已知基因笔录ID，然后可以将其转换为Gene name。ccdsInfo表将CCDS ID转换为ENSEMBL成绩单或RefSeq成绩单，然后您可以将其进一步转换为Gene name。

为非人类物种创建自己的基因定义数据库

除人类基因组外，其他物种也可以处理。但是，ANNOVAR没有为其他基因定义提供内置的mRNA FASTA文件，因此用户必须自己构建它。

要进一步了解这一点，请尝试处理黑猩猩的基因组：

[kai@beta ~/]$ annotate_variation.pl -downdb -buildver panTro2 gene chimpdb
NOTICE: Downloading annotation database ftp://hgdownload.cse.ucsc.edu/goldenPath/panTro2/database/refGene.txt.gz ... OK 
NOTICE: Downloading annotation database ftp://hgdownload.cse.ucsc.edu/goldenPath/panTro2/database/refLink.txt.gz ... OK 
NOTICE: Downloading annotation database http://www.openbioinformatics.org/annovar/download/panTro2_refGeneMrna.fa.gz ... Failed
NOTICE: Uncompressing downloaded files
NOTICE: Finished downloading annotation files for panTro2 build version, with files saved at the 'chimpdb' directory
WARNING: Some files cannot be downloaded, including http://www.openbioinformatics.org/annovar/download/panTro2_refGeneMrna.fa.gz
--------------------------------IMPORTANT---------------------------------
--------------------------------------------------------------------------
NOTICE: the FASTA file http://www.openbioinformatics.org/annovar/download/panTro2_refGeneMrna.fa.gz is not available to download but can be generated by the ANNOVAR software. PLEASE RUN THE FOLLOWING TWO COMMANDS CONSECUTIVELY TO GENERATE THE FASTA FILES:
annotate_variation.pl --buildver panTro2 --downdb seq chimpdb/panTro2_seq
retrieve_seq_from_fasta.pl chimpdb/panTro2_refGene.txt -seqdir chimpdb/panTro2_seq -format refGene -outfile chimpdb/panTro2_refGeneMrna.fa
--------------------------------------------------------------------------
--------------------------------------------------------------------------

The above command will run, but will print out some warning message: the FASTA sequences are not provided in ANNOVAR website so users need to build them. Just follow the exact intructions and run the two commands:

[kai@beta ~/]$ annotate_variation.pl --buildver panTro2 --downdb seq chimpdb/panTro2_seq
NOTICE: Downloading annotation database ftp://hgdownload.cse.ucsc.edu/goldenPath/panTro2/bigZips/chromFa.zip ... Failed
NOTICE: Downloading annotation database ftp://hgdownload.cse.ucsc.edu/goldenPath/panTro2/bigZips/chromFa.tar.gz ... OK 
NOTICE: Uncompressing downloaded files
NOTICE: Finished downloading annotation files for panTro2 build version, with files saved at the 'chimpdb/panTro2_seq' directory

[kai@beta ~/]$ retrieve_seq_from_fasta.pl chimpdb/panTro2_refGene.txt -seqdir chimpdb/panTro2_seq -format refGene -outfile chimpdb/panTro2_refGeneMrna.fa
NOTICE: Finished reading 1 sequences from chimpdb/panTro2_seq/12/chr12_random.fa
NOTICE: Finished reading 1 sequences from chimpdb/panTro2_seq/22/chr22.fa
NOTICE: Finished reading 1 sequences from chimpdb/panTro2_seq/14/chr14.fa
......
......
NOTICE: Finished writting FASTA for 1337 genomic regions to chimpdb/panTro2_refGeneMrna.fa.

因此，在运行上述命令后，黑猩猩基因组的基因注释数据库将是完整，准确和最新的。

练习：尝试对rheMac2（猕猴）运行上述相同的过程，并查看其与panTro2有何不同。UCSC没有为不同的基因组使用相同的文件命名约定或目录结构规则，这使程序员的生活变得更加复杂。ANNOVAR可以处理许多基因组，但是将存在另一个ANNOVAR无法自动检索序列的基因组。如果是这种情况，请报告给我，我将进行调查并添加功能。

练习：尝试对sacCer2（酵母）运行与上面相同的过程，并观察它们之间的不同。

练习：尝试对bosTau6（cow）运行与上面相同的过程。请注意，截至2012年4月，UCSC尚未将bosTau6基因组序列的FASTA文件拆分为单个染色体。因此，用户需要在retrieve_seq_from_fasta.pl命令中使用“ -seqfile bosTau6.fa”，而不是“ -seqdir cowdb / bosTau6_seq”。同样，尝试对micMur1（鼠标狐猴）运行相同的过程，并注意使用-seqfile而不是-seqdir。

练习：尝试对rn5（大鼠）或dm6（果蝇）运行相同的步骤。同样，用户需要提供FASTA文件而不是FASTA目录。

仅当特定物种或特定构建体的UCSC中存在基于基因的注释时，以上过程才有效。例如，如果您想在猪上使用ANNOVAR，由于RefSeq基因和UCSC基因不适用于猪，因此您必须使用RefSeq基因和UCSC基因，annotate_variation.pl --downdb -buildver susScr2 ensgene pigdb并将其-dbtype ensgene用于基于基因的注释。

那么新物种的GFF3文件呢？

有时，用户会自己对新物种进行测序，因此无法从UCSC或Ensembl或任何地方获得基因组构建。一般而言，您可以使用几种基因预测工具之一从您构建的重叠群中编译基因结构，并且通常采用GFF3格式。从Ensembl下载或由用户编译的GFF3或GTF文件需要转换为GenePred格式。可以通过http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64/上的gff3ToGenePred或gtfToGenePred工具执行转换。根据我们的经验，有时Ensembl的某些GFF3文件无法正确转换。因此，对于无法从UCSC注释数据库中获得的非人类物种，我们建议使用GTF文件生成基因定义文件。

在本文中，我们举了一个示例来对拟南芥进行注释，其过程复制如下：

a. 访问http://plants.ensembl.org/info/website/ftp/index.html将该植物的GTF文件和基因组FASTA文件下载到名为atdb的文件夹中。

mkdir atdb
cd atdb
wget ftp://ftp.ensemblgenomes.org/pub/release-27/plants/fasta/arabidopsis_thaliana/dna/Arabidopsis_thaliana.TAIR10.27.dna.genome.fa.gz
wget ftp://ftp.ensemblgenomes.org/pub/release-27/plants/gtf/arabidopsis_thaliana/Arabidopsis_thaliana.TAIR10.27.gtf.gz

b. 请解压缩两个文件：

gunzip Arabidopsis_thaliana.TAIR10.27.dna.genome.fa.gz 
gunzip Arabidopsis_thaliana.TAIR10.27.gtf.gz

c. 请使用gtfToGenePred工具将GTF文件转换为GenePred文件：

gtfToGenePred -genePredExt Arabidopsis_thaliana.TAIR10.27.gtf AT_refGene.txt`

d. 请使用我们提供的脚本生成转录本FASTA文件：

perl retrieve_seq_from_fasta.pl --format refGene --seqfile Arabidopsis_thaliana.TAIR10.27.dna.genome.fa AT_refGene.txt --out AT_refGeneMrna.fa

完成此步骤后，就准备好了基于基因的注释所需的注释数据库文件。现在，您可以注释给定的VCF文件。请注意，--buildver参数应设置为AT。

最近，一个用户分享了他的经验来注释大麦基因组。通过下面的大麦示例，在ENSEMBL PLANTS中发布的组装序列和注释中使用ANNOVAR建立数据库的程序

a. 在ANNOVAR目录中，通过键入以下内容来创建数据库目录：

mkdir barleydb

b. 转到barleydb并从ensembl plants网站下载两个文件（一个组装好的序列和一个注释）

cd barleydb
wget ftp://ftp.ensemblgenomes.org/pub/plants/release-39/fasta/hordeum_vulgare/dna/Hordeum_vulgare.Hv_IBSC_PGSB_v2.dna.toplevel.fa.gz
wget ftp://ftp.ensemblgenomes.org/pub/plants/release-39/gff3/hordeum_vulgare/Hordeum_vulgare.Hv_IBSC_PGSB_v2.39.gff3.gz

c。安装gff3ToGenePred。请注意，这是一个二进制文件。

请注意，假设您在服务器中安装了anaconda，可以采用以下方法： conda install -c bioconda ucsc-gff3togenepred

请注意，另一种方法是使用Google进行搜索，然后从网络上获取它，并将其放在ANNOVAR目录中： hgdownload.cse.ucsc.edu/admin/exe/linux.x86_64/gff3ToGenePred

d。解压缩两个下载的大麦H * .gz文件并将gff3转换为GenePred格式

gunzip Hordeum*.gz
./gff3ToGenePred Hordeum_vulgare.Hv_IBSC_PGSB_v2.39.gff3 HV39-refGene0.txt

e. 使用nl在HV39_refGene0.txt中添加一个随机的第一个字段（将行号添加到每行）

nl HV39_refGene0.txt >HV39_refGene.txt

f. 返回到ANNOVAR目录并在下面执行

cd ..
retrieve_seq_from_fasta.pl barleydb/HV39_refGene.txt -seqfile barleydb/Hordeum_vulgare.Hv_IBSC_PGSB_v2.dna.toplevel.fa -format ensGene -outfile barleydb/HV39_refGeneMrna.fa

如果ANNOVAR目录未放在.bash_profile或.bashsrc的PATH中，则必须在命令前面使用perl：

perl retrieve_seq_from_fasta.pl barleydb/HV39_refGene.txt -seqfile barleydb/Hordeum_vulgare.Hv_IBSC_PGSB_v2.dna.toplevel.fa -format ensGene -outfile barleydb/HV39_refGeneMrna.fa

了解ANNOVAR如何解决基因定义中的罕见问题

有时，refGene或knownGene注释本身包含错误。下面是一个清晰的示例：

gene_shot

从上图可以看出，refGene对终止密码子进行了错误注释（以致阅读框中缺少一个核苷酸T）。但是，UCSC基因正确地识别了最后一个终止密码子。ANNOVAR会识别此类错误，对于refGene注释，它会将位于hr17：3141680的变异注释为“ UNKNOWN”作为外显子功能（因为密码子不能只有2个核苷酸）。

对于某些变异，ANNOVAR可能会与其他竞争性注释软件或流程相比生成不同的注释。这可能是由于多种原因，例如使用不同的基因定义。尽管如此，在Genome Browser中人工检查变异并确认ANNOVAR在注释中是否出错始终是一个好主意。不幸的是，这仅适用于UCSC基因（请参见上面的示例），但是对于大多数基因而言，UCSC基因与refGene注释完全一致。

下面说明一个有趣的示例。另一个软件将“ 3 17028503 17028503 A G”注释为同义突变synonymous，但ANNOVAR通过refGene注释将其注释为“非同义突变”non-synonymous。事实证明，refGene在此区域提供了两个转录本注释，并且相同的突变可以是同义也可以是非同义。在ANNOVAR批注中，非同义突变会覆盖同义突变（同样，请阅读上面的优先级表以了解这一点），因此生成的是非同义突变。如果用户在命令行中使用参数-separate，则ANNOVAR将在输出文件中打印两个注释。

gene_shot

下面给出了另一个有趣的示例。使用RefSeq注释，ANNOVAR将突变“ chr12 6945846 6945846 A C”注释为stop-lost。在基因组浏览器中对其进行检查表明，UCSC基因注释了该区域中的多个转录本，因此，根据用户想要使用的实际基因定义，该突变可以是stop-lost，3-UTR或内含子的。这是一种极为罕见的情况，但是用户在解释数据时应牢记这一点，尤其是在找到潜在的候选变异之后。与其他注释系统相比，ANNOVAR在这方面表现出色，因为它可以灵活选择基因定义系统。

gene_shot

另一个有趣的示例如下所示。两个基因（DGCR14和TSSK2）彼此重叠，但是一个基因的编码区是另一个基因的UTR。对于其中一个基因，Ensembl提供了三个肽段/转录物注释，但是其中一个（ENST00000383058）在方向上与另一个完全相反，带有一个额外的氨基酸。因此，突变本身无法确定转录本的方向（由基因注释系统和基因组构建决定），这就是为什么ANNOVAR不提供方向注释的原因。要了解方向，用户将必须执行grep ENST00000383058 hg18_ensGene.txt，其中包括有关基因注释和构建注释的信息。

gene_shot

有用的教程

2020年夏天，我教了一个在线培训课程，该课程为绝对的初学者提供了一些资料，包括使用个人Windows和Mac膝上型计算机而不是Linux服务器的初学者。本教程要求您在个人计算机上安装conda，并且可以在此处进行访问。一些用户可能会发现它很有用。但是，如果您已经在使用计算集群，并且已经熟悉Linux，则无需遵循本教程，而只需阅读以下内容。

table_annovar.pl

对于初学者，使用ANNOVAR的最简单方法是使用该table_annovar.pl程序。该程序采用一个输入变异文件（例如VCF文件），并生成一个由制表符分隔的输出文件，其中包含许多列，每列代表一组注释。此外，如果输入是VCF文件，则程序还将生成一个新的输出VCF文件，其INFO字段中会填充注释信息。

假设我们已经下载了ANNOVAR软件包并用于tar xvfz annovar.latest.tar.gz解压缩该软件包。您将看到该bin/目录包含几个带有.pl后缀的Perl程序。（请注意，如果您已经将ANNOVAR路径添加到系统可执行文件路径中，则可以输入annotate_variation.pl而不是输入perl annotate_variation.pl）。首先，我们需要使用来下载适当的数据库文件annotate_variation.pl，然后我们将运行table_annovar.pl程序来注释example/ex1.avinput文件中的变体。

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -downdb -webfrom annovar refGene humandb/

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -downdb cytoBand humandb/

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -downdb -webfrom annovar exac03 humandb/ 

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -downdb -webfrom annovar avsnp147 humandb/ 

[kaiwang@biocluster ~/]$ annotate_variation.pl -buildver hg19 -downdb -webfrom annovar dbnsfp30a humandb/

[kaiwang@biocluster ~/]$ table_annovar.pl example/ex1.avinput humandb/ -buildver hg19 -out myanno -remove -protocol refGene,cytoBand,exac03,avsnp147,dbnsfp30a -operation gx,r,f,f,f -nastring . -csvout -polish -xref example/gene_xref.txt

逐一运行以上命令。前几个命令将适当的数据库下载到humandb/目录中。最终命令运行TABLE_ANNOVAR，使用ExAC版本0.3（称为exac03），dbNFSP版本3.0a（称为dbnsfp30a），左侧标准化的dbSNP版本147（称为avsnp147），并删除所有临时文件，并生成输出文件myanno.hg19_multianno.txt。没有任何注释的字段将由“.”填充。在Excel中打开输出文件，然后查看其中包含的内容。我生成输出文件可以在这里下载：ex1.hg19_multianno.csv。前几列的屏幕截图如下所示：

table_annovar

输出文件包含多个列。前几列是您的输入列。以下每列对应于用户在命令行中指定的“协议”（protocol）之一。Func.refGene，Gene.refGene，GeneDetail.refGene，ExonicFunc.refGene，AAChange.refGene列包含有关突变如何影响基因结构的各种注释。Xref.refGene列包含该基因的交叉引用（cross-reference ）。在这种情况下，已知的遗传疾病是否是由该基因的缺陷引起的（此信息在命令行中的example/gene_xref.txt文件）。在接下来的几列中，ExAC 列表示所有样本中的等位基因频率，以及Exome Aggregation Consortium数据集中的sub-populations，而avsnp147表示dbSNP版本147中的SNP标识符。其他列包含几种广泛使用的工具获得的非同义突变的预测分数*，包括SIFT分数，PolyPhen2 HDIV分数，PolyPhen2 HVAR分数，LRT分数，MutationTaster分数，MutationAssessor分数，FATHMM分数，GERP ++分数，CADD分数，DANN分数，PhyloP分数和SiPhy分数等等。

我们可以更详细地检查命令行。-operation参数告诉ANNOVAR每种协议的操作：

g表示基于基因的，
gx表示（gene-based with cross-reference annotation）带有交叉引用注释的基于基因的（来自-xref变量），
r表示基于区域和f表示基于过滤器。

如果不提供xref文件，则该操作只能是g。您将在其他网页中找到有关基于基因/区域/过滤器的注释的详细信息。有时，用户需要制表符分隔的文件，而不是逗号分隔的文件。删除上面命令的-csvout参数就可以了。

在上面的命令中，我们使用-xreffile参数为基因提供注释。如果文件包含标题行，则可以为基因提供多条注释（而不只是一列）。为了说明这一点，我们可以检查example/gene_fullxref.txt文件的前两行（包括标题行）：

[kaiwang@biocluster ~/project/annotate_variation]$ head -n 2 example/gene_fullxref.txt
#Gene_name      pLi     pRec    pNull   Gene_full_name  Function_description    Disease_description     Tissue_specificity(Uniprot)     Expression(egenetics)  Expression(GNF/Atlas)    P(HI)   P(rec)  RVIS    RVIS_percentile GDI     GDI-Phred
A1BG    9.0649236354772e-05     0.786086131023045       0.2138232197406 alpha-1-B glycoprotein  .       .       TISSUE SPECIFICITY: Plasma.;    unclassifiable (Anatomical System);amygdala;prostate;lung;islets of Langerhans;liver;spleen;germinal center;brain;thymus;       fetal liver;liver;fetal lung;trigeminal ganglion;       0.07384 0.31615 -0.466531444    23.51380042     79.3774 1.88274

标题行以开头#。然后，交叉引用文件（ cross-reference file ）包含15种类型的基因注释。您可以运行上面相同的命令，但将-xreffile从gene_xref.txt变为gene_fullxref.txt，结果文件可以从以下网址下载这里。下面显示了文件的一部分，为用户提供了示例：

table_annovar_fullxref

table_annovar.pl可以直接支持VCF文件的输入和输出（注释将被写入输出VCF文件的INFO字段）。

[kaiwang@biocluster ~/]$ table_annovar.pl example/ex2.vcf humandb/ -buildver hg19 -out myanno -remove -protocol refGene,cytoBand,exac03,avsnp147,dbnsfp30a -operation g,r,f,f,f -nastring . -vcfinput -polish

您可以在此处下载输出文件：ex2.hg19_multianno.vcf。此外，制表符分隔的输出文件也可以作为ex2.hg19_multianno.txt获得，该文件包含不同格式的相似信息。您可以在文本编辑器中打开新的VCF文件，并检查文件中所做的更改：VCF文件中的INFO字段现在包含所需的annotations ，以字符串ANNOVAR_DATE开头，以符号ALLELE_END结尾。如果在同一基因座中有多个等位基因，您将在INFO字段中看到多个这样的符号。屏幕截图如下所示：

table_vcf

希望在完成上述练习之后，您现在可以更好地了解什么是ANNOVAR，并且可以开始享受对变异进行注释的过程。

annotate_variation.pl

annotate_variation.pl程序是ANNOVAR的核心程序。我们可以转到ANNOVAR目录，然后一个接一个地运行以下三个命令。

annotate_variation.pl -geneanno -dbtype refGene -buildver hg19 example/ex1.avinput humandb/

annotate_variation.pl -regionanno -dbtype cytoBand -buildver hg19 example/ex1.avinput humandb/ 

annotate_variation.pl -filter -dbtype exac03 -buildver hg19 example/ex1.avinput humandb/

请注意，这三个命令对应于基于基因，基于区域和基于过滤器的注释。

第一个命令注释ex1.avinput文件中的12个变异，并将它们分类为基因间，内含子，非同义SNP，移码删除，大规模复制等。ex1.avinput文件每行一个变异。注释会生成ex1.avinput.variant_function和ex1.avinput.exonic_variant_function两个输出文件。检查example/目录中的两个输出文件以查看它们包含的内容：在variant_function文件中，第一列和第二列注释对基因结构和受变异影响的基因，而其他列从输入文件复制。在exonic_variant_function 文件中，第一，第二和第三列注释输入文件中的变异行号，变异对编码序列的影响以及受影响的基因/转录本，而其他列则从输入文件中复制。

接下来，程序在ex1.avinput文件中注释变异并识别这些变异的细胞遗传学带（cytogenetic band）。注释过程应花费几秒钟。检查输出文件ex1.avinput.hg19_cytoBand以查看其包含的内容。第一列显示cytoBand，第二列显示注释结果，其他列从输入文件复制。

接下来，程序将识别出一个变异的子集，这些变异在ex1.avinputexac03数据库中未观察到（保存在ex1.avinput.hg19_exac03_filtered中）和在等位基因频率下观察到的变异（保存在ex1.avinput.hg19_exac03_dropped文件中）。

技术说明：默认情况下，ANNOVAR在hg18（人类基因组NCBI build 36）坐标上标注变异。由于输入文件位于hg19坐标中，因此我们在上面的每个命令中都添加了-buildver hg19。同样，如果您是从人类GRCh38坐标生成variant calls的，则每个命令添加-buildver hg38；如果变异文件来自fly，则您使用的每个命令添加-buildver dm3；如果您的变异文件来自鼠，请您使用的每个命令添加-buildver mm9……

上面的命令代表了一组有关ANNOVAR如何帮助研究人员注释从高通量测序数据产生的遗传变异的基本示例。如果您有兴趣，请单击左侧的横幅以了解输入格式，并查看ANNOVAR可以为您的研究做些什么...

额外资源

以下参考资料提供了有关使用ANNOVAR和wANNOVAR的分步协议（step-by-step protocol）

Yang H，Wang K. 使用ANNOVAR和wANNOVAR Nature Protocols进行基因组变异注释和优先级排序，10：1556-1566，2015

VCF文件

table_annovar.pl程序可以获取VCF文件并对其进行注释（-vcfinput参数）。如今，VCF已经成为大多数研究人员使用的黄金标准格式。有关处理VCF文件的其他建议，请参阅本文的“ VCF处理指南 ”。

ANNOVAR输入文件

annotate_variation.pl程序需要一种简单的基于文本的格式，我们将其称为ANNOVAR输入格式。在此文件中，每一行对应一个变异。在每行上，前五个以空格或制表符分隔的列代表chromosome, start position, end position, the reference nucleotides and the observed nucleotides。可以提供其他列，并以相同的形式打印出来。为方便起见，如果此信息不易获得，则用户可以使用“ 0”填写参考核苷酸。通过使用“ –”代表无效核苷酸，可以通过这种简单的文件格式轻松表示Insertions, deletions 或 block substitutions。下面一个示例（此示例包含在ex1.avinput文件）作为有关变异的注释，并带有额外的列。默认情况下，使用基于1的坐标系。

ANNOVAR软件包包含一些示例输入文件。例如，ex1.avinput文件的内容如下：

[kaiwang@biocluster ~/]$ cat example/ex1.avinput
1 948921 948921 T C comments: rs15842, a SNP in 5' UTR of ISG15
1 1404001 1404001 G T comments: rs149123833, a SNP in 3' UTR of ATAD3C
1 5935162 5935162 A T comments: rs1287637, a splice site variant in NPHP4
1 162736463 162736463 C T comments: rs1000050, a SNP in Illumina SNP arrays
1 84875173 84875173 C T comments: rs6576700 or SNP_A-1780419, a SNP in Affymetrix SNP arrays
1 13211293 13211294 TC - comments: rs59770105, a 2-bp deletion
1 11403596 11403596 - AT comments: rs35561142, a 2-bp insertion
1 105492231 105492231 A ATAAA comments: rs10552169, a block substitution
1 67705958 67705958 G A comments: rs11209026 (R381Q), a SNP in IL23R associated with Crohn's disease
2 234183368 234183368 A G comments: rs2241880 (T300A), a SNP in the ATG16L1 associated with Crohn's disease
16 50745926 50745926 C T comments: rs2066844 (R702W), a non-synonymous SNP in NOD2
16 50756540 50756540 G C comments: rs2066845 (G908R), a non-synonymous SNP in NOD2
16 50763778 50763778 - C comments: rs2066847 (c.3016_3017insC), a frameshift SNP in NOD2
13 20763686 20763686 G - comments: rs1801002 (del35G), a frameshift mutation in GJB2, associated with hearing loss
13 20797176 21105944 0 - comments: a 342kb deletion encompassing GJB6, associated with hearing loss

前五个空格或制表符分隔的字段是染色体（“ chr”前缀是可选的），“开始”，“结束”，“参考等位基因”，“替代等位基因”。其余的列是完全可选的。

上面的例子包含几个遗传变异。

第一个变异是单核苷酸变异，参考基因组中的C被T取代。
第三个变异是2 bp缺失，观察到的核苷酸由“-”表示。
第四个变异是2 bp的插入，因为参考基因组中的参考核苷酸由“-”表示。
最后一个变异是大规模缺失，但是参考等位基因由“ 0”表示，从而无需在此行上包括参考核苷酸。

另一个示例如下所示。请注意，前五列符合上述规范，而其他所有列完全是可选的，用户可以在此处放置任何内容。

[kaiwang@biocluster ~/]$ cat ex3.avinput
7 92570705 92570705 T C 7 3 43 D G SAMD9 1.56
7 98870495 98870495 G A 26 16 62 R C PTCD1 3.06
7 99835402 99835402 C T 13 6 46 P L PILRA 1.75
7 100122289 100122289 - CCT 5 3 60 EQ ERQ GIGYF1 3.98
7 100209410 100209410 G A 15 8 53 R H ZAN 1.81
7 100473466 100473466 A G 38 13 34 T A MUC17 0.60
7 105066159 105066160 TC - 19 3 16 E X ATXN7L1 4.92
7 113306419 113306419 C T 15 6 40 S N PPP1R3A 1.05
7 115411632 115411632 C T 14 5 36 D N TFEC -0.45
7 119702880 119702880 T C 29 3 10 C R KCND2 5.00
7 120216091 120216091 C G 20 10 50 A P TSPAN12 5.32
7 120555712 120555712 T C 10 3 30 L S C7orf58 3.00
7 128099699 128099699 C G 7 5 71 I M FAM71F2 0.42
7 128221650 128221650 T C 13 3 23 L P CCDC136 1.26

在某些情况下，用户可能只想指定位置，而不要指定实际核苷酸。在这种情况下，可以使用“ 0”填写第4和第5列。ANNOVAR仍可以在此输入文件上运行，但是显然氨基酸变化没有输出。另外，假定观察到的氨基酸与野生型等位基因的长度相等（由每行的开始和结束位置指定）。

如果ANNOVAR遇到无效的输入行，它将把无效的行写入到--outfile参数指定的名为<outfile>.invalid_inputwhere 的<outfile>文件中。如果所有输入行均为有效格式，则此输出文件将不存在。因此，即使输入文件包含空行或无效格式，ANNOVAR仍可以继续下一个输入行。

下载包包含几个示例输入文件。用户可以签出查看它们。

格式转换

该convert2annovar.pl脚本可以将其他“genotype calling”格式转换为ANNOVAR格式。目前，该程序可以处理Samtools genotype-calling pileup格式，Illumina CASAVA格式，SOLiD GFF genotype-calling格式，Complete Genomics variant 格式，SOAPsnp格式，MAQ格式和VCF格式。此外，该程序还可以从dbSNP标识符列表，转录本标识符或基因组区域生成ANNOVAR输入文件。应该强调的是，当今VCF是使用最广泛的文件格式，并且不再使用其他大多数格式。

-VCF4格式

该-format vcf4参数可用于将VCF文件转换为ANNOVAR输入格式。ANNOVAR软件包应在example /目录中包含一个示例VCF文件。我们可以使用此文件作为示例：

[kaiwang@biocluster ~/]$ cat example/ex2.vcf 
##fileformat=VCFv4.0
##fileDate=20090805
##source=myImputationProgramV3.1
##reference=1000GenomesPilot-NCBI36
##phasing=partial
##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">
##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">
##FILTER=<ID=q10,Description="Quality below 10">
##FILTER=<ID=s50,Description="Less than 50% of samples have data">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">
##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA00001 NA00002 NA00003
16 50745926 rs2066844 C T 80 PASS NS=3;DP=14;AF=0.5;DB;H2 GT:GQ:DP:HQ 0|0:48:1:51,51 1|0:48:8:51,51 1/1:43:5:.,.
20 14370 rs6054257 G A 29 PASS NS=3;DP=14;AF=0.5;DB;H2 GT:GQ:DP:HQ 0|0:48:1:51,51 1|0:48:8:51,51 1/1:43:5:.,.
20 17330 . T A 3 q10 NS=3;DP=11;AF=0.017 GT:GQ:DP:HQ 0|0:49:3:58,50 0|1:3:5:65,3 0/0:41:3
20 1110696 rs6040355 A G,T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27 2|1:2:0:18,2 2/2:35:4
20 1230237 . T G 47 PASS NS=3;DP=13;AA=T GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:51,51 0/0:61:2
20 1230288 . T . 50 PASS NS=3;DP=13;AA=T GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:51,51 0/0:61:2
20 1234567 microsat1 GTCT G,GTACT 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4 0/2:17:2 1/1:40:3

此VCF文件的详细说明在此处给出。您会看到文件中包含七个基因座以及许多注释行。

现在让我们进行转换：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf > ex2.avinput
WARNING to old ANNOVAR users: this program no longer does line-to-line conversion for multi-sample VCF files. If you want to include all variants in output, use '-format vcf4old' or use '-format vcf4 -allsample -withfreq' instead.
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 2 SNPs (1 transitions and 1 transversions) and 1 indels/substitutions for 1 sample (but input contains 3 samples)
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.avinput 
20 1110696 1110696 A G het 67 6
20 1110696 1110696 A T het 67 6
20 1234568 1234570 TCT - het 50 4

上面的命令将ex2.vcf作为输入文件，ex2.avinput作为输出文件。另外3列是zygosity status（接合性状态）, genotype quality 和 read depth。

如果您仔细阅读屏幕消息，则看到该VCF文件中仅处理了3个样本中的1个。默认情况下，VCF文件中的仅第一个样本将被写入输出文件。输入包含七个基因座，但是对于第一个样本，它们中的许多不具有非参考基因型，这就是为什么输出仅包含3个变异的原因。

除了使用重定向之外，我们还可以执行以下操作将输出写入ex2.avinput文件中：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2.avinput

现在，如果我们要将所有3个样本都写入输出文件（作为三个单独的输出文件），则可以添加-allsample参数

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2 -allsample
NOTICE: output files will be written to ex2.<samplename>.avinput
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 10 SNPs (6 transitions and 4 transversions) and 3 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.NA00001.avinput 
20 1110696 1110696 A G het 67 6
20 1110696 1110696 A T het 67 6
20 1234568 1234570 TCT - het 50 4

[kaiwang@biocluster ~/]$ cat ex2.NA00002.avinput 
16 50745926 50745926 C T het 80 8
20 14370 14370 G A het 29 8
20 17330 17330 T A het 3 5
20 1110696 1110696 A G het 67 0
20 1110696 1110696 A T het 67 0
20 1234567 1234570 GTCT GTACT het 50 2

[kaiwang@biocluster ~/]$ cat ex2.NA00003.avinput 
16 50745926 50745926 C T hom 80 5
20 14370 14370 G A hom 29 5
20 1110696 1110696 A T hom 67 4
20 1234568 1234570 TCT - hom 50 3

显然，输入的VCF文件具有七个基因座。其中之一具有未知的替代等位基因，因此将其从输出中排除。另外两个都有两个备选等位基因，因此总共有9个变异，包括6个SNP和3个indel。生成三个输出文件。

请注意NA0002中的chr20：1110696。该基因座有两个备选等位基因，因此存在两个输出线，每个输出线用于一个突变。它们都是杂合子：

20 1110696 1110696 A G het 67 0
20 1110696 1110696 A T het 67 0

如果需要输出文件包含所有信息作为VCF文件，则可以添加-includeinfo参数：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2 -allsample -includeinfo
NOTICE: output files will be written to ex2.<samplename>.avinput
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 10 SNPs (6 transitions and 4 transversions) and 3 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.NA00001.avinput
20 1110696 1110696 A G 20 1110696 rs6040355 A G 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1110696 1110696 A T 20 1110696 rs6040355 A T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1234568 1234570 TCT - 20 1234567 microsat1 GTCT G 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4

如果在输出行中还需要接合性，质量和reads覆盖率信息，请添加-withzyg参数：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2 -allsample -include -withzyg
NOTICE: output files will be written to ex2.<samplename>.avinput
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 10 SNPs (6 transitions and 4 transversions) and 3 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.NA00001.avinput 
20 1110696 1110696 A G het 67 6 20 1110696 rs6040355 A G 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1110696 1110696 A T het 67 6 20 1110696 rs6040355 A T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1234568 1234570 TCT - het 50 4 20 1234567 microsat1 GTCT G 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4

现在，有一个非常重要的论点：-withfreq。当-withfreq被设置时，它会打印出每个SNP的等位基因频率在VCF文件，基于文件内的所有样本。因为我们没有整体看待所有样本，所以这里不考虑单个基因型，因此输出文件应包含输入文件中的所有基因座：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2.avinput -allsample -withfreq
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 17 SNPs (8 transitions and 9 transversions) and 6 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.avinput 
16 50745926 50745926 C T 0.5 80 5
20 14370 14370 G A 0.5 29 5
20 17330 17330 T A 0.1667 3 3
20 1110696 1110696 A G 0.5 67 0
20 1110696 1110696 A T 0.6667 67 4
20 1230237 1230237 T G 0 47 2
20 1230288 1230288 T 0 0 50 2
20 1234568 1234570 TCT - 0.75 50 3
20 1234567 1234570 GTCT GTACT 0.5 50 2

实际上，您可能想要添加，-includeinfo以便所有样本的所有基因型记录都包含在最终输出文件中：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2.avinput -allsample -withfreq -include
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 17 SNPs (8 transitions and 9 transversions) and 6 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.avinput
16 50745926 50745926 C T 0.5 80 5 16 50745926 rs2066844 C T 80 PASS NS=3;DP=14;AF=0.5;DB;H2 GT:GQ:DP:HQ 0|0:48:1:51,51 1|0:48:8:51,51 1/1:43:5:.,.
20 14370 14370 G A 0.5 29 5 20 14370 rs6054257 G A 29 PASS NS=3;DP=14;AF=0.5;DB;H2 GT:GQ:DP:HQ 0|0:48:1:51,51 1|0:48:8:51,51 1/1:43:5:.,.
20 17330 17330 T A 0.1667 3 3 20 17330 . T A 3 q10 NS=3;DP=11;AF=0.017 GT:GQ:DP:HQ 0|0:49:3:58,50 0|1:3:5:65,3 0/0:41:3
20 1110696 1110696 A G 0.5 67 0 20 1110696 rs6040355 A G,T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27 2|1:2:0:18,2 2/2:35:4
20 1110696 1110696 A T 0.6667 67 4 20 1110696 rs6040355 A G,T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27 2|1:2:0:18,2 2/2:35:4
20 1230237 1230237 T G 0 47 2 20 1230237 . T G 47 PASS NS=3;DP=13;AA=T GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:51,51 0/0:61:2
20 1230288 1230288 T 0 0 50 2 20 1230288 . T . 50 PASS NS=3;DP=13;AA=T GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:51,51 0/0:61:2
20 1234568 1234570 TCT - 0.75 50 3 20 1234567 microsat1 GTCT G,GTACT 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4 0/2:17:2 1/1:40:3
20 1234567 1234570 GTCT GTACT 0.5 50 2 20 1234567 microsat1 GTCT G,GTACT 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4 0/2:17:2 1/1:40:3

同样，请注意，对于chr20：1110696有两行，因为在同一位点有两个等位基因。同样，同一VCF记录有两行（chr20：1234568和chr20：1234567），因为同一位点有两个等位基因（indels）。

如果您需要VCF标头信息，只需添加-comment参数：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format vcf4 example/ex2.vcf -outfile ex2 -allsample -include -comment
NOTICE: output files will be written to ex2.<samplename>.avinput
NOTICE: Finished reading 25 lines from VCF file
NOTICE: A total of 7 locus in VCF file passed QC threshold, representing 6 SNPs (3 transitions and 3 transversions) and 3 indels/substitutions
NOTICE: Finished writing 10 SNPs (6 transitions and 4 transversions) and 3 indels/substitutions for 3 samples
WARNING: Skipped 1 invalid alternative alleles found in input file

[kaiwang@biocluster ~/]$ cat ex2.NA00001.avinput
##fileformat=VCFv4.0
##fileDate=20090805
##source=myImputationProgramV3.1
##reference=1000GenomesPilot-NCBI36
##phasing=partial
##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">
##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">
##FILTER=<ID=q10,Description="Quality below 10">
##FILTER=<ID=s50,Description="Less than 50% of samples have data">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">
##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA00001
20 1110696 1110696 A G 20 1110696 rs6040355 A G 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1110696 1110696 A T 20 1110696 rs6040355 A T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27
20 1234568 1234570 TCT - 20 1234567 microsat1 GTCT G 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 0/1:35:4

请注意，输出文件仅包含该特定样本的VCF信息（换句话说，输出文件中仅样本NA00001的GT：GQ：DP：HQ信息ex2.NA00001.avinput）。

因此，在实践中，您可以convert2annovar.pl通过将ANNOVAR输入文件反向转换为特定样本的VCF文件，将庞大的VCF文件分成许多部分，每个部分针对一个特定样本。从某种意义上说，它类似于vcf-subsetVCFtools中的程序，但是在这种情况下，ANNOVAR会更加高效，尤其是在处理大型VCF文件（如超过1TB的文件）时。见下文：

[kaiwang@biocluster ~/]$ grep -P '^#' example1.NA00001.avinput > example1.NA00002.vcf
[kaiwang@biocluster ~/]$ grep -v -P '^#' example1.NA00001.avinput | cut -f 6- >> example1.NA00002.vcf
[kaiwang@biocluster ~/]$ grep -P '^#' example1.NA00002.avinput > example1.NA00002.vcf
[kaiwang@biocluster ~/]$ grep -v -P '^#' example1.NA00002.avinput | cut -f 6- >> example1.NA00002.vcf
[kaiwang@biocluster ~/]$ grep -P '^#' example1.NA00003.avinput > example1.NA00003.vcf
[kaiwang@biocluster ~/]$ grep -v -P '^#' example1.NA00003.avinput | cut -f 6- >> example1.NA00003.vcf

[kaiwang@biocluster ~/]$ cat example1.NA00003.vcf 
##fileformat=VCFv4.0
##fileDate=20090805
##source=myImputationProgramV3.1
##reference=1000GenomesPilot-NCBI36
##phasing=partial
##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">
##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">
##FILTER=<ID=q10,Description="Quality below 10">
##FILTER=<ID=s50,Description="Less than 50% of samples have data">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">
##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA00003
20 14370 rs6054257 G A 29 PASS NS=3;DP=14;AF=0.5;DB;H2 GT:GQ:DP:HQ 1/1:43:5:.,.
20 1110696 rs6040355 A T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 2/2:35:4
20 1234567 microsat1 GTCT G 50 PASS NS=3;DP=9;AA=G GT:GQ:DP 1/1:40:3

-dbSNP标识符

许多用户都有dbSNP rs标识符的列表，并想注释这些SNP的功能。这可以通过convert2annovar.pl以下-format rsid参数实现：

[kaiwang@biocluster ~/]$ cat example/snplist.txt 
rs74487784
rs41534544
rs4308095
rs12345678

[kaiwang@biocluster ~/]$ convert2annovar.pl -format rsid example/snplist.txt -dbsnpfile humandb/hg19_snp138.txt > snnplist.avinput 
NOTICE: Scanning dbSNP file humandb/hg19_snp138.txt...
NOTICE: input file contains 4 rs identifiers, output file contains information for 4 rs identifiers
WARNING: 1 rs identifiers have multiple records (due to multiple mapping) and they are all written to output

[kaiwang@biocluster ~/]$ cat snplist.avinput 
chr2 186229004 186229004 C T rs4308095
chr7 6026775 6026775 T C rs41534544
chr7 6777183 6777183 G A rs41534544
chr9 3901666 3901666 T C rs12345678
chr22 24325095 24325095 A G rs74487784

如上所示，新文件的前五列为chr，start，end，ref，alt，第六列为dbSNP标识符。LOG消息告诉我们1个SNP（rs41534544）具有与基因组的多个映射，因此，它在输出文件中具有两个条目。

-基因组区域中的所有变体

假设我有兴趣注释3 bp基因组区域chr1：2000001-2000003中的所有SNP，1 bp插入和1 bp缺失。convert2annovar.pl现在可以轻松完成此操作：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format region -seqdir humandb/hg19_seq/ chr1:2000001-2000003
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr1.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
1 2000001 2000001 A C
1 2000001 2000001 A G
1 2000001 2000001 A T
1 2000002 2000002 T A
1 2000002 2000002 T C
1 2000002 2000002 T G
1 2000003 2000003 C A
1 2000003 2000003 C G
1 2000003 2000003 C T

要为此添加x-bp插入和删除，请使用-inssize x和-delsize x，其中x是整数。

[kaiwang@biocluster ~/]$ convert2annovar.pl -format region -seqdir humandb/hg19_seq/ chr1:2000001-2000003 -inssize 1 -delsize 2
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr1.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
1 2000001 2000001 A C
1 2000001 2000001 A G
1 2000001 2000001 A T
1 2000002 2000002 T A
1 2000002 2000002 T C
1 2000002 2000002 T G
1 2000003 2000003 C A
1 2000003 2000003 C G
1 2000003 2000003 C T
1 2000001 2000001 AT -
1 2000002 2000002 TC -
1 2000001 2000001 - A
1 2000001 2000001 - C
1 2000001 2000001 - G
1 2000001 2000001 - T
1 2000002 2000002 - A
1 2000002 2000002 - C
1 2000002 2000002 - G
1 2000002 2000002 - T
1 2000003 2000003 - A
1 2000003 2000003 - C
1 2000003 2000003 - G
1 2000003 2000003 - T

如果您只想删除2个bp，而不想要SNV和插入，请使用-subsize 0 -delsize 2。默认情况下，-subsize设置为1表示始终需要SNV。

-成绩单中的所有变体

与-format region上述类似，用户可以为转录本的外显子中所有可能的变体（加上剪接变体）生成ANNOVAR输入文件。让我们以NM_022162为例：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format transcript NM_022162 -gene humandb/hg19_refGene.txt -seqdir humandb/hg19_seq/ > NM_022162.avinput
NOTICE: 12 regions will be analyzed for possible mutations
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout
NOTICE: Reading region from STDIN ... Done with 1 regions from 1 chromosomes
NOTICE: Finished reading 1 sequences from humandb/hg19_seq/chr16.fa
NOTICE: Finished writting FASTA for 1 genomic regions to stdout

默认情况下，还会处理剪接位点（每个外显子外部2 bp）。如果您不希望它们，或者想更改拼接的定义，请使用-splicing_threshold参数。

类似于-format region以上，-delsize并且-dupsize和-subsize全部的支持，这样就可以自定义输出文件。

使用的一项重大实用功能-format region是用户可以将蛋白质或cDNA突变命名法反向转换回基因组坐标。为此，只需要为一个脚本生成所有可能的突变，然后注释生成的* .avinput文件，然后在该文件中搜索匹配项即可。

-SAMtools堆积格式

该部分现在已过时，实际上samtools现在使用mpileup，而不是“旧的”堆积。

请注意，有许多不同的堆积格式，但是在这里，我们要处理（自2011年起已过时）“基因型调用”堆积，其中一列中包含变体调用。Samtools网站上提供了更详细的描述。生成“基因型调用”堆积文件的示例如下所示：

samtools pileup -vcf ref.fa aln.bam > raw.pileup

该命令生成包含MAQ中实现的模型的共识调用的堆积文件（当然，还有许多其他指定的SNP调用者可用，用户可以自由选择）。由SamTools生成的示例基因型调用堆积格式如下所示：

chr1 556674 G G 54 0 60 16 a,.....,...,.... (B%A+%7B;0;%=B<:
chr1 556675 C C 55 0 60 16 ,,..A..,...,.... CB%%5%,A/+,%....
chr1 556676 C C 59 0 60 16 g,.....,...,.... .B%%.%.?.=/%...1
chr1 556677 G G 75 0 60 16 ,$,.....,...,.... .B%%9%5A6?)%;?:<
chr1 556678 G K 60 60 60 24 ,$.....,...,....^~t^~t^~t^~t^~t^~t^~t^~t^~t B%%B%<A;AA%??<=??;BA%B89
chr1 556679 C C 61 0 60 23 .....a...a....,,,,,,,,, %%1%&?*:2%*&)(89/1A@B@@
chr1 556680 G K 88 93 60 23 ..A..,..A,....ttttttttt %%)%7B:B0%55:7=>>A@B?B;
chr1 556681 C C 102 0 60 25 .$....,...,....,,,,,,,,,^~,^~. %%3%.B*4.%.34.6./B=?@@>5.
chr1 556682 A A 70 0 60 24 ...C,...,....,,,,,,,,,,. %:%(B:A4%7A?;A><<999=<<'
chr1 556683 G G 99 0 60 24 ....,...,....,,,,,,,,,,. %A%3B@%?%C?AB@BB/./-1A7?

这些列是染色体，基于1的坐标，参考碱基，共有碱基（核苷酸的IUPAC命名法），共有质量，SNP质量，覆盖该位点的最大阅读质量，覆盖该位点的阅读数，各个碱基和碱基品质。

该convert2annovar.pl程序可以将堆积文件格式转换为ANNOVAR输入文件。默认情况下，-snpqual 20将强加该参数，以便仅处理质量得分> = 20的SNP并将其写入输出文件。输出的varlist文件包含ANNOVAR格式的被称为突变（非突变显然不在输出文件中）。

在2011年的Januaray版本的ANNOVAR中，用于处理堆积文件的格式已经非常成熟/固定。请注意，前五列符合标准ANNOVAR输入格式，而第六列和后续列则提供有关等位基因的信息。

[kaiwang@biocluster ~/]$ convert2annovar.pl 84060.pileup -coverage 10 | head
NOTICE: the default --format argument is set as 'pileup'
NOTICE: the default --snpqual argument for pileup format is set as 20
NOTICE: Column 6-9 in output are heterozygosity status, SNP quality, total reads, reads with mutation
1 20139 20140 CA - het 563 53 10
1 59374 59374 A G hom 129 39 37
1 798677 798677 T G het 30 26 4
1 798785 798785 G A hom 132 38 37
1 798791 798791 C T hom 156 46 45
1 799544 799544 G A het 35 39 7
1 799550 799550 G C het 64 38 8
1 799595 799595 T C het 24 28 5
1 861034 861034 A C het 46 14 4

上面的NOTICE行告诉用户输出中第6-9列的含义。在第一行中，我们看到一个Indel的深度覆盖率为53，其中有10个支持Indel。在第二行中，我们看到一个深度覆盖为39的SNP，其中37个支持替代等位基因（G）。第6列之后的这些附加数字可帮助用户确定变量调用是否可靠。

该-fraction参数可用于过滤出所有读段中备选等位基因百分比过低的变体。例如，如果我们假设所有变体调用必须得到覆盖网站的至少40％的读取的支持，则可以使用：

[kaiwang@biocluster ~/]$ convert2annovar.pl 84060.pileup -coverage 10 -fraction 0.4 | head
NOTICE: the default --format argument is set as 'pileup'
NOTICE: the default --snpqual argument for pileup format is set as 20
NOTICE: Column 6-9 in output are heterozygosity status, SNP quality, total reads, reads with mutation
1 59374 59374 A G hom 129 39 37
1 798785 798785 G A hom 132 38 37
1 798791 798791 C T hom 156 46 45
1 871781 871781 G A het 198 13 8
1 873954 873954 C G het 211 16 7
1 873964 873964 A C het 172 17 8
1 877423 877423 A C hom 84 19 19
1 877664 877664 A G hom 63 12 12
1 878502 878502 T C hom 69 14 14
1 881483 881483 C A het 228 24 10

通过比较两个输出文件可以看出，indel的第一行不再在输出中，因为10/53 <40％。

一些其他有用的参数包括：--altcov，它指定替代等位基因的最小覆盖率（--coverage指定所有读取的覆盖率，无论它们支持参考等位基因还是替代等位基因）；--maxcoverage，它指定了要打印此变体的最大覆盖级别；--includeinfo，它指定输入行中的所有信息都应通过将它们附加在打印列之后而包含在输出行中。

程序完成后，它将打印出一些统计信息。通常，对于人类的全基因组测序，杂合子：纯合子的比例应为2：1左右，过渡：转化的比例应为2：1。

高级说明：当染色体为“ M”时，ANNOVAR不会打印出“ hom”或“ het”，而是会打印出介于0和1之间的数字，表明支持替代等位基因的读数的百分比。-chrmt如果您的路线中线粒体未标注为M，请使用参数。

-完整的基因组学格式

完整的基因组公司为客户提供许多基因分型调用文件。其中是一个var * ASM.tsv文件，如下所示。

[kai@beta ~/]$ head -n 20 var-GS000000088-ASM.tsv
#BUILD 1.5.0.5
#GENERATED_AT 2009-Nov-03 19:52:21.722927
#GENERATED_BY dbsnptool
#TYPE VAR-ANNOTATION
#VAR_ANN_SET /Proj/Pipeline/Production_Data/REF/HUMAN-F_06-REF/dbSNP.csv
#VAR_ANN_TYPE dbSNP
#VERSION 0.3
>locus ploidy haplotype chromosome begin end varType reference alleleSeq totalScore hapLink xRef
1 2 all chr1 0 959 no-call = ?
2 2 all chr1 959 972 = = =
3 2 all chr1 972 1001 no-call = ?
4 2 all chr1 1001 1008 = = =
5 2 all chr1 1008 1114 no-call = ?
6 2 all chr1 1114 1125 = = =
7 2 all chr1 1125 1191 no-call = ?
8 2 all chr1 1191 1225 = = =
9 2 all chr1 1225 1258 no-call = ?
10 2 all chr1 1258 1267 = = =
12 2 all chr1 1267 1275 no-call = ?

该convert2annovar.pl程序可以使用-format cg参数将该文件转换为ANNOVAR格式。输出文件如下所示：

[kai@beta ~/]$ head var-GS000000088-ASM.tsv.snp 
1 28095 28095 A G snp 67 dbsnp:rs806727 hom
1 31844 31844 A G snp 133 dbsnp:rs806721 hom
1 37155 37155 T G snp 51 dbsnp:rs2691275 het 
1 44449 44449 T C snp 74 het 
1 45679 45679 G A snp 191 dbsnp:rs3020699 het 
1 45713 45713 C G snp 191 het 
1 45789 45789 T C snp 138 dbsnp:rs3020698 hom
1 46662 46662 T C snp 69 dbsnp:rs2691309 het 
1 47109 47109 C G snp 56 dbsnp:rs2691313 het 
1 47815 47815 A C snp 67 dbsnp:rs2691334 hom

下面是一个示例命令行会话：

[kaiwang@biocluster ~/]$ convert2annovar.pl -format cg -out GS000000455.query var-GS000000455-ASM.tsv
NOTICE: Converting variants from var-GS000000455-ASM.tsv
NOTICE: Done with 25667914 lines

[kaiwang@biocluster ~/]$ wc -l GS000000455.query 
3728645 GS000000455.query

在此示例中，处理了来自Complete Genomics数据的var * ASM.tsv文件中的2560万行，并将370万个变体以ANNOVAR输入格式写入了输出文件。

-GFF3-SOLiD调用格式

由于当今很少有人处理SOLiD数据，因此本节也被删除。SOLiD BioScope生成GFF3格式的基因型调用，在此我们将其称为GFF3-SOLiD。（此输入文件不应与GFF3注释数据库混淆，因为它们具有不同的用途。在这里，我们仅处理输入文件。）例如，SOLiD提供以下格式的SNP变体调用：

[kaiwang@biocluster ~/]$ head -n 20 var/Yoruban_snp_18x.gff 
##gff-version 3
##solid-gff-version 0.3
##source-version 2
##type DNA
##date 2009-03-13
##time 0:0:0
##feature-ontology http://song.cvs.sourceforge.net/*checkout*/song/ontology/sofa.obo?revision=1.141
##reference-file 
##input-files /data/results3/yoruban_strikes_back/files_for_Aaron/NA18507_18x_SOLiD_SNP_calls_hg18_dbSNP129_annotated.txt
##run-path 
1 AB_SOLiD SNP caller SNP 997 997 1 . . coverage=3;ref_base=A; ref_score=0.3359;ref_confi=0.9528;ref_single=0/0;ref_paired=1/1;consen_base=G; consen_score=0.6641;consen_confi=0.9420;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 1371 1371 1 . . coverage=2;ref_base=A; ref_score=0.0000;ref_confi=0.0000;ref_single=0/0;ref_paired=0/0;consen_base=G; consen_score=1.0000;consen_confi=0.8717;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 2061 2061 1 . . coverage=2;ref_base=G; ref_score=0.0000;ref_confi=0.0000;ref_single=0/0;ref_paired=0/0;consen_base=C; consen_score=1.0000;consen_confi=0.9138;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 4770 4770 1 . . coverage=2;ref_base=A; ref_score=0.0000;ref_confi=0.0000;ref_single=0/0;ref_paired=0/0;consen_base=G; consen_score=1.0000;consen_confi=0.8699;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 4793 4793 1 . . coverage=16;ref_base=A; ref_score=0.0689;ref_confi=0.9384;ref_single=0/0;ref_paired=1/1;consen_base=G; consen_score=0.6858;consen_confi=0.8494;consen_single=0/0;consen_paired=11/10
1 AB_SOLiD SNP caller SNP 5074 5074 1 . . coverage=4;ref_base=T; ref_score=0.5165;ref_confi=0.9355;ref_single=2/2;ref_paired=0/0;consen_base=K; consen_score=0.4835;consen_confi=0.8759;consen_single=1/1;consen_paired=1/1
1 AB_SOLiD SNP caller SNP 6241 6241 1 . . coverage=5;ref_base=T; ref_score=0.4450;ref_confi=0.9383;ref_single=0/0;ref_paired=2/2;consen_base=Y; consen_score=0.3646;consen_confi=0.7688;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 9089 9089 1 . . coverage=5;ref_base=T; ref_score=0.4260;ref_confi=0.9483;ref_single=1/1;ref_paired=1/1;consen_base=W; consen_score=0.3639;consen_confi=0.8100;consen_single=1/1;consen_paired=1/1
1 AB_SOLiD SNP caller SNP 9131 9131 1 . . coverage=8;ref_base=C; ref_score=0.7547;ref_confi=0.9306;ref_single=3/3;ref_paired=3/3;consen_base=Y; consen_score=0.2453;consen_confi=0.9073;consen_single=0/0;consen_paired=2/2
1 AB_SOLiD SNP caller SNP 18426 18426 1 . . coverage=3;ref_base=A; ref_score=0.0000;ref_confi=0.0000;ref_single=0/0;ref_paired=0/0;consen_base=G; consen_score=1.0000;consen_confi=0.8163;consen_single=0/0;consen_paired=3/3

可以使用-format gff3-solid参数完成转换。

[kaiwang@biocluster ~/]$ convert2annovar.pl var/Yoruban_snp_18x.gff -format gff3-solid | head
1 997 997 A G hom
1 1371 1371 A G hom
1 2061 2061 G C hom
1 4770 4770 A G hom
1 4793 4793 A G hom
1 5074 5074 T G het
1 6241 6241 T C het
1 9089 9089 T A het
1 9131 9131 C T het
1 18426 18426 A G hom

添加--includeinfo参数将打印出带有调用的详细属性的附加列。

-SOAPsnp调用格式

短寡核苷酸分析软件包（SOAP）套件是由BGI开发的，而SOAPsnp是生成变体调用的组件。下面给出了基因型调用文件的示例：

chr10 84026 G R 55 A 32 9 9 G 29 3 5 14 0.275000 1.42857 1 81
chr10 84541 C M 45 C 27 5 5 A 25 3 4 9 0.285714 1.11111 1 4
chr10 284953 A G 76 G 33 26 26 A 0 0 0 26 1.00000 1.00000 1 9472
chr10 313283 A R 99 G 32 15 15 A 28 8 8 23 0.162302 1.00000 1 28330
chr10 363048 T Y 99 T 26 12 12 C 30 10 10 22 0.461435 1.00000 1 14012
chr10 377060 G A 55 A 33 11 11 G 0 0 0 11 1.00000 1.00000 1 7654
chr10 384714 G A 74 A 33 18 18 G 0 0 0 18 1.00000 1.00000 1 552
chr10 418503 A G 34 G 34 4 4 A 0 0 0 4 1.00000 1.00000 1 7377
chr10 434997 C Y 73 T 26 5 5 C 26 4 4 9 0.682540 1.00000 1 64
chr10 435061 C T 36 T 30 16 16 C 0 0 0 16 1.00000 1.00000 1 64

该convert2annovar.pl程序可以处理这种格式，使用的-format soapsnp参数。下面给出了输出文件的示例：

10 84026 84026 G A het
10 84541 84541 C A het
10 284953 284953 A G hom
10 313283 313283 A G het
10 363048 363048 T C het
10 377060 377060 G A hom
10 384714 384714 G A hom
10 418503 418503 A G hom
10 434997 434997 C T het
10 435061 435061 C T hom

请注意，如果使用--includeinfo参数，则来自输入文件的所有信息将包含在输出文件中。

-MAQ呼叫格式

该convert2annovar.pl程序可以处理这种格式，使用的-format maq参数。SNP和插入/缺失均可被正确处理。

-CASAVA呼叫格式

由于CASAVA调用文件本身不包含染色体信息convert2annovar.pl，因此程序可以使用-format casava参数并按参数指定染色体来处理这种格式--chr。SNP和插入/缺失均可被正确处理。

https://doc-openbio.readthedocs.io/projects/annovar/en/latest/user-guide/gene/

ANNOVAR | 变异注释【下】

ANNOVAR | 变异注释【下】

结果文件的解读

输出文件1：*.variant_function

输出文件2（refSeq gene annotation）

包括RefSeq基因版本号

注释线粒体变异

切换到UCSC / Ensembl基因注释

切换到hg38 Ensembl基因注释

切换到GENCODE / CCDS基因注释

为非人类物种创建自己的基因定义数据库

那么新物种的GFF3文件呢？

了解ANNOVAR如何解决基因定义中的罕见问题

有用的教程

table_annovar.pl

annotate_variation.pl

额外资源

VCF文件

ANNOVAR输入文件

格式转换

-VCF4格式

-dbSNP标识符

-基因组区域中的所有变体

-成绩单中的所有变体

-SAMtools堆积格式

-完整的基因组学格式

-GFF3-SOLiD调用格式

-SOAPsnp调用格式

-MAQ呼叫格式

-CASAVA呼叫格式