search

grenaud / glactools

command-line tools for the management of genotype likelihoods and allele counts
http://grenaud.github.io/glactools/
GNU General Public License v3.0
29 stars 2 forks source link

vcfm2acf issue #26

Closed jaesoonpark closed 2 months ago

jaesoonpark commented 2 years ago

Hello.

I have a problem using vcfm2acf.

glactools vcfm2acf --onlyGT --fai human_g1k_v37.fasta.fai chr1.recode.vcf > chr1.acf.gz

error message libgab.h base2int() Invalid base *

Please let me know if any improvements can be made.

Thanks.

grenaud commented 2 years ago

Could you check the ALT field? ? do you have a * in there?

jaesoonpark commented 2 years ago

Thanks for your comment. I have solved the aforementioned issues.

However, I had another issue in using EPO file.

The command line was as followed: glactools/glactools vcfm2acf --onlyGT --epo all.epo.gz --fai human_g1k_v37.fasta.fai test > test.acf.gz

Then the error message was as followed: Error, missing data in the EPO file

I wonder if you have any clue what this error message means?

Thanks for your help in advance.

grenaud commented 2 years ago

Is it possible that there is a different encoding of chromosomes? "1" vs "chr1"

jaesoonpark commented 2 years ago

Great! Removing 'chr' from my VCF file helps and ACF file seems to be made partially. (Previously, ACF file was not made at all.) However, the error message pops out after a while. Could you suggest any other reason for this message?

Thanks in advance.

grenaud commented 2 years ago

just to be sure, are you sure your variants were mapped to hg19/v37? Are there chromosomes that are not 1,2,3..22? did you remove any _random chr? can you do a cut -f 1 |uniq -c on column 1?

Priambodobagus commented 1 year ago

Hello I got an issue I run vcfm2acf with a merged of multiple vcf files generated by haplotype caller in GATK. and my reference genome was at the scaffold level.

my command: ./glactools vcfm2acf --onlyGT --fai ReferenceGenome.fasta.fai my.vcf > my.acf

but I got this error: the PL does not have 3,6, or 10 fields, has 1 fields

Does anyone know how to fix this problem?

Best, BAGUS

grenaud commented 1 year ago

Can you please post a few of the lines that are causing the problem?

M-Asim-Javed commented 3 months ago

Hello, I encountered the same error while running the below command

glactools vcfm2acf --onlyGT --fai pb3A_genome.fa.fai pb_filtered_biallelic_snps.vcf > pb_filterd_bialleic_snps.acf.gz

Only this below Error and no other error

SimpleVCF: for line = chromosome1 181 0 121 the PL field does not have 3,6 or 10 fields, has 2 fields

Here're a few lines from my vcf file!

#CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  ab_01   ab_11   ab_15   ab_21   ab_5    qc_01   qc_40   qc_6    qc_7    sa_2    sa_21   sa_9    th20_10 th20_5  th76_14 th76_2  th76_6  th76_8
chromosome1 181 .   A   C   33.12   PASS    AC=1;AF=0.083;AN=12;DP=178;FS=0;MLEAC=1;MLEAF=0.083;MQ=38.24;QD=33.12;SOR=1.609 GT:AD:DP:GQ:PL  0:8,0:8:99:0,121    0:2,0:2:45:0,45 0:11,0:11:45:0,45   0:6,0:6:90:0,90 0:4,0:4:99:0,133    .:32,0:32:.:0,0 .:22,0:22:.:0,0 .:50,0:50:.:0,0 1:0,1:1:45:45,0 0:5,0:5:42:0,42 .:0,0:0:.:0,0   .:4,0:4:.:0,0   0:10,0:10:99:0,296  0:3,0:3:99:0,125    0:5,0:5:90:0,90 0:3,0:3:99:0,99 .:0,0:0:.:0,0   0:10,0:10:45:0,45
chromosome1 890 .   C   A   82.38   PASS    AC=3;AF=0.273;AN=11;DP=187;FS=0;MLEAC=4;MLEAF=0.364;MQ=35.41;QD=27.46;SOR=2.833 GT:AD:DP:GQ:PL  .:20,0:20:.:0,0 1:0,1:1:38:38,0 .:27,0:27:.:0,0 1:0,1:1:32:32,0 0:4,0:4:99:0,135    0:50,0:50:99:0,350  0:3,0:3:99:0,134    0:4,0:4:99:0,135    0:3,0:3:99:0,118    .:6,0:6:.:0,0   .:12,0:12:.:0,0 .:12,0:12:.:0,0 0:3,0:3:99:0,125    0:3,0:3:99:0,125    1:0,1:1:31:31,0 0:3,0:3:99:0,99 .:15,0:15:.:0,0 .:19,0:19:.:0,0
chromosome1 9838    .   C   A   25606.5 PASS    AC=10;AF=0.556;AN=18;BaseQRankSum=-0.807;DP=982;FS=0;MLEAC=10;MLEAF=0.556;MQ=31.05;MQRankSum=1.51;QD=28.23;ReadPosRankSum=0.523;SOR=0.67    GT:AD:DP:GQ:PL  1:4,89:93:99:2480,0 1:3,64:67:99:1825,0 1:5,151:156:99:4239,0   1:0,107:107:99:3209,0   0:3,0:3:99:0,113    0:7,0:7:99:0,270    0:8,0:8:99:0,358    0:35,0:35:99:0,1478 0:12,0:12:99:0,509  1:0,36:36:99:1087,0 1:1,59:60:99:1651,0 1:0,45:45:99:1285,0 0:3,0:3:99:0,113    0:4,0:4:99:0,135    1:0,124:124:99:3670,0   0:3,0:3:99:0,125    1:1,81:82:99:2302,0 1:1,136:137:99:3884,0
chromosome1 9943    .   G   T   3748.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=0;DP=280;FS=44.722;MLEAC=4;MLEAF=0.222;MQ=31.37;MQRankSum=-0.942;QD=26.97;ReadPosRankSum=-0.796;SOR=3.445  GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:19,27:46:99:540,0 1:10,15:25:99:225,0 1:1,15:16:99:630,0  1:0,52:52:99:2385,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9968    .   G   A   3729.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.251;DP=318;FS=15.179;MLEAC=4;MLEAF=0.222;MQ=31.1;MQRankSum=0.202;QD=21.07;ReadPosRankSum=0.773;SOR=1.877    GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:28,35:63:99:360,0 1:15,21:36:99:270,0 1:2,16:18:99:630,0  1:2,58:60:99:2501,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9970    .   C   T   3819.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=0.306;DP=318;FS=15.179;MLEAC=4;MLEAF=0.222;MQ=31.1;MQRankSum=0.205;QD=21.58;ReadPosRankSum=0.773;SOR=1.877 GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:28,36:64:99:405,0 1:15,21:36:99:315,0 1:2,16:18:99:630,0  1:2,57:59:99:2501,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9979    .   A   C   7708.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.33;DP=324;FS=1.678;MLEAC=4;MLEAF=0.222;MQ=31.12;MQRankSum=0.677;QD=29.56;ReadPosRankSum=1.31;SOR=0.491   GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:1,67:68:99:3060,0 1:1,38:39:99:1620,0 1:3,16:19:99:585,0  1:2,55:57:99:2475,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9980    .   C   G   7708.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.289;DP=326;FS=1.653;MLEAC=4;MLEAF=0.222;MQ=31.06;MQRankSum=0.668;QD=30.62;ReadPosRankSum=1.3;SOR=0.365  GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:1,68:69:99:3060,0 1:2,38:40:99:1620,0 1:3,16:19:99:585,0  1:2,55:57:99:2475,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9983    .   G   C   4063.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.79;DP=325;FS=8.518;MLEAC=4;MLEAF=0.222;MQ=31.09;MQRankSum=0.165;QD=22.09;ReadPosRankSum=0.951;SOR=1.382  GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:28,40:68:99:675,0 1:17,23:40:99:360,0 1:3,16:19:99:585,0  1:2,55:57:99:2475,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9989    .   T   C   7978.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.48;DP=333;FS=3.579;MLEAC=4;MLEAF=0.222;MQ=30.9;MQRankSum=0.659;QD=28.17;ReadPosRankSum=0.99;SOR=0.254   GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:1,69:70:99:3105,0 1:2,40:42:99:1755,0 1:3,17:20:99:630,0  1:2,58:60:99:2520,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9993    .   G   C   8293.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.555;DP=331;FS=7.64;MLEAC=4;MLEAF=0.222;MQ=30.96;MQRankSum=0.896;QD=26.8;ReadPosRankSum=1.88;SOR=0.146   GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:2,67:69:99:3195,0 1:2,41:43:99:1890,0 1:3,17:20:99:630,0  1:3,55:58:99:2610,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9995    .   C   G   8158.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.542;DP=340;FS=7.72;MLEAC=4;MLEAF=0.222;MQ=30.89;MQRankSum=0.861;QD=26;ReadPosRankSum=1.57;SOR=0.128 GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:2,70:72:99:3105,0 1:2,43:45:99:1890,0 1:3,17:20:99:630,0  1:3,59:62:99:2565,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
chromosome1 9999    .   A   T   5688.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.42;DP=348;FS=0;MLEAC=4;MLEAF=0.222;MQ=30.86;MQRankSum=0.164;QD=27.48;ReadPosRankSum=-0.108;SOR=0.672 GT:AD:DP:GQ:PL  0:6,0:6:99:0,135    0:3,0:3:99:0,110    0:6,0:6:99:0,135    0:4,0:4:99:0,135    0:3,0:3:99:0,113    1:29,46:75:99:1154,0    1:17,29:46:99:900,0 1:3,17:20:99:630,0  1:3,63:66:99:3036,0 0:10,0:10:99:0,101  0:3,0:3:99:0,125    0:8,0:8:99:0,239    0:3,0:3:99:0,113    0:4,0:4:99:0,135    0:5,0:5:99:0,125    0:3,0:3:99:0,125    0:5,0:5:99:0,135    0:3,0:3:99:0,99
grenaud commented 3 months ago

