search

brentp / slivar

genetic variant expressions, annotation, and filtering for great good.
MIT License
248 stars 23 forks source link

Unable to access split-vep annotation in INFO #119

Closed MattWellie closed 2 years ago

MattWellie commented 2 years ago

Passing my VCF through Split-VEP to expose some CSQ fields, I can't work out the syntax for Slivar to see them:

error from duktape: unknown attribute:gnomAD_AF for expression:INFO.impactful & **variant**.gnomAD_AF < 0.01 & variant.FILTER == 'PASS' & variant.ALT[0] != '*'

error from duktape: unknown attribute:gnomAD_AF for expression:INFO.impactful & **INFO**.gnomAD_AF < 0.01 & variant.FILTER == 'PASS' & variant.ALT[0] != '*'

error from duktape: ReferenceError: identifier 'gnomAD_AF' undefined for expression:INFO.impactful & gnomAD_AF < 0.01 & variant.FILTER == 'PASS' & variant.ALT[0] != '*'

example variant (tons of additional CSQ removed)

chrX    154830834   rs12009271  A   G   2445.09 PASS    AC=6;AF=0.014;AN=432;AS_BaseQRankSum=1;AS_FS=0;AS_FilterStatus=PASS;AS_InbreedingCoeff=0.5457;AS_MQ=60;AS_MQRankSum=0;AS_QD=23.07;AS_ReadPosRankSum=0.9;AS_SOR=0.616;AS_VQSLOD=10.1229;AS_culprit=AS_MQRankSum;BaseQRankSum=1.82;DB;DP=5321;ExcessHet=0.0005;FS=0;InbreedingCoeff=0.5457;MLEAC=6;MLEAF=0.014;MQ=60;MQRankSum=0;POSITIVE_TRAIN_SITE;QD=23.07;QUALapprox=2550;RAW_GT_COUNT=0;ReadPosRankSum=0.938;SOR=0.616;VarDP=106;CSQ=G|downstream_gene_variant|MODIFIER|F8|ENSG00000185010|Transcript|ENST00000330287|protein_coding||||||||||rs12009271|1|4954|-1||SNV|1|HGNC|HGNC:3546||||1||CCDS44026.1|ENSP00000327895|P00451.260||UPI000006FF7B|P00451-2|1||||||0.1253|0.4437|0.0477|0|0.0026|0.0014|0.3763|0.0008365|0.02681|0.3964|0.0198|0|0|0|0.0008267|0.01055|0.0007907|0.4437|AFR|||||||||||||,;Consequence=missense_variant;gnomAD_AF=0.02681;Feature=ENST00000369529;Feature_type=Transcript;CANONICAL=YES;Gene=ENSG00000203870;MANE_SELECT=NM_001162936.4;MANE_PLUS_CLINICAL=.;impactful;genic;highest_impact_order=14

relevant header section

##INFO=<ID=CSQ,Number=.,Type=String,Description="Consequence annotations from Ensembl VEP. Format: Allele|Consequence|IMPACT|SYMBOL|Gene|Feature_type|Feature|BIOTYPE|EXON|INTRON|HGVSc|HGVSp|cDNA_position|CDS_position|Protein_position|Amino_acids|Codons|Existing_variation|ALLELE_NUM|DISTANCE|STRAND|FLAGS|VARIANT_CLASS|MINIMISED|SYMBOL_SOURCE|HGNC_ID|CANONICAL|MANE_SELECT|MANE_PLUS_CLINICAL|TSL|APPRIS|CCDS|ENSP|SWISSPROT|TREMBL|UNIPARC|UNIPROT_ISOFORM|GENE_PHENO|SIFT|PolyPhen|DOMAINS|miRNA|HGVS_OFFSET|AF|AFR_AF|AMR_AF|EAS_AF|EUR_AF|SAS_AF|AA_AF|EA_AF|gnomAD_AF|gnomAD_AFR_AF|gnomAD_AMR_AF|gnomAD_ASJ_AF|gnomAD_EAS_AF|gnomAD_FIN_AF|gnomAD_NFE_AF|gnomAD_OTH_AF|gnomAD_SAS_AF|MAX_AF|MAX_AF_POPS|CLIN_SIG|SOMATIC|PHENO|PUBMED|MOTIF_NAME|MOTIF_POS|HIGH_INF_POS|MOTIF_SCORE_CHANGE|TRANSCRIPTION_FACTORS|LoF|LoF_filter|LoF_flags|LoF_info">
##LoF=Loss-of-function annotation (HC = High Confidence; LC = Low Confidence)
##LoF_filter=Reason for LoF not being HC
##LoF_flags=Possible warning flags for LoF
##LoF_info=Info used for LoF annotation
##INFO=<ID=Consequence,Number=.,Type=String,Description="The Consequence field from INFO/CSQ">
##INFO=<ID=gnomAD_AF,Number=.,Type=Float,Description="The gnomAD_AF field from INFO/CSQ">
##INFO=<ID=Feature,Number=.,Type=String,Description="The Feature field from INFO/CSQ">
##INFO=<ID=Feature_type,Number=.,Type=String,Description="The Feature_type field from INFO/CSQ">
##INFO=<ID=CANONICAL,Number=.,Type=String,Description="The CANONICAL field from INFO/CSQ">
##INFO=<ID=Gene,Number=.,Type=String,Description="The Gene field from INFO/CSQ">
##INFO=<ID=MANE_SELECT,Number=.,Type=String,Description="The MANE_SELECT field from INFO/CSQ">
##INFO=<ID=MANE_PLUS_CLINICAL,Number=.,Type=String,Description="The MANE_PLUS_CLINICAL field from INFO/CSQ">
##bcftools_split-vepVersion=1.14+htslib-1.14
brentp commented 2 years ago

