text: user1099: 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.
  
  user1099: 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,
  
original_text: 1099: <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>
   1099: <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>
  