一文读懂Manta的结果文件

2019-04-16 本文已影响13人徐广惠_6f76

Manta输出VCF文件

Manta运行完毕后，将在$ {MANTA_ANALYSIS_PATH}/results/variants目录下输出一组VCF格式的结果文件。

如果用户使用的是germline的检测模式，结果文件将包括： diploidSV.vcf.gz，candidateSV.vcf.gz 和 candidateSmallIndels.vcf.gz。

Germline 模式
如果用户使用的是somatic检测模式中的Tumor-Normal，结果文件将包括： somaticSV.vcf.gz，diploidSV.vcf.gz，candidateSV.vcf.gz和 candidateSmallIndels.vcf.gz 。
如果用户使用的是somatic检测模式中的Tumor-Only，结果文件将包括： tumorSV.vcf.gz，candidateSV.vcf.gz 和 candidateSmallIndels.vcf.gz。

Somatic 模式

无论是diploidSV.vcf.gz，somaticSV.vcf.gz还是tumorSV.vcf.gz，他们描述sv的规则是一致的，只是在记录的信息上略有不同。如，

基因型判定信息：somaticSV.vcf.gz和tumorSV.vcf.gz不包含基因型判定的相关信息，例如GT, GQ, PL等；
打分信息：diploidSV.vcf.gz的胚系突变打分展示在QUAL中，somaticSV.vcf.gz中的体细胞变异打分展示在FORMAT的SOMATICSCORE中，而tumorSV.vcf.gz中不包含打分信息，需要自己通过PR和SR信息进行筛选，得到较为可靠的SV。

使用gzip -d -c *.file.gz > *.file命令可生成解压缩的VCF文件。

输出VCF中记录的SV类型

片段缺失（Deletion)

对于大的片段缺失，在VCF中ALT一列会有<DEL>的标志，ID中将以MantaDEL开头，使用grep "<DEL>" diploidSV.vcf命令可以直接将这一类的变异提取出来。CHROM和POS中记录的是该Deletion在参考基因组上的起始位置，FORMAT中END记录的是Deletion在参考基因组上的终止位置，SVLEN记录的是缺失片段的长度。
FORMAT中的PR和SR记录的是支持REF和ALT基因型的Paired Reads数和Split Reads数。
在diploidSV.vcf中还会在FORMAT中包含基因型相关的信息，如GT，GQ， PL等（不懂这些概念？请参考：https://software.broadinstitute.org/gatk/documentation/article.php?id=1268）。

# diploidSV.vcf
CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Sample1
1       15552819        MantaDEL:1225:0:1:0:0:0 G       <DEL>   442     PASS    END=15563511;SVTYPE=DEL;SVLEN=-10692;SVINSLEN=2;SVINSSEQ=TA     GT:FT:GQ:PL:PR:SR       0/1:PASS:334:492,0,331:15,11:15,8
# somaticSV.vcf
CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Normal     Tumor
3       122354925       MantaDEL:7341:0:1:0:0:0 A       <DEL>   .       MinSomaticScore END=185785642;SVTYPE=DEL;SVLEN=-63430717;SVINSLEN=1;SVINSSEQ=C;SOMATIC;SOMATICSCORE=16     PR:SR   47,0:149,0      283,2:975,4

片段插入（Insertions with incomplete insert sequence assembly）

对于大的片段插入，Manta会在CHROM和POS中记录DNA片段的插入位置，并在ALT中加入<INS>的标志，ID中将以MantaINS开头。这里插入的“DNA片段”，个人理解指的是外源的DNA片段，即无法比对到参考基因组，或者无法比对到参考基因组唯一位置。因此，Manta只能通过断点附近的reads得到插入片段两端的序列，但无法将整个插入片段的序列组装起来（如果有不同意见，欢迎留言讨论）。INFO中的LEFT_SVINSSEQ和RIGHT_SVINSSEQ给出了插入片段左右两端的序列信息。

CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Sample1
1       11031132        MantaINS:5:22234:22234:0:3:0    A       <INS>   999     PASS    END=11031132;SVTYPE=INS;CIPOS=0,20;CIEND=0,20;HOMLEN=20;HOMSEQ=GAGGCAGAGGCTGCAGTGAG;LEFT_SVINSSEQ=GAGGCAGAGGCTGCAGTGAGTCCAGCCTGGGGGACAGAGTGAGACCCTGTCTCAAAAAGAAAAAAAAAACAGCATAGGCACTGGTGTCAGTAGGCATCTGGGTTTGAATCCCACCTCTGTTGTGTGTATGTGTGTGTGTGTGTGTGTGTACCTGTTGCTTAGTTTCAGTTTATTTCTGTGAGTTGATTGTATGATAATGATGGTGATGATAGTAATAATAGTGATGGTAGTAGAGGGATGATATTGATGGTGATGGTGGTGATGATGATGTGAATGGTGGTGATGATAGTGATGGTGGTGATGGTGGTGATGATGATGGTGATGGTGACAATCATGGTAGTGATGGTCACAGTGATGATGGTGCTGGTGATGGTGGTGATGATGGTGTTAATGGTGGTGAT;RIGHT_SVINSSEQ=GACATGGATTATGGGATACTCACGTGTACTTTAAAAAATACAGGCTGGGGCCGAGCACGGTGGCTCACGCCTGTAACCCCAGCACTTTGGGAGGCCGAGGCGGGTGGATCACGAGGTCAGGAGTTCAAGACCAGCCTGGCCAACATGGCGAAACCCCATCTCTACTAAACATACAAAAATTAGCAGGGCATGGTGGTGTGTACCTGTAATCCCAGCTACCCAGGAGGCTGAGGCAGGAGAATCACTGGAACCCGG     GT:FT:GQ:PL:PR:SR       1/1:PASS:136:999,139,0:0,12:0,37

小的插入和缺失（Small indels）

Manta中，符合以下几个条件的插入或缺失会被归类于small indels:

该突变可以完全用插入序列和缺失序列来表示。
插入序列或缺失序列的长度小于1000bp。
有精确的变异的断点和插入/缺失序列。

虽然这些小的indels的ID也以MantaDEL或MantaINS开头，但在VCF中的表示方式和前述的DEL和INS不同，Manta将这些变异的完整的插入/缺失序列给在了REF或ALT中。并且会在INFO中增加CIGAR标签，对此类变异进行描述。

CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Sample1
1       2954348 MantaINS:244:0:0:0:0:0  A       ACCTGGGTCTCGTCTGCCACGGATTGCTCTCCGTGCTCCCCAGAGCGAGGTGCAGATGCCAGGGACCCTCTC        999     PASS    END=2954348;SVTYPE=INS;SVLEN=71;CIGAR=1M71I;CIPOS=0,17;HOMLEN=17;HOMSEQ=CCTGGGTCTCGTCTGCC  GT:FT:GQ:PL:PR:SR       1/1:PASS:88:999,91,0:0,0:0,33
1       1302326 MantaDEL:98:0:0:0:1:1   GAATGAGTGGATTGGTGAGTGAATTGGTGAGTTGAATTGGTGTGTGTAGTGGATGAGTGTGGATGAATGTGAATTGGCGAGTATGGATGTGTGAATTGGTGAGTGTGAATGTGTGGATTGGTGAGTGAATTGGTGAGTTGAATTGGTGTGTGTAGTGTGGATGAGTGTGAATTGGCGAGTGTGGATGAGTGTGAATTGGTGAGTGTG    GCAGTGTGAA      904     PASS    END=1302532;SVTYPE=DEL;SVLEN=-206;CIGAR=1M9I206D        GT:FT:GQ:PL:PR:SR  1/1:PASS:61:957,64,0:0,1:0,24

串联重复（Tandem Duplicate）

Manta没有办法检测散在重复（Dispersed duplications），但可以检出串联重复（Tandem Duplicate）。