Yes, normally you get 3 values for PL, homo ref, het, homo alt. How did you get "0,45" for example?

M-Asim-Javed commented 3 months ago

Thank you for the immediate response.

I used the GATK pipeline to call variants and combined gvcf files, then selected SNPs only. The species under study is haploid. I don't know how to resolve the issue. I also have an example VCF file from a group working on population genomics studies, and their VCF file also shows two PL values. Could you please guide me on how to resolve this issue or suggest any other way to handle this PL error?

Here's the example VCF file...

#CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  032i    053i    117 12.1.014    12.1.037    12.1.078    12.1.181    12.1.205    12.1.225    12.1.234    12.1.241    12.1.291    127 169 204 205 37  B71 BR32    BTBa-B1 BTBa-B2 BTGP1-b BTGS6e  BTGS6f  BTGS6g  BTGS6h  BTJP4-1 BTJP4-11    BTJP4-12    BTJP4-15    BTJP4-16    BTJP4-18    BTJP4-2 BTJP4-3 BTJP4-5 BTJP4-6 BTJP4-9 BTMP-S13-1  BTMP-S13-2  BdBar16-1   Br116.5 Br118.2 Br48    Br7 Br80    P29 PY0925  PY25.1  PY35.3  PY36.1  PY5003  PY5033  PY5035  PY6017  PY6025  PY6045  PY6047  T25 ZMW18-06    ZMW18-08    ZMW18-10    ZMW18-11    ZMW19-09    ZMW19-13    ZMW19-17    ZMW20-03    ZMW20-04    ZMW20-07    ZMW20-14    ZMW20-15    ZMW20-16
1   168 .   A   G   268.53  PASS    AC=1;AF=0.333;AN=2;DP=18;FS=0;MLEAC=8;MLEAF=1;MQ=20.21;QD=26.85;SOR=3.258   GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,10:10:99:249,0  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   195 .   C   T   819.72  PASS    AC=1;AF=0.143;AN=6;DP=26;FS=0;MLEAC=6;MLEAF=0.857;MQ=22.63;QD=25.36;SOR=1.863   GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:27:0,27 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:16:0,16 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,18:18:99:810,0  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:21:0,21 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:41:0,41 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   205 .   C   T   906.3   PASS    AC=1;AF=0.111;AN=8;DP=29;FS=0;MLEAC=5;MLEAF=0.556;MQ=22.48;QD=28.73;SOR=1.402   GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:42:0,42 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:16:0,16 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,20:20:99:900,0  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:38:0,38 .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:42:0,42 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:32:0,32 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   258 .   C   T   1176.77 PASS    AC=1;AF=0.111;AN=9;DP=39;FS=0;MLEAC=5;MLEAF=0.556;MQ=21.93;QD=30.62;SOR=0.693   GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:32:0,32 .:0,0:0:.:0,0   0:3,0:3:76:0,76 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,26:26:99:1170,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:6:0,6   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:45:0,45 .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:23:0,23 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   267 .   G   A   1223.07 PASS    AC=1;AF=0.125;AN=8;DP=41;FS=0;MLEAC=6;MLEAF=0.75;MQ=21.86;QD=28.17;SOR=0.765    GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:16:0,16 .:0,0:0:.:0,0   0:3,0:3:80:0,80 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,27:27:99:1215,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:40:0,40 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   0:1,0:1:40:0,40 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:49:0,49 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:11:0,11 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   268 .   G   A   1224.79 PASS    AC=1;AF=0.143;AN=7;DP=41;FS=0;MLEAC=6;MLEAF=0.857;MQ=21.86;QD=26.8;SOR=0.765    GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   0:3,0:3:55:0,55 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,27:27:99:1215,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:49:0,49 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:12:0,12 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   306 .   C   T   1207.81 PASS    AC=2;AF=0.182;AN=10;DP=48;FS=0;MLEAC=9;MLEAF=0.818;MQ=21.86;QD=26;SOR=0.765 GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   0:2,0:2:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:45:0,45 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,26:26:99:1150,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:38:0,38 .:0,0:0:.:0,0   0:1,0:1:41:0,41 .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 0:1,0:1:21:0,21 .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:49:0,49 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:42:0,42 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,1:1:45:45,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   318 .   G   A   1180.54 PASS    AC=2;AF=0.167;AN=11;DP=47;FS=0;MLEAC=9;MLEAF=0.75;MQ=21.93;QD=30.02;SOR=0.693   GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:42:0,42 .:0,0:0:.:0,0   0:2,0:2:74:0,74 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:45:0,45 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,25:25:99:1125,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:2,0:2:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:19:0,19 .:0,0:0:.:0,0   0:1,0:1:32:0,32 .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 0:1,0:1:21:0,21 .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:45:0,45 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:22:0,22 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,1:1:45:45,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
1   321 .   C   T   1134.69 PASS    AC=2;AF=0.154;AN=11;DP=47;FS=0;MLEAC=9;MLEAF=0.692;MQ=22.1;MQRankSum=-1.939;QD=31.98;SOR=0.223  GT:AD:DP:GQ:PL  .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:16:0,16 .:0,0:0:.:0,0   0:2,0:2:62:0,62 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:1,25:26:99:1080,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:2,0:2:79:0,79 .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   0:2,0:2:2:0,2   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:39:0,39 .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 0:1,0:1:21:0,21 .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:43:0,43 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:1,0:1:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   0:1,0:1:35:0,35 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   1:0,1:1:45:45,0 .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0   .:0,0:0:.:0,0
grenaud commented 3 months ago

