VCF file, calculate depth per sample
1
1
Entering edit mode
2.3 years ago
pixie@bioinfo ★ 1.5k

Hello,

I have GBS data with a fairly large raw vcf file. I wanted to generate some quick custom stats. I wanted to find the sample depth for each sample by doing a row-sum across all the SNPs. For example, here it would be S1=875 and S2=1448 (i guess). It will be best if I can do this using some command line tools like bcftools because it will be part of a snakemake or bash pipeline. Also, I don't want to run into any memory issues here. Any other suggestions ?

Thanks

#CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  S1  S2
SL3.0ch01   83929063    .   A   G   999 PASS    DP=1667;VDB=0.0;RPB=0.999232;MQB=1.0;BQB=0.969803;SGB=293.924;MQ0F=0.0;AF1=0.517147;AC1=94.0;MQ=60;FQ=999.0;PV4=1.0,0.128142,1.0,1.0;G3=0.340041,0.286409,0.373549;HWE=3.7029e-05;DP4=758,0,909,0   GT:PL:DP:SP:ADF:ADR:AD  1/1:216,66,0:22:0:0,22:0,0:0,22 1/1:214,57,0:19:0:0,19:0,0:0,19
SL3.0ch01   97631153    .   A   T   999 PASS    DP=9192;VDB=0.0;RPB=0.996356;MQB=1.0;BQB=0.997647;SGB=960.891;MQ0F=0.0;AF1=0.477998;AC1=87.0;MQ=60;FQ=999.0;PV4=1.0,1.0,1.0,1.0;G3=0.362637,0.318703,0.318659;HWE=0.000491945;DP4=4817,0,4374,0 GT:PL:DP:SP:ADF:ADR:AD  0/1:198,0,193:122:0:59,63:0,0:59,63 1/1:255,255,0:115:0:0,115:0,0:0,115
SL3.0ch01   83952645    .   T   C   999 PASS    DP=7200;VDB=0.0;RPB=0.99006;MQB=0.999827;BQB=0.562521;SGB=1279.27;MQ0F=0.0;AF1=0.516484;AC1=94.0;MQ=60;FQ=999.0;PV4=1.0,1.61487e-06,0.180924,1.0;G3=0.340659,0.285715,0.373626;HWE=3.34968e-05;DP4=3269,0,3931,0    GT:PL:DP:SP:ADF:ADR:AD  1/1:195,175,0:104:0:4,100:0,0:4,100 1/1:234,216,0:85:0:1,84:0,0:1,84
SL3.0ch01   83936249    .   C   T   999 PASS    DP=14208;VDB=0.0;RPB=0.963987;MQB=0.998652;BQB=0.87956;SGB=676.78;MQ0F=0.0;AF1=0.516484;AC1=94.0;MQ=60;FQ=999.0;PV4=1.0,0.0640157,0.286743,1.0;G3=0.340659,0.285714,0.373626;HWE=3.34953e-05;DP4=6722,0,7486,0  GT:PL:DP:SP:ADF:ADR:AD  1/1:255,255,0:144:0:0,144:0,0:0,144 1/1:255,255,0:190:0:0,190:0,0:0,190
SL3.0ch01   96485259    .   C   A   999 PASS    DP=13689;VDB=0.0;RPB=0.973135;MQB=0.999874;BQB=0.964526;SGB=1322.73;MQ0F=0.0;AF1=0.532967;AC1=97.0;MQ=60;FQ=999.0;PV4=1.0,1.0,1.0,1.0;G3=0.296703,0.340659,0.362637;HWE=0.00240154;DP4=6373,0,7316,0    GT:PL:DP:SP:ADF:ADR:AD  1/1:255,255,0:180:0:0,180:0,0:0,180 1/1:255,255,0:168:0:0,168:0,0:0,168
SL3.0ch01   83969274    .   T   G   999 PASS    DP=84605;VDB=0.0;RPB=0.972479;MQB=0.999999;BQB=0.924048;SGB=-1202.62;MQ0F=0.0;AF1=0.494505;AC1=90.0;MQ=60;FQ=999.0;PV4=1.0,1.0,0.446962,0.375381;G3=0.362637,0.285714,0.351648;HWE=3.25784e-05;DP4=40264,0,44341,0  GT:PL:DP:SP:ADF:ADR:AD  1/1:255,255,0:1112:0:5,1107:0,0:5,1107  1/1:255,255,0:1083:0:3,1079:0,0:3,1079
SL3.0ch01   98281320    .   T   A   999 PASS    DP=713;VDB=0.0;RPB=0.349379;MQB=1.0;BQB=0.996143;SGB=-24.2855;MQ0F=0.0;AF1=0.279049;AC1=51.0;MQ=60;FQ=999.0;PV4=1.0,1.0,1.0,0.271776;G3=0.574036,0.261202,0.164762;HWE=0.00356772;DP4=592,0,121,0   GT:PL:DP:SP:ADF:ADR:AD  0/1:53,0,124:7:0:5,2:0,0:5,2    0/0:0,18,155:6:0:6,0:0,0:6,0
vcf bcftools • 617 views
ADD COMMENT
1
Entering edit mode
2.3 years ago
bcftools query -l in.vcf | while read SAMPLE; do echo -n "${SAMPLE}: " && bcftools view --samples  "${SAMPLE}" -O u in.vcf| bcftools query -f '[%DP\n]' | paste -s -d '+' | bc; done


SAMPLE1: 2423
SAMPLE2: 3838
SAMPLE3: 2659
SAMPLE4: 3249
SAMPLE5: 3812
ADD COMMENT

Login before adding your answer.

Traffic: 3297 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6