text: user1515: I'm doing some analysis and I need to subset a large VCF file (~8GB gziped) given a bed interval and identify within a subset of rsid. 
  Unfortunately, both my normal choices to do this analysis (snpSift and bedtools) are taking way to long or failing due to memory issues in my local computer and a remote server.
  Do you guys know any other options or suggestions to speed up this process?
  Follow the commands I use:
  bedtools intersect -a <myvcf>.vcf.gz -b <myinterval>.bed -wa | \
      java -Xmx10g -jar snpSift.jar filter --set <myrsid>.txt "ID in SET[0]"
  
  or
  gzcat <myvcf>.vcf.gz | \
      java -Xmx10g -jar snpSift.jar intervals <mybed>.bed | \
      java -Xmx10g -jar snpSift.jar filter --set <myrsid>.txt "ID in SET[0]"
  
  The bedtools command usually fail due to unknown reason and SnpSift runs over 6 hours even given 10GB of ram. My local machine have 8GB of RAM, but the server have 32GB.
  
  user1310: Here is another way to do it in 2019:
  $ bcftools view -T regions.bed -i 'ID=@<myrsid>.txt' input.vcf > output.vcf
  
  Things can be sped up if we compress and index the input.vcf. Then bcftools will have near-immediate random access to the position.
  $ bgzip -c input.vcf > input.vcf.gz
  $ tabix input.vcf.gz
  $ bcftools view -T regions.bed -i 'ID=@<myrsid>.txt' input.vcf.gz > output.vcf
  
  Internally bcftools works with bcf. Whenever the input and output is vcf it converts vcf>bcf>vcf. When dealing with large files, this produces a large overhead. Converting vcf to bcf first can reduce runtime for further downstream analyses dramatically. In my experience 2-3x faster.
  $ bcftools view -Ob input.vcf > input.bcf
  $ bcftools index input.bcf
  $ bcftools view -Ob -T regions.bed -i 'ID=@<myrsid>.txt' input.bcf > output.bcf
  
  
  user1515: I was tinkering with the command and was able to complete the execution inverting the order of rsid and bed intervals filtering. The command is as follows:
  gzcat <myvcf>.vcf.gz | \
      java -Xmx10g -jar snpSift.jar filter --set <myrsid>.txt "ID in SET[0]" | \
      java -Xmx10g -jar snpSift.jar intervals <mybed>.bed
  
  Maybe my rsid subset is smaller than my BED intervals and thus the process is quicker and more efficient.
  Anyway thanks everyone
  
  user71: this:
  bedtools intersect -a <myvcf>.vcf.gz -b <myinterval>.bed -wa | \
      java -Xmx10g -jar snpSift.jar filter --set <myrsid>.txt "ID in SET[0]"
  
  can be replaced with one GATK SelectVariants https://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_gatk_tools_walkers_variantutils_SelectVariants.php
  java -jar gatk.jar -T SelectVariants -R ref.fa \
    -L myinterval.bed --keepIDs myrsid.txt -V myvcf.vcf.gz -o out.vcf
  
  and furthermore, you can also work in parallel by splitting your bed file (e.g. by chromosomes), and merge the results.
  
  user776: Generically, with BEDOPS:
  $ vcf2bed < <(gunzip -c snps.vcf) | bedops -e 1 - myRegions.bed > answer.bed
  
  Or:
  $ vcf2bed < <(gunzip -c snps.vcf) | bedmap --echo --echo-map-id --delim '\t' myRegions.bed - > answer.bed
  
  Etc.
  BEDOPS supports input streams; use them where you can.
  That said, if you're going to query the VCF data frequently, convert once and use it in downstream operations.
  $ vcf2bed < <(gunzip -c snps.vcf) > snps.bed
  
  Then:
  $ bedops -e 1 snps.bed myRegions.bed > answer.bed
  
  Etc.
  If you have a computational cluster, you can further parallelize work by splitting SNPs on chromosomes with bedextract:
  $ for chr in `bedextract --list-chr snps.bed`; do echo ${chr}; bedextract ${chr} snps.bed > snps.${chr}.bed; done
  
  Then submit a batch job for each chromosome:
  $ bedmap --chrom ${chr} --echo --echo-map-id --delim '\t' myRegions.bed snps.${chr}.bed > answer.${chr}.bed
  
  Using bedextract and adding --chrom ${chr} to your bedmap call will yield large speed improvements over other toolkits, as these options focus work on the chromosome of interest and nothing else.
  