is this haploid data?

M-Asim-Javed commented 3 months ago

Yes, the above dataset is haploid, and my data is also haploid.

grenaud commented 2 months ago

I see! Normally I would modify the code to handle this but can you either code it yourself or ask chatgpt/claude to transform the GT from 0 = 0/0 1 = 1/1

Now the tricky part is getting the PL fields, I would transform PL1,PL2 to PL1,0,PL2 this should make the program unable to distinguish between homo/hetero and keep a single dominant allele, can you test it on a few lines and check if the output is satisfactory?

M-Asim-Javed commented 2 months ago

Thank you for getting back to me.

I transformed the VCF file with both GT and PL fields. Now, the VCF file looks like this below. I am not sure if I transformed it correctly, but it worked with Glactools using the vcfm2acf and then acf2betascan. However, the values I received seem odd, especially the haploid genome count, which should be 18 in total, but Glactools gave me 36, double of each.

I ran a Python script to count the allele frequency at each SNP site as generated by Glactools to run betascan, and it gave me 18 haploid counts for each SNP site and also generated half of each allele count. I don't know which one is correct; it seems like Glactools is giving me double of each.

I would appreciate your help and assistance.

#CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  ab_01   ab_11   ab_15   ab_21   ab_5    qc_01   qc_40   qc_6    qc_7    sa_2    sa_21   sa_9    th20_10 th20_5  th76_14 th76_2  th76_6  th76_8
chromosome1 181 .   A   C   33.12   PASS    AC=1;AF=0.083;AN=12;DP=178;FS=0;MLEAC=1;MLEAF=0.083;MQ=38.24;QD=33.12;SOR=1.609 GT:AD:DP:GQ:PL  0/0:8,0:8:99:0,0,121    0/0:2,0:2:45:0,45,. 0/0:11,0:11:45:0,45,.   0/0:6,0:6:90:0,90,. 0/0:4,0:4:99:0,133,.    ./.:32,0:32:.:0,0,. ./.:22,0:22:.:0,0,. ./.:50,0:50:.:0,0,. 1/1:0,1:1:45:45,0,. 0/0:5,0:5:42:0,42,. ./.:0,0:0:.:0,0,.   ./.:4,0:4:.:0,0,.   0/0:10,0:10:99:0,296,.  0/0:3,0:3:99:0,125,.    0/0:5,0:5:90:0,90,. 0/0:3,0:3:99:0,99,. ./.:0,0:0:.:0,0,.   0/0:10,0:10:45:0,45,.
chromosome1 890 .   C   A   82.38   PASS    AC=3;AF=0.273;AN=11;DP=187;FS=0;MLEAC=4;MLEAF=0.364;MQ=35.41;QD=27.46;SOR=2.833 GT:AD:DP:GQ:PL  ./.:20,0:20:.:0,0,0 1/1:0,1:1:38:38,0,. ./.:27,0:27:.:0,0,. 1/1:0,1:1:32:32,0,. 0/0:4,0:4:99:0,135,.    0/0:50,0:50:99:0,350,.  0/0:3,0:3:99:0,134,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,118,.    ./.:6,0:6:.:0,0,.   ./.:12,0:12:.:0,0,. ./.:12,0:12:.:0,0,. 0/0:3,0:3:99:0,125,.    0/0:3,0:3:99:0,125,.    1/1:0,1:1:31:31,0,. 0/0:3,0:3:99:0,99,. ./.:15,0:15:.:0,0,. ./.:19,0:19:.:0,0,.
chromosome1 9838    .   C   A   25606.5 PASS    AC=10;AF=0.556;AN=18;BaseQRankSum=-0.807;DP=982;FS=0;MLEAC=10;MLEAF=0.556;MQ=31.05;MQRankSum=1.51;QD=28.23;ReadPosRankSum=0.523;SOR=0.67    GT:AD:DP:GQ:PL  1/1:4,89:93:99:2480,0,0 1/1:3,64:67:99:1825,0,. 1/1:5,151:156:99:4239,0,.   1/1:0,107:107:99:3209,0,.   0/0:3,0:3:99:0,113,.    0/0:7,0:7:99:0,270,.    0/0:8,0:8:99:0,358,.    0/0:35,0:35:99:0,1478,. 0/0:12,0:12:99:0,509,.  1/1:0,36:36:99:1087,0,. 1/1:1,59:60:99:1651,0,. 1/1:0,45:45:99:1285,0,. 0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    1/1:0,124:124:99:3670,0,.   0/0:3,0:3:99:0,125,.    1/1:1,81:82:99:2302,0,. 1/1:1,136:137:99:3884,0,.
chromosome1 9943    .   G   T   3748.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=0;DP=280;FS=44.722;MLEAC=4;MLEAF=0.222;MQ=31.37;MQRankSum=-0.942;QD=26.97;ReadPosRankSum=-0.796;SOR=3.445  GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:19,27:46:99:540,0,. 1/1:10,15:25:99:225,0,. 1/1:1,15:16:99:630,0,.  1/1:0,52:52:99:2385,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9968    .   G   A   3729.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.251;DP=318;FS=15.179;MLEAC=4;MLEAF=0.222;MQ=31.1;MQRankSum=0.202;QD=21.07;ReadPosRankSum=0.773;SOR=1.877    GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:28,35:63:99:360,0,. 1/1:15,21:36:99:270,0,. 1/1:2,16:18:99:630,0,.  1/1:2,58:60:99:2501,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9970    .   C   T   3819.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=0.306;DP=318;FS=15.179;MLEAC=4;MLEAF=0.222;MQ=31.1;MQRankSum=0.205;QD=21.58;ReadPosRankSum=0.773;SOR=1.877 GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:28,36:64:99:405,0,. 1/1:15,21:36:99:315,0,. 1/1:2,16:18:99:630,0,.  1/1:2,57:59:99:2501,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9979    .   A   C   7708.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.33;DP=324;FS=1.678;MLEAC=4;MLEAF=0.222;MQ=31.12;MQRankSum=0.677;QD=29.56;ReadPosRankSum=1.31;SOR=0.491   GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:1,67:68:99:3060,0,. 1/1:1,38:39:99:1620,0,. 1/1:3,16:19:99:585,0,.  1/1:2,55:57:99:2475,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9980    .   C   G   7708.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.289;DP=326;FS=1.653;MLEAC=4;MLEAF=0.222;MQ=31.06;MQRankSum=0.668;QD=30.62;ReadPosRankSum=1.3;SOR=0.365  GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:1,68:69:99:3060,0,. 1/1:2,38:40:99:1620,0,. 1/1:3,16:19:99:585,0,.  1/1:2,55:57:99:2475,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9983    .   G   C   4063.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.79;DP=325;FS=8.518;MLEAC=4;MLEAF=0.222;MQ=31.09;MQRankSum=0.165;QD=22.09;ReadPosRankSum=0.951;SOR=1.382  GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:28,40:68:99:675,0,. 1/1:17,23:40:99:360,0,. 1/1:3,16:19:99:585,0,.  1/1:2,55:57:99:2475,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9989    .   T   C   7978.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.48;DP=333;FS=3.579;MLEAC=4;MLEAF=0.222;MQ=30.9;MQRankSum=0.659;QD=28.17;ReadPosRankSum=0.99;SOR=0.254   GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:1,69:70:99:3105,0,. 1/1:2,40:42:99:1755,0,. 1/1:3,17:20:99:630,0,.  1/1:2,58:60:99:2520,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9993    .   G   C   8293.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.555;DP=331;FS=7.64;MLEAC=4;MLEAF=0.222;MQ=30.96;MQRankSum=0.896;QD=26.8;ReadPosRankSum=1.88;SOR=0.146   GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:2,67:69:99:3195,0,. 1/1:2,41:43:99:1890,0,. 1/1:3,17:20:99:630,0,.  1/1:3,55:58:99:2610,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9995    .   C   G   8158.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=-0.542;DP=340;FS=7.72;MLEAC=4;MLEAF=0.222;MQ=30.89;MQRankSum=0.861;QD=26;ReadPosRankSum=1.57;SOR=0.128 GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:2,70:72:99:3105,0,. 1/1:2,43:45:99:1890,0,. 1/1:3,17:20:99:630,0,.  1/1:3,59:62:99:2565,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.
chromosome1 9999    .   A   T   5688.88 PASS    AC=4;AF=0.222;AN=18;BaseQRankSum=1.42;DP=348;FS=0;MLEAC=4;MLEAF=0.222;MQ=30.86;MQRankSum=0.164;QD=27.48;ReadPosRankSum=-0.108;SOR=0.672 GT:AD:DP:GQ:PL  0/0:6,0:6:99:0,0,135    0/0:3,0:3:99:0,110,.    0/0:6,0:6:99:0,135,.    0/0:4,0:4:99:0,135,.    0/0:3,0:3:99:0,113,.    1/1:29,46:75:99:1154,0,.    1/1:17,29:46:99:900,0,. 1/1:3,17:20:99:630,0,.  1/1:3,63:66:99:3036,0,. 0/0:10,0:10:99:0,101,.  0/0:3,0:3:99:0,125,.    0/0:8,0:8:99:0,239,.    0/0:3,0:3:99:0,113,.    0/0:4,0:4:99:0,135,.    0/0:5,0:5:99:0,125,.    0/0:3,0:3:99:0,125,.    0/0:5,0:5:99:0,135,.    0/0:3,0:3:99:0,99,.

