Hi, I am running the PharmCAT VCF Preprocessor 3.4.0. Here, it converts an explicitly missing genotype
(./.) to homozygous reference (0/0) while completing the ALT alleles
from pharmcat_positions.vcf.
This is related to Issue #218, but the important difference is that the input
genotype in this example is missing, not already 0/0. The direct
preprocessor output therefore changes the genotype call.
My environment has:
- PharmCAT VCF Preprocessor: 3.4.0
- bcftools: 1.22
- htslib: 1.22.1
- Python: 3.10.6
- Input: ordinary single-sample VCF, not a block gVCF
None of the following options were used:
--missing-to-ref
--absent-to-ref
--unspecified-to-ref
Here I have some example data:
The official PharmCAT 3.4.0 position catalog contains:
chr6 18130687 rs1142345 T C,G . PASS PX=TPMT GT 0/0
Example input vcf:
##fileformat=VCFv4.2
##contig=<ID=chr6>
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT SAMPLE
chr6 18130687 . T C . PASS . GT ./.
I used this command:
pharmcat_vcf_preprocessor \
-vcf missing_gt.vcf \
-o output \
-refFna GRCh38_reference.fa \
-refVcf pharmcat_positions.vcf.bgz
And my output is:
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT SAMPLE
chr6 18130687 rs1142345 T C,G . PASS PX=TPMT GT 0/0
The missing genotype ./. has become 0/0.
Completing the ALT representation from C to C,G is understandable, but it should not promote an missing genotype to a confident homozygous-reference genotype, am I correct?
The output genotype should remain missing:
GT=./.
This is also consistent with the PharmCAT VCF requirements, which state that reference (0/0) and missing (./.) positions should be specified explicitly and that missing positions should not be interpreted as
reference.
The behavior appears to originate from the catalog ALT completion block:
https://github.com/PharmGKB/PharmCAT/blob/v3.4.0/preprocessor/pcat/utilities.py#L1206-L1215
If upstream QC marks a site as missing, this can turn that unknown position into reference. The Named Allele Matcher can then use information that was not present in the input.
Could you therefore please confirm whether the ./. -> 0/0 conversion is intended?
Best,
Lukas
Hi, I am running the PharmCAT VCF Preprocessor 3.4.0. Here, it converts an explicitly missing genotype
(
./.) to homozygous reference (0/0) while completing the ALT allelesfrom
pharmcat_positions.vcf.This is related to Issue #218, but the important difference is that the input
genotype in this example is missing, not already
0/0. The directpreprocessor output therefore changes the genotype call.
My environment has:
None of the following options were used:
--missing-to-ref--absent-to-ref--unspecified-to-refHere I have some example data:
The official PharmCAT 3.4.0 position catalog contains:
Example input vcf:
I used this command:
And my output is:
The missing genotype
./.has become0/0.Completing the ALT representation from C to C,G is understandable, but it should not promote an missing genotype to a confident homozygous-reference genotype, am I correct?
The output genotype should remain missing:
GT=./.This is also consistent with the PharmCAT VCF requirements, which state that reference (
0/0) and missing (./.) positions should be specified explicitly and that missing positions should not be interpreted asreference.
The behavior appears to originate from the catalog ALT completion block:
https://github.com/PharmGKB/PharmCAT/blob/v3.4.0/preprocessor/pcat/utilities.py#L1206-L1215
If upstream QC marks a site as missing, this can turn that unknown position into reference. The Named Allele Matcher can then use information that was not present in the input.
Could you therefore please confirm whether the
./.->0/0conversion is intended?Best,
Lukas