Manta输出VCF文件
Manta运行完毕后,将在$ {MANTA_ANALYSIS_PATH}/results/variants
目录下输出一组VCF格式的结果文件。
-
如果用户使用的是germline的检测模式,结果文件将包括:
diploidSV.vcf.gz
,candidateSV.vcf.gz
和candidateSmallIndels.vcf.gz
。
如果用户使用的是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
。
无论是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)。
在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这种情况。
染色体倒位(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