Test-drive the GATK tools and Best Practices pipelines on Terra


Check out this blog post to learn how you can get started with GATK and try out the pipelines in preconfigured workspaces (with a user-friendly interface!) without having to install anything.

JEXL and selectVariants

louwslouws San FranciscoMember

Hi,

I am trying to use select variants along with the JEXL syntax. In this specific case I wish to select those variants that

(1)Were heterozygous and,

(2)The ALT allele had more than 8 reads and,

(3) The alternate allele frequency exceeds 20% or as I define in a calculation: (ALT_AD)/(ALT_AD+REF_AD) >= 0.20

My JEXL looks like this:

(vc.getGenotype(“SAMP1").isHet()) && (vc.getGenotype(“SAMP1").getAD().1 >= 8) && (((vc.getGenotype(“SAMP1").getAD().1) / ((vc.getGenotype(“SAMP1").getAD().0)+(vc.getGenotype("SAMP1").getAD().1))) >= 0.2)

In this case, I get no results ( although GATK executes succesfully without errors). However when I exclude the third clause I get the expected results namely those variants that are (1)heterozygous and (2) The ALT_AD exceeds 8 reads. Like so:

(vc.getGenotype(“SAMP1").isHet()) && (vc.getGenotype(“SAMP1").getAD().1 >= 8)

To troubleshoot this issue I thought it might be a decimal vs. integer issue and I added a decimal value to my AD calculation, like so:

(vc.getGenotype(“SAMP1").isHet()) && (vc.getGenotype(“SAMP1").getAD().1 >= 8) && (((vc.getGenotype(“SAMP1").getAD().1 + 0.000000001) / ((vc.getGenotype(“SAMP1").getAD().0)+(vc.getGenotype("SAMP1").getAD().1))) >= 0.2)

The result was that select variants proceeded to select calls as if that third clause did not exist. I am not sure what I might be overlooking, but any insight much appreciated.

Comments

  • louwslouws San FranciscoMember

    Hi all,

    I found the issue, it is an integer vs. decimal issue. For those interested I modified my JEXL by converting the array value into a float, like so:

    (((vc.getGenotype(“SAMP1").getAD().1.floatValue()) / ((vc.getGenotype(“SAMP1").getAD().0.floatValue())+(vc.getGenotype("SAMP1").getAD().1,floatValue()))) >= 0.2)

    I hope that will help others in the future.

  • pdexheimerpdexheimer Member ✭✭✭✭

    Nice catch! That solution was not at all obvious to me

Sign In or Register to comment.