0

我有一个 vcf 文件，我正在尝试从这些列中提取信息：

#CHROM  POS   REF     ALT

但是，只有当 SAMPLE-1 列包含字符串 DeNovo（不是 DeNovoSV）并且 SAMPLE-1、SAMPLE-2 和 SAMPLE-3 都包含 PASS 时，我才想提取这些。

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  SAMPLE-1  SAMPLE-2  SAMPLE-3
chr1    10230   .       AC      A       186.90  .       AC=4;AF=0.667;AN=6;DP=77;FS=0.000;MQ=26.38;MQRankSum=0.436;QD=3.89;ReadPosRankSum=0.000;SOR=1.162       GT:AD:AF:DP:GQ:FT:F1R2:F2R1:PL:GP:PP
:DN 0/1:6,12:0.667:18:32:PASS:2,3:4,9:69,0,30:3.9669e+01,2.7888e-03,3.2724e+01:295,0,215:Inherited  0/1:5,15:0.750:20:11:PASS:3,6:2,9:60,0,8:3.0340e+01,3.6694e-01,1.0964e+01:172,0,137:.   1/1:0,10
:1.000:10:26:PASS:0,2:0,8:93,29,0:6.1212e+01,2.6342e+01,1.0101e-02:378,0,183:.
chr1    61871   .       C       CT      60.27   .       AC=3;AF=0.500;AN=6;DP=29;FS=11.290;MQ=33.00;MQRankSum=-0.423;QD=2.51;ReadPosRankSum=0.705;SOR=0.478     GT:AD:AF:DP:GQ:FT:F1R2:F2R1:PL:GP:PP:DPL:DN:DQ  0/0:5,0:0.000:5:15:PASS:.:.:0,15,182:.:0,7,93:0,15,100:DeNovo:2.9227e-07        1/1:0,2:1.000:2:5:PASS:0,1:0,1:42,6,0:2.4787e+01,4.7870e+00,1.7754e+00:29,0,9:40,6,0:.  0/1:15,7:0.318:22:26:PASS:6,3:9,4:43,0,41:2.6538e+01,9.8206e-03,4.4010e+01:65,0,250:74,0,234:.
chr1    66369   .       TA      T       116.77  .       AC=2;AF=0.500;AN=4;DP=56;FS=10.138;MQ=173.59;MQRankSum=1.468;QD=4.32;ReadPosRankSum=0.929;SOR=0.367     GT:AD:AF:DP:GQ:FT:F1R2:F2R1:PL:GP:PP    ./.:11,5:0.312:16:0:LowGQ:.:.   0/1:8,4:0.333:12:40:PASS:3,2:5,2:71,0,43:4.1762e+01,4.0824e-04,4.5625e+01:71,0,43       0/1:8,7:0.467:15:45:PASS:4,4:4,3:77,0,47:4.7400e+01,1.2244e-04,5.0000e+01:77,0,47
chr1    934273  .       G       C       8.67    .       AC=1;AF=0.167;AN=6;DP=26;FS=0.000;MQ=19.17;MQRankSum=-1.179;QD=0.96;ReadPosRankSum=1.666;SOR=0.223      GT:AD:AF:DP:GQ:FT:F1R2:F2R1:PL:GP:PP
:DPL:DN:DQ  0/1:7,2:0.222:9:11:PASS:5,1:2,1:45,0,32:1.0868e+01,3.7242e-01,3.5372e+01:45,0,60:45,0,32:DeNovoSV:4.3945e-09    0/0:6,3:0.333:9:0:LowGQ:.:.:0,0,140:.:55,0,191:46,0,186:.       0/0:8,0:
0.000:8:23:PASS:.:.:0,23,190:.:0,25,195:0,23,190:.
chr1    934274  .       G       C       8.68    .       AC=1;AF=0.167;AN=6;DP=26;FS=0.000;MQ=19.17;MQRankSum=-1.179;QD=0.96;ReadPosRankSum=1.666;SOR=0.223      GT:AD:AF:DP:GQ:FT:F1R2:F2R1:PL:GP:PP
:DPL:DN:DQ  0/1:7,2:0.222:9:11:PASS:5,1:2,1:45,0,32:1.0868e+01,3.7242e-01,3.5372e+01:45,0,60:45,0,32:DeNovoSV:4.3945e-09    0/0:6,3:0.333:9:0:PASS:.:.:0,0,140:.:55,0,191:46,0,186:.       0/0:8,0:
0.000:8:23:PASS:.:.:0,23,190:.:0,25,195:0,23,190:.

我试过使用 bcftools，见下文。

bcftools query -f '%CHROM %POS %REF %ALT\n' file.vcf | head -3
chr1 10230 AC A
chr1 61871 C CT
chr1 66369 TA T

有没有办法使用 bcftools，或者将它与 awk 结合以获得我正在寻找的 vcf 文件格式的输出？

非常感谢

4

1 回答 1

1

FWIW，我不熟悉，bcftools但如果打算将bcftools输出输入，awk那么我们可以在awk......

假设：

  • 没有字段（在文件中）包含空格
  • 总是要为列拉数据#CHROM///POSREFALT
  • 这 4 列将始终位于相同的字段/位置（#CHROM = $1/// ）POS = $2REF = $4ALT = $5
  • DeNovo总是出现在一对冒号之间（即，我们:DeNovo:在 field 中寻找$10）
  • 三个测试字段始终显示在相同的字段/位置 ( SAMPLE_1 = $10/ SAMPLE_2 = $11/ SAMPLE_3 = $12)
  • PASS总是出现在一对冒号之间（即，我们:PASS:在字段$10/ $11/中查找$12）

一个awk想法：

awk '
$10 ~ ":DeNovo:" &&
$10 ~ ":PASS:"   &&
$11 ~ ":PASS:"   &&
$12 ~ ":PASS:"      { print $1,$2,$4,$5 }
' file.vcf

这会产生：

chr1 61871 C CT

注意：这种方法的明显（？）缺点是我们对列引用进行了硬编码；我们当然可以修改awk代码以提供更动态的界面（复制bcftools功能？），但我不确定增加的复杂性是否值得付出努力，ymmv。

假设 OP 需要的功能bcftools（例如，指定一组可变的列），一个想法是修改当前bcftools调用以包含SAMPLE_X列，然后将该输出通过管道传输到awk; 就像是：

bcftools query -f '%CHROM %POS %REF %ALT %SAMPLE_1 %SAMPLE_2 %SAMPLE_3\n' file.vcf | awk '
$5 ~ ":DeNovo:" &&
$5 ~ ":PASS:"   &&
$6 ~ ":PASS:"   &&
$7 ~ ":PASS:"      { print $1,$2,$3,$4 }'

笔记：

  • 进一步的改进将包括参数化DeNovo和PASS字符串
  • 回到消除bcftools和完成整个事情的想法awk......可能需要两个列列表（要显示的列，要测试的列）以及一个字符串列表来比较每个“测试”列反对......可行但现在我们正在进行相当大量的编码工作（与我们迄今为止所做的相比）
于 2022-02-13T23:12:04.927 回答