Hi, with this vcf:

##fileformat=VCFv4.2
##FILTER=<ID=PASS,Description="All filters passed">
##contig=<ID=chrX,length=249250621>
##INFO=<ID=CSQ,Number=.,Type=String,Description="Consequence annotations from Ensembl VEP. Format: Allele|Consequence|IMPACT|SYMBOL|Gene|Feature_type|Feature|BIOTYPE|EXON|INTRON|HGVSc|HGVSp|cDNA_position|CDS_position|Protein_position|Amino_acids|Codons|Existing_variation|ALLELE_NUM|DISTANCE|STRAND|FLAGS|VARIANT_CLASS|MINIMISED|SYMBOL_SOURCE|HGNC_ID|CANONICAL|MANE_SELECT|MANE_PLUS_CLINICAL|TSL|APPRIS|CCDS|ENSP|SWISSPROT|TREMBL|UNIPARC|UNIPROT_ISOFORM|GENE_PHENO|SIFT|PolyPhen|DOMAINS|miRNA|HGVS_OFFSET|AF|AFR_AF|AMR_AF|EAS_AF|EUR_AF|SAS_AF|AA_AF|EA_AF|gnomAD_AF|gnomAD_AFR_AF|gnomAD_AMR_AF|gnomAD_ASJ_AF|gnomAD_EAS_AF|gnomAD_FIN_AF|gnomAD_NFE_AF|gnomAD_OTH_AF|gnomAD_SAS_AF|MAX_AF|MAX_AF_POPS|CLIN_SIG|SOMATIC|PHENO|PUBMED|MOTIF_NAME|MOTIF_POS|HIGH_INF_POS|MOTIF_SCORE_CHANGE|TRANSCRIPTION_FACTORS|LoF|LoF_filter|LoF_flags|LoF_info">
##LoF=Loss-of-function annotation (HC = High Confidence; LC = Low Confidence)
##LoF_filter=Reason for LoF not being HC
##LoF_flags=Possible warning flags for LoF
##LoF_info=Info used for LoF annotation
##INFO=<ID=Consequence,Number=.,Type=String,Description="The Consequence field from INFO/CSQ">
##INFO=<ID=gnomAD_AF,Number=.,Type=Float,Description="The gnomAD_AF field from INFO/CSQ">
##INFO=<ID=Feature,Number=.,Type=String,Description="The Feature field from INFO/CSQ">
##INFO=<ID=Feature_type,Number=.,Type=String,Description="The Feature_type field from INFO/CSQ">
##INFO=<ID=CANONICAL,Number=.,Type=String,Description="The CANONICAL field from INFO/CSQ">
##INFO=<ID=Gene,Number=.,Type=String,Description="The Gene field from INFO/CSQ">
##INFO=<ID=MANE_SELECT,Number=.,Type=String,Description="The MANE_SELECT field from INFO/CSQ">
##INFO=<ID=MANE_PLUS_CLINICAL,Number=.,Type=String,Description="The MANE_PLUS_CLINICAL field from INFO/CSQ">
#CHROM  POS ID  REF ALT QUAL    FILTER  INFO
chrX    154830834   rs12009271  A   G   2445.09 PASS    CSQ=G|downstream_gene_variant|MODIFIER|F8|ENSG00000185010|Transcript|ENST00000330287|protein_coding||||||||||rs12009271|1|4954|-1||SNV|1|HGNC|HGNC:3546||||1||CCDS44026.1|ENSP00000327895|P00451.260||UPI000006FF7B|P00451-2|1||||||0.1253|0.4437|0.0477|0|0.0026|0.0014|0.3763|0.0008365|0.02681|0.3964|0.0198|0|0|0|0.0008267|0.01055|0.0007907|0.4437|AFR|||||||||||||;Consequence=missense_variant;gnomAD_AF=0.02681;Feature=ENST00000369529;Feature_type=Transcript;CANONICAL=YES;Gene=ENSG00000203870;MANE_SELECT=NM_001162936.4;MANE_PLUS_CLINICAL=.;impactful;genic;highest_impact_order=14

and this command:

slivar expr --info "INFO.gnomAD_AF > 0.01" -v x.vcf

I see no error. Can you share a full, valid VCF that has the error?

MattWellie commented 2 years ago

@brentp thanks for the reply! I think this syntax was working for me, but delivering no variants... so I kept tweaking the syntax into something incorrect. I've spotted what looks like a problem elsewhere in the pipeline, so Slivar wasn't the culprit. Facepalms all round!

brentp commented 2 years ago

Glad you figured it out!