Here's python script i used to calculate allele frquency on chromosome20:

import pysam
import numpy as np

def parse_vcf(vcf_file, output_file, contig, start_pos, end_pos, haploid_genomes):
    # Open the VCF file
    vcf_in = pysam.VariantFile(vcf_file)

    # Open the output file
    with open(output_file, 'w') as out_f:
        for record in vcf_in.fetch(contig, start_pos, end_pos + 1):
            pos = record.pos
            allele_counts = np.zeros(2)

            for sample in record.samples:
                genotype = record.samples[sample]['GT']
                if genotype is not None:
                    allele_counts[0] += genotype.count(0)
                    allele_counts[1] += genotype.count(1)

            # Calculate minor allele frequency (MAF)
            minor_allele_count = min(allele_counts)
            minor_allele_freq = minor_allele_count / haploid_genomes

            # Write the data to the output file
            out_f.write(f"{pos}\t{int(minor_allele_count)}\t{haploid_genomes}\n")

if __name__ == "__main__":
    vcf_file = "pb_snps_for_genotype.vcf.gz"  # Replace with your compressed VCF file path
    output_file = "output_chr20.txt"  # Replace with your desired output file path
    contig = "chromosome20"  # Use the correct contig name
    start_pos = 311
    end_pos = 625161
    haploid_genomes = 18

    # Run the function to parse the VCF file
    parse_vcf(vcf_file, output_file, contig, start_pos, end_pos, haploid_genomes)
grenaud commented 2 months ago

Dear Mohammed, a quick question before, isn't BetaScan configured for balancing selection and considers heterozygous advantage?

M-Asim-Javed commented 2 months ago

Hello grenaud,

Yes, you are right about balancing selection. I mixed the Python script with the Glactools output. I ran Betascan with the Glactools output, and it worked fine. I was only concerned about the third column with the double individual count. Thank you for your help and patience.

grenaud commented 2 months ago

sure but I would question what do the results even mean for haploid data. I would suggest asking Katie Siewert-Rocks directly.