VariantFiltration and SelectVariants

kssrkssr Posts: 6Member
edited January 2013 in Ask the GATK team

I am trying to filter variant calls which have "GQ>=20.0".

GATK SelectVariants, gives no error but gives only the header in the output file

java -Xmx2g -jar ~/GenomeAnalysisTKLite-2.1-8-gbb7f038/GenomeAnalysisTKLite.jar -R xxx -T SelectVariants --variant xxx.var.flt.vcf -o xxx.vcf -select "GQ >= 20.0"

So, I tried using VariantFiltration followed by SelectVariants. The variant filtration seems to work fine adding FT tag to the format field. And then I am trying to get records having FT tag using the following commands

java -Xmx2g -jar ~/GenomeAnalysisTKLite-2.1-8-gbb7f038/GenomeAnalysisTKLite.jar -R xxx -T VariantFiltration --variant xxx.var.flt.vcf -o xxx_filtered.vcf --genotypeFilterExpression "GQ >= 20.0" --genotypeFilterName "qual_1_filters"

java -Xmx4g -jar ~/GenomeAnalysisTKLite-2.1-8-gbb7f038/GenomeAnalysisTKLite.jar -T SelectVariants -R xxx --variant xxx_filtered.vcf -select 'vc.hasAttribute("FT")' -o xxx_qual20.vcf 

but I only get header in the output vcf file.

I am not sure if this is the right approach. Any help would be appreciated.

Post edited by Geraldine_VdAuwera on

Answers

  • Geraldine_VdAuweraGeraldine_VdAuwera Posts: 7,163Administrator, GATK Developer admin

    The simplest approach is the one in your first command; that should work. Have you looked at the article on JEXL expressions to find out what may be messing it up?

    http://gatkforums.broadinstitute.org/discussion/1255

    Geraldine Van der Auwera, PhD

  • kssrkssr Posts: 6Member

    I did check out the article on JEXL expressions.GQ is not a missing field, so all my records should have it.I don't find anything wrong with the expression "GQ >= 20.0".

  • pdexheimerpdexheimer Posts: 417Member, GSA Collaborator ✭✭✭✭

    I'm pretty sure SelectVariants can only select on variant-level annotations (in the INFO field) - GQ is specified at the sample level (FORMAT field). I think your second approach is close, but the filters added by "VariantFiltration --genotypeFilter" are sample-level (they're added to the FORMAT field), so you still can't select on them.

    To my knowledge, there's no way to do what you want in the GATK. The appropriate place would be in SelectVariants, but there would also have to be thresholding logic in there (do you select if you see it in 1 sample? 10%? 90%?). I think your best bet is probably to write a Perl filter

  • kssrkssr Posts: 6Member

    I am working on single sample data.So,I don't see a way to set thresholds. I guess the best bet would be to write a script.Thanks for the help.

  • SharonCoxSharonCox Posts: 7Member

    Dear Geraldine....I also am having the same problem. I have a .vcf file containing three samples obtained with the CombineVariants. I checked out the final .vcf and in the INFO field I have "set=variant" , I am trying to SelectVariants on the basis of this expression with the following command line:

    java -Xmx4g -jar /home/cox/GenomeAnalysisTK-2.3-9-ge5ebf34/GenomeAnalysisTK.jar -R hg19.fa -T SelectVariants --variant input.vcf -select “set=variant” -o output.vcf

    This gives no error but gives only the header in the output file .....what am I doing wrong?

  • Geraldine_VdAuweraGeraldine_VdAuwera Posts: 7,163Administrator, GATK Developer admin

    Hi Sharon,

    Your JEXL expression should be "set == 'variant'".

    Geraldine Van der Auwera, PhD

Sign In or Register to comment.