original_text: 1515: <p>I'm doing some analysis and I need to subset a large VCF file (~8GB gziped) given a bed interval and identify within a subset of rsid. </p>
  
  <p>Unfortunately, both my normal choices to do this analysis (<code>snpSift</code> and <code>bedtools</code>) are taking way to long or failing due to memory issues in my local computer and a remote server.</p>
  
  <p>Do you guys know any other options or suggestions to speed up this process?</p>
  
  <p>Follow the commands I use:</p>
  
  <pre><code>bedtools intersect -a &lt;myvcf&gt;.vcf.gz -b &lt;myinterval&gt;.bed -wa | \
      java -Xmx10g -jar snpSift.jar filter --set &lt;myrsid&gt;.txt "ID in SET[0]"
  </code></pre>
  
  <p>or</p>
  
  <pre><code>gzcat &lt;myvcf&gt;.vcf.gz | \
      java -Xmx10g -jar snpSift.jar intervals &lt;mybed&gt;.bed | \
      java -Xmx10g -jar snpSift.jar filter --set &lt;myrsid&gt;.txt "ID in SET[0]"
  </code></pre>
  
  <p>The <code>bedtools</code> command usually fail due to unknown reason and <code>SnpSift</code> runs over 6 hours even given 10GB of ram. My local machine have 8GB of RAM, but the server have 32GB.</p>
   1310: <p>Here is another way to do it in 2019:</p>
  
  <pre><code>$ bcftools view -T regions.bed -i 'ID=@&lt;myrsid&gt;.txt' input.vcf &gt; output.vcf
  </code></pre>
  
  <p>Things can be sped up if we compress and index the <code>input.vcf</code>. Then <code>bcftools</code> will have near-immediate random access to the position.</p>
  
  <pre><code><span class="math-container">$ bgzip -c input.vcf &gt; input.vcf.gz
  $</span> tabix input.vcf.gz
  $ bcftools view -T regions.bed -i 'ID=@&lt;myrsid&gt;.txt' input.vcf.gz &gt; output.vcf
  </code></pre>
  
  <p>Internally <code>bcftools</code> works with <code>bcf</code>. Whenever the input and output is <code>vcf</code> it converts <code>vcf</code>><code>bcf</code>><code>vcf</code>. When dealing with large files, this produces a large overhead. Converting <code>vcf</code> to <code>bcf</code> first can reduce runtime for further downstream analyses dramatically. In my experience 2-3x faster.</p>
  
  <pre><code><span class="math-container">$ bcftools view -Ob input.vcf &gt; input.bcf
  $</span> bcftools index input.bcf
  $ bcftools view -Ob -T regions.bed -i 'ID=@&lt;myrsid&gt;.txt' input.bcf &gt; output.bcf
  </code></pre>
   1515: <p>I was tinkering with the command and was able to complete the execution inverting the order of rsid and bed intervals filtering. The command is as follows:</p>
  
  <pre><code>gzcat &lt;myvcf&gt;.vcf.gz | \
      java -Xmx10g -jar snpSift.jar filter --set &lt;myrsid&gt;.txt "ID in SET[0]" | \
      java -Xmx10g -jar snpSift.jar intervals &lt;mybed&gt;.bed
  </code></pre>
  
  <p>Maybe my rsid subset is smaller than my BED intervals and thus the process is quicker and more efficient.</p>
  
  <p>Anyway thanks everyone</p>
   71: <p>this:</p>
  
  <pre><code>bedtools intersect -a &lt;myvcf&gt;.vcf.gz -b &lt;myinterval&gt;.bed -wa | \
      java -Xmx10g -jar snpSift.jar filter --set &lt;myrsid&gt;.txt "ID in SET[0]"
  </code></pre>
  
  <p>can be replaced with one GATK SelectVariants <a href="https://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_gatk_tools_walkers_variantutils_SelectVariants.php" rel="nofollow noreferrer">https://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_gatk_tools_walkers_variantutils_SelectVariants.php</a></p>
  
  <pre><code>java -jar gatk.jar -T SelectVariants -R ref.fa \
    -L myinterval.bed --keepIDs myrsid.txt -V myvcf.vcf.gz -o out.vcf
  </code></pre>
  
  <p>and furthermore, you can also work in parallel by splitting your bed file (e.g. by chromosomes), and merge the results.</p>
   776: <p>Generically, with <a href="http://bedops.readthedocs.io/en/latest/" rel="nofollow noreferrer">BEDOPS</a>:</p>
  
  <pre><code>$ vcf2bed &lt; &lt;(gunzip -c snps.vcf) | bedops -e 1 - myRegions.bed &gt; answer.bed
  </code></pre>
  
  <p>Or:</p>
  
  <pre><code>$ vcf2bed &lt; &lt;(gunzip -c snps.vcf) | bedmap --echo --echo-map-id --delim '\t' myRegions.bed - &gt; answer.bed
  </code></pre>
  
  <p>Etc.</p>
  
  <p>BEDOPS supports input streams; use them where you can.</p>
  
  <p>That said, if you're going to query the VCF data frequently, convert once and use it in downstream operations.</p>
  
  <pre><code>$ vcf2bed &lt; &lt;(gunzip -c snps.vcf) &gt; snps.bed
  </code></pre>
  
  <p>Then:</p>
  
  <pre><code>$ bedops -e 1 snps.bed myRegions.bed &gt; answer.bed
  </code></pre>
  
  <p>Etc.</p>
  
  <p>If you have a computational cluster, you can further parallelize work by splitting SNPs on chromosomes with <code>bedextract</code>:</p>
  
  <pre><code>$ for chr in `bedextract --list-chr snps.bed`; do echo ${chr}; bedextract ${chr} snps.bed &gt; snps.${chr}.bed; done
  </code></pre>
  
  <p>Then submit a batch job for each chromosome:</p>
  
  <pre><code>$ bedmap --chrom ${chr} --echo --echo-map-id --delim '\t' myRegions.bed snps.${chr}.bed &gt; answer.${chr}.bed
  </code></pre>
  
  <p>Using <code>bedextract</code> and adding <code>--chrom ${chr}</code> to your <code>bedmap</code> call will yield large speed improvements over other toolkits, as these options focus work on the chromosome of interest and nothing else.</p>
  