串联重复（Tandem Duplicate）
在VCF结果文件中，串联重复的ID以MantaDUP:TANDEM开头，CHROM和POS记录开始位置，END记录结束为止。如下：

1       1413234 MantaDUP:TANDEM:123:0:1:0:0:0   C       <DUP:TANDEM>    514     PASS    END=1413364;SVTYPE=DUP;SVLEN=130;SVINSLEN=3;SVINSSEQ=TGT        GT:FT:GQ:PL:PR:SR 0/1:PASS:504:564,0,501:14,1:41,19

染色体易位（Translocation）

Manta对于染色体间易位和染色体内易位不做特殊区分，ID都以MantaBND开头，BND即breakend的缩写。在CHROM、POS中展示第一个断点位置，在ALT中展示第二个断点位置，例如：A]10:115172011]、[12:70547434[C。通过比较第一个断点和第二个断点的染色体，可以判断是染色体间易位还是染色体内易位（这里为了描述方便，使用了“第一个断点”、“第二个断点”的描述方式，事实上两个断点并没有顺序之分）。
值得注意的是，ALT中方括号的方向在判断融合基因中有重要的作用。...]...]指易位序列在第一个断点位置的3'端，[...[...指易位序列在第一个断点位置的5‘端，如下图。

...]...]指易位序列在第一个断点位置的右侧

[...[...指易位序列在第一个断点位置的左侧
Manta会用两条记录（也就是两个BND）来描述一个易位产生的新的连接点，这两条记录互为MATE关系，在FORMATA的MATEID标签可找到另一条记录的ID。如下，

CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Sample1
# example1
1       180903258       MantaBND:13113:0:1:0:0:0:0      C       C[3:48537167[   314     PASS    SVTYPE=BND;MATEID=MantaBND:13113:0:1:0:0:0:1;CIPOS=0,3;HOMLEN=3;HOMSEQ=GCA;BND_DEPTH=30;MATE_BND_DEPTH=31  GT:FT:GQ:PL:PR:SR       0/1:PASS:314:364,0,423:20,9:24,6
3       48537167        MantaBND:13113:0:1:0:0:0:1      G       ]1:180903258]G  314     PASS    SVTYPE=BND;MATEID=MantaBND:13113:0:1:0:0:0:0;CIPOS=0,3;HOMLEN=3;HOMSEQ=CAC;BND_DEPTH=31;MATE_BND_DEPTH=30  GT:FT:GQ:PL:PR:SR       0/1:PASS:314:364,0,423:20,9:24,6
# example2
12      34017350        MantaBND:114233:0:1:0:0:0:0     C       C]13:48856953]  58      PASS    SVTYPE=BND;MATEID=MantaBND:114233:0:1:0:0:0:1;IMPRECISE;CIPOS=-318,319;BND_DEPTH=35;MATE_BND_DEPTH=34      GT:FT:GQ:PL:PR  0/1:PASS:58:108,0,286:20,9
13      48856953        MantaBND:114233:0:1:0:0:0:1     A       A]12:34017350]  58      PASS    SVTYPE=BND;MATEID=MantaBND:114233:0:1:0:0:0:0;IMPRECISE;CIPOS=-287,288;BND_DEPTH=34;MATE_BND_DEPTH=35      GT:FT:GQ:PL:PR  0/1:PASS:58:108,0,286:20,9

染色体片段在易位的过程中，可能会平移并连接到另一段染色体上(见下图 variant a），也可能翻转之后再连接到另一段染色体上（见下图 variant b）。具体看上面的两个例子，其中example1的两条记录ALT中的方括号方向不一样，它对应的是variant a这种情况；example2的两条记录中方括号方向一致，对应的是variant b这种情况。

Translocation

染色体倒位（Inversion）

在默认情况下，Manta会用4条BND记录来表述一个倒位事件，并且这四条记录拥有相同的EVENT标签。下面是官网上给的例子：

chr1    17124941        MantaBND:1445:0:1:1:3:0:0       T       [chr1:234919886[T       999     PASS    SVTYPE=BND;MATEID=MantaBND:1445:0:1:1:3:0:1;CIPOS=0,1;HOMLEN=1;HOMSEQ=T;INV5;EVENT=MantaBND:1445:0:1:0:0:0:0;JUNCTION_QUAL=254;BND_DEPTH=107;MATE_BND_DEPTH=100 GT:FT:GQ:PL:PR:SR       0/1:PASS:999:999,0,999:65,8:15,51
chr1    17124948        MantaBND:1445:0:1:0:0:0:0       T       T]chr1:234919824]       999     PASS    SVTYPE=BND;MATEID=MantaBND:1445:0:1:0:0:0:1;INV3;EVENT=MantaBND:1445:0:1:0:0:0:0;JUNCTION_QUAL=999;BND_DEPTH=109;MATE_BND_DEPTH=83      GT:FT:GQ:PL:PR:SR       0/1:PASS:999:999,0,999:60,2:0,46
chr1    234919824       MantaBND:1445:0:1:0:0:0:1       G       G]chr1:17124948]        999     PASS    SVTYPE=BND;MATEID=MantaBND:1445:0:1:0:0:0:0;INV3;EVENT=MantaBND:1445:0:1:0:0:0:0;JUNCTION_QUAL=999;BND_DEPTH=83;MATE_BND_DEPTH=109      GT:FT:GQ:PL:PR:SR       0/1:PASS:999:999,0,999:60,2:0,46
chr1    234919885       MantaBND:1445:0:1:1:3:0:1       A       [chr1:17124942[A        999     PASS    SVTYPE=BND;MATEID=MantaBND:1445:0:1:1:3:0:0;CIPOS=0,1;HOMLEN=1;HOMSEQ=A;INV5;EVENT=MantaBND:1445:0:1:0:0:0:0;JUNCTION_QUAL=254;BND_DEPTH=100;MATE_BND_DEPTH=107 GT:FT:GQ:PL:PR:SR       0/1:PASS:999:999,0,999:65,8:15,51

但开发者另外提供了一个脚本$MANTA_INSTALL_FOLDER/libexec/convertInversion.py可以将BND记录的Inversion转换成另一种形式（见下），并以MantaINV作为ID的开头，每条记录表述一个新的连接点的信息，位置信息记录在CHROM和POS中。一条标准的Inversion应该有两连接点的记录，并且拥有相同的EVENT标签。

chr1    17124940        MantaINV:1445:0:1:1:3:0 C       <INV>   999     PASS    END=234919885;SVTYPE=INV;SVLEN=217794945;CIPOS=0,1;CIEND=-1,0;HOMLEN=1;HOMSEQ=T;EVENT=MantaINV:1445:0:1:0:0:0;JUNCTION_QUAL=254;INV5    GT:FT:GQ:PL:PR:SR       0/1:PASS:999:999,0,999:65,8:15,51
chr1    17124948        MantaINV:1445:0:1:0:0:0 T       <INV>   999     PASS    END=234919824;SVTYPE=INV;SVLEN=217794876;EVENT=MantaINV:1445:0:1:0:0:0;JUNCTION_QUAL=999;INV3   GT:FT:GQ:PL:PR:SR        0/1:PASS:999:999,0,999:60,2:0,46

另外，在Inversion的记录中，INFO中还提供了INV3、INV5两个标签，INV3指发生倒位的序列位于此记录报道的连接点的3'端，INV5指发生倒位的序列位于此记录报道的连接点的5'端。在IGV中，INV5标签对应的是"RR" reads，INV3标签对应的是"LL"reads（可参考我的另一篇笔记）。
需要注意的是，在实际应用中得到的VCF完成格式转换后，存在很多虽然标注为Manta:INV，但只有一条记录情况，因此实际上并不是一个标准的Inversion事件。

写在后面

不同的SV检测软件都有自己的一套描述规则，有很多细节值得琢磨，以后有新的体会再慢慢补充。

参考

https://github.com/Illumina/manta/blob/master/docs/userGuide/README.md