<p>I am following the guidelines given in this links for variant selection for some specific cases that are not as much detailed in the GATK site.</p>
<ul>
<li><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants</a></li>
<li><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions</a></li>
</ul>
<p>I just had to do some filtering of my own data and tested several methods (to see what the outcome would be). I am sharing those here for others and myself, if it is ever needed again.</p>
<h2><strong>Tutorials .... !!!</strong></h2>
<p><strong>- method to call the sites that are all homozygous Variants (1/1, 2/2 genotype)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() == 6'
# outcome: no sites will have ambigous (./.) and/or homRef (0/0). All other homVar are allowed
</code></pre>
<p><strong>- in the above script if we want the site where at least 4 samples are homozygous ref or variants</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() &gt;= 4'
# outcome: atleast 4 samples are homVar, rest can be any genotype (./., 0/0, 0/1, etc.)
</code></pre>
<p><strong>- at least 1 sample is homVar (1/1, 2/2 or 3/3)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() != 0'
</code></pre>
<p><strong>- all the samples are homRef (0/0 genotypes)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomRefCount() == 6'
</code></pre>
<h2>Note: 6 chr = 3 samples (when samples are diploid)</h2>
<p><strong>- sites that are all HomVar or nocalls (./.)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() == vc.getCalledChrCount()/2'
# outcome: no 0/0, 0/1
# Note: samples with GT = ./. are not counted as called chr
</code></pre>
<h2>Useful after samples are sub-sampled from multisample vcfs</h2>
<p><strong>- sites where no genotypes were called for all the samples in the selected vcf</strong>
<strong>-  or say, select vcf sites that are not called in any sample (i.e don't have any genotypes, only contain ./.)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() == 0'
# outcome = vcf sites (lines) that are all no calls (GT = ./.) on all samples are selected
</code></pre>
<p><strong>- samples are either homRef (0/0) or no call (./.)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomRefCount() == vc.getCalledChrCount()/2'
# outcome: vcf sites where all the samples are GT = ./. (no call) or GT = 0/0
</code></pre>
<p><strong>- site where 6 chromosomes (exact) or 3 samples were called for diploid genome</strong></p>
<pre><code>java -jar /home/everestial007/GenomeAnalysisTK-3.7/GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() == 6'
# Note of Caution: 6 chr = 3 samples (when samples are diploid)
</code></pre>
<p><strong>- site where at least 3 samples (diploid) are called</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() &gt;=6'
</code></pre>
<p><strong>- if I have 6-sample vcf</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() ==12'
# 12 haploid chr = 6 samples (diploid)
# outcome: sites that were called in all samples, there will be no GT = ./.
</code></pre>
<p><strong>- select vcf sites where no sample has hetVar (no GT = 0/1, 0/2), only GT = ./., 0/0, 1/1, 2/2</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHetCount() == 0'
</code></pre>
<p><strong>- select vcf sites where no sample has hetVar (no Gt = 0/1, 0/2), and at least one sample is HomVar (GT = 1/1)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHetCount() == 0 &amp;&amp; vc.getHomVarCount() &gt;= 1'
</code></pre>
<p><strong>- site where at least 1 sample is homVar (GT = 1/1 and/or 2/2 etc.)</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() &gt;= 1'
</code></pre>
<p><strong>- where a site is only homRef for a particular sample and is only SNP</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getGenotype(&quot;ms01e&quot;).isHomRef()' -O test_variants.ms01e.Chr_2_3.vcf
</code></pre>
<p><strong>- where a site is only homVar for a particular sample and is only SNP</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getGenotype(&quot;ms01e&quot;).isHomVar()' -O test_variants.ms01e.Chr_2_3.vcf
</code></pre>
<p><strong>- where a site is not homRef for a particular sample and is only SNP. This will also select empty sites</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select '! vc.getGenotype(&quot;ms01e&quot;).isHomRef()' -O test_variants.ms01e.Chr_2_3.vcf
</code></pre>
<p><strong>NOTE: following genotype call methods won't work:</strong></p>
<pre><code># to get only heterozygous sites
-select 'vc.getGenotype(&quot;ms01e&quot;).isHetVar() 

# to get any variant sites
-select 'vc.getGenotype(&quot;ms01e&quot;).isVar() 

# won't work if getGenotype method is combined
-select 'vc.getGenotype(&quot;ms01e&quot;).isHomVar() || vc.getGenotype(&quot;ms01e&quot;.isHomRef()'  
</code></pre>
<p><strong>So, I figured out the following methods where count helps with selecting particular site</strong></p>
<p><strong>- where a site is a heterozygous variant site for a particular sample</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select 'vc.getHetCount() == 1'  -O test_variants.ms01e.vcf

NOTE: 'vc.getGenotype(&quot;ms01e&quot;).isHetVar()' wont work so generated the above expression.
</code></pre>
<p><strong>- where a site is a variant site (either homVar or hetVar) for a particular sample and is only SNP</strong></p>
<pre><code>java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getHetCount() == 1 || vc.getHomVarCount() == 1'  -O test_variants.ms01e.vcf

NOTE: 'vc.getGenotype(&quot;NA12878&quot;).isHomRef()' || vc.getGenotype(&quot;NA12878&quot;).isHomVar()' wont work together, so generated this above expression.
</code></pre>
<p>Other updates to follow.
Thanks,</p>

CC BY-SA 4.0
2022-02-26T03:49:42.750Z
Restored because this is an answer to the question [i.e. a user has answered their own question].
gringer
CC BY-SA 4.0
2022-02-25T00:30:56.340Z
2022-02-28T14:23:10.657Z
2022-02-28T14:23:10.657Z
everestial
<p>How to do variants selection in some corner cases using GATK and JEXL expressions?</p>
<p>I am following the guidelines given in this links for variant selection for some specific cases are not as much detailed.</p>
<ul>
<li><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants</a></li>
<li><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions</a></li>
</ul>
<p>I just had to do some filtering of my own data and tested several methods (to see what the outcome would be). I am sharing those here for others and myself, if it is ever needed again.</p>

<variant-calling><gatk><variant-filtration>
GATK variant filtration using "SelectVariants" and use of JEXL queries
# GATK variant filtration using "SelectVariants" and use of JEXL queries

(Asked by: username_0 on 2022-02-25)

How to do variants selection in some corner cases using GATK and JEXL expressions?

I am following the guidelines given in this links for variant selection for some specific cases are not as much detailed.

* <https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants>
* <https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions>

I just had to do some filtering of my own data and tested several methods (to see what the outcome would be). I am sharing those here for others and myself, if it is ever needed again.

**Tags:** variant-calling, gatk, variant-filtration

---

## 1 answers

### Answer by username_0 on 2022-02-25. Score: 1

I am following the guidelines given in this links for variant selection for some specific cases that are not as much detailed in the GATK site.

https://gatk.broadinstitute.org/hc/en-us/articles/360036430432-SelectVariants
https://gatk.broadinstitute.org/hc/en-us/articles/360035891011-JEXL-filtering-expressions

I just had to do some filtering of my own data and tested several methods (to see what the outcome would be). I am sharing those here for others and myself, if it is ever needed again.
Tutorials .... !!!
- method to call the sites that are all homozygous Variants (1/1, 2/2 genotype)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() == 6'
# outcome: no sites will have ambigous (./.) and/or homRef (0/0). All other homVar are allowed

- in the above script if we want the site where at least 4 samples are homozygous ref or variants
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() >= 4'
# outcome: atleast 4 samples are homVar, rest can be any genotype (./., 0/0, 0/1, etc.)

- at least 1 sample is homVar (1/1, 2/2 or 3/3)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() != 0'

- all the samples are homRef (0/0 genotypes)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomRefCount() == 6'

Note: 6 chr = 3 samples (when samples are diploid)
- sites that are all HomVar or nocalls (./.)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() == vc.getCalledChrCount()/2'
# outcome: no 0/0, 0/1
# Note: samples with GT = ./. are not counted as called chr

Useful after samples are sub-sampled from multisample vcfs
- sites where no genotypes were called for all the samples in the selected vcf
-  or say, select vcf sites that are not called in any sample (i.e don't have any genotypes, only contain ./.)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() == 0'
# outcome = vcf sites (lines) that are all no calls (GT = ./.) on all samples are selected

- samples are either homRef (0/0) or no call (./.)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomRefCount() == vc.getCalledChrCount()/2'
# outcome: vcf sites where all the samples are GT = ./. (no call) or GT = 0/0

- site where 6 chromosomes (exact) or 3 samples were called for diploid genome
java -jar /home/everestial007/GenomeAnalysisTK-3.7/GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() == 6'
# Note of Caution: 6 chr = 3 samples (when samples are diploid)

- site where at least 3 samples (diploid) are called
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() >=6'

- if I have 6-sample vcf
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getCalledChrCount() ==12'
# 12 haploid chr = 6 samples (diploid)
# outcome: sites that were called in all samples, there will be no GT = ./.

- select vcf sites where no sample has hetVar (no GT = 0/1, 0/2), only GT = ./., 0/0, 1/1, 2/2
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHetCount() == 0'

- select vcf sites where no sample has hetVar (no Gt = 0/1, 0/2), and at least one sample is HomVar (GT = 1/1)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHetCount() == 0 && vc.getHomVarCount() >= 1'

- site where at least 1 sample is homVar (GT = 1/1 and/or 2/2 etc.)
java -jar GenomeAnalysisTK.jar -T SelectVariants -R lyrata_genome.fa -V MY.phased_variants.Final.vcf -o My.filtered_01.vcf -select 'vc.getHomVarCount() >= 1'

- where a site is only homRef for a particular sample and is only SNP
java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getGenotype("ms01e").isHomRef()' -O test_variants.ms01e.Chr_2_3.vcf

- where a site is only homVar for a particular sample and is only SNP
java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getGenotype("ms01e").isHomVar()' -O test_variants.ms01e.Chr_2_3.vcf

- where a site is not homRef for a particular sample and is only SNP. This will also select empty sites
java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select '! vc.getGenotype("ms01e").isHomRef()' -O test_variants.ms01e.Chr_2_3.vcf

NOTE: following genotype call methods won't work:
# to get only heterozygous sites
-select 'vc.getGenotype("ms01e").isHetVar() 

# to get any variant sites
-select 'vc.getGenotype("ms01e").isVar() 

# won't work if getGenotype method is combined
-select 'vc.getGenotype("ms01e").isHomVar() || vc.getGenotype("ms01e".isHomRef()'  

So, I figured out the following methods where count helps with selecting particular site
- where a site is a heterozygous variant site for a particular sample
java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select 'vc.getHetCount() == 1'  -O test_variants.ms01e.vcf

NOTE: 'vc.getGenotype("ms01e").isHetVar()' wont work so generated the above expression.

- where a site is a variant site (either homVar or hetVar) for a particular sample and is only SNP
java -jar GenomeAnalysisTK.jar SelectVariants -R lyrata_genome.fa -V phased.MySpF1.vcf -sn ms01e -select-type SNP -select 'vc.getHetCount() == 1 || vc.getHomVarCount() == 1'  -O test_variants.ms01e.vcf

NOTE: 'vc.getGenotype("NA12878").isHomRef()' || vc.getGenotype("NA12878").isHomVar()' wont work together, so generated this above expression.

Other updates to follow.
Thanks,

> **Comments:**
> **username_1**: Restored because this is an answer to the question [i.e. a user has answered their own question].  

---

