Answers.Body: <p>You can use <code>grep</code> for this. The problem is the <code>-w</code> flag which means "match entire words", but VCF header lines start with <code>##</code>, not <code>#</code>. Here are some alternative <code>grep</code> commands:</p>
  
  <pre><code>grep -E '^(#|2[[:space:]])' my.vcf 
  </code></pre>
  
  <p>Or, if you have GNU <code>grep</code> (default on Linux):</p>
  
  <pre><code>grep -P '^(#|2\t)' my.vcf 
  </code></pre>
  
  <p>To match either <code>2</code> or <code>chr2</code>:</p>
  
  <pre><code>grep -E '^(#|(chr2|2)[[:space:]])' my.vcf 
  </code></pre>
  
  <p>or</p>
  
  <pre><code>grep -P '^(#|(chr){0,1}2\t)' my.vcf 
  </code></pre>
  
  <p>Alternatively, you can use <code>awk</code>:</p>
  
  <pre><code>awk '/^#/ || $1=="2" || $2=="chr2"' my.vcf 
  </code></pre>
  
  <p>Finally, if you want to use <code>vcftools</code>, you need to also use the <code>--recode</code> option:</p>
  
  <blockquote>
    <p>--recode</p>
    
    <p>--recode-bcf</p>
    
    <p>These  options  are used to generate a new file in either VCF or BCF from the input VCF or BCF file after applying   the
    filtering options specified by the user. The output file has the
    suffix  ".recode.vcf"  or  ".recode.bcf".  By   default,  the  INFO
    fields are removed from the output file, as the INFO values may be
    invalidated by the recoding</p>
  </blockquote>
  
  <p>(e.g. the total depth may need to be recalculated if individuals are removed). This behavior may be  overriden  by
    the following options. By default, BCF files are written out as BGZF compressed files.</p>
  
  <p>But that will remove the INFO fields as explained above, so just use normal text-parsing tools like the ones I suggested.</p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-25T20:25:01.847Z
Answers.Text: grep (at least the version on Linux) can use multiple `-e` commands as expressions to search for; lines matching any expression will be displayed.
Answers.UserDisplayName: gringer
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-25T20:45:26.240Z
Answers.Text: @gringer yes, that's [POSIX](http://pubs.opengroup.org/onlinepubs/9699919799/utilities/grep.html). Personally, I prefer ERE or PCRE syntax over multiple `-e`, but both work.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-26T17:12:18.180Z
Answers.Text: These are both great answers. I will probably mark @gringer as the accepted as there are less opportunities for readers to make mistakes, but personally, I'm using these grep commands. Thank you!
Answers.UserDisplayName: EB2127
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-26T17:33:46.797Z
Answers.Text: @EB2127 you're very welcome, and please feel very free to mark whichever you prefer. I personally prefer rolling my own parser scripts, but with formats as complex as vcf, it is indeed usually better to use dedicated tools. I would go for gringer's answer as well in your place :)
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-25T19:13:22.367Z
Answers.LastActivityDate: 2018-01-25T19:13:22.367Z
Answers.UserDisplayName: terdon
Answers.Body: <p>You can do this more easily with the GATK tool SelectVariants, which allows you to subset variants based on various criteria. Here's the doc page:</p>
  
  <p><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360041849831-SelectVariants" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360041849831-SelectVariants</a></p>
  
  <p>For this case (subsetting by chromosome) you would simply do:</p>
  
  <pre><code>gatk SelectVariants \
       -R reference.fasta \
       -V input.vcf \
       -L chr2 \
       -O output.chr2.vcf
  </code></pre>
  
  <p>If you wanted to do this for several chromosomes, or regions within chromosomes, you would give the tool a file of intervals with <code>-L</code> instead of the chromosome name. Suitable interval file formats are described here: </p>
  
  <p><a href="https://gatk.broadinstitute.org/hc/en-us/articles/360035531852-Intervals-and-interval-lists" rel="nofollow noreferrer">https://gatk.broadinstitute.org/hc/en-us/articles/360035531852-Intervals-and-interval-lists</a></p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-09-03T06:02:23.973Z
Answers.Text: In the folder of the command ( ./gatk ), the Fasta index file and the Fasta dict file of the reference.fasta should also be included.
Answers.UserDisplayName: envs_h_gang_5
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-04-23T15:56:08.440Z
Answers.LastActivityDate: 2020-04-23T15:56:08.440Z
Answers.UserDisplayName: Geraldine_VdAuwera
Answers.Body: <p>This seems to work for me with <code>bcftools filter</code> and the <code>-r</code> or <code>-R</code> argument.</p>
  
  <blockquote>
  <pre><code>   -r, --regions chr|chr:pos|chr:from-to|chr:from-[,...]
  </code></pre>
    
    <p>Comma-separated list of regions, see also <strong>-R</strong>, <strong>--regions-file</strong>. Note that <strong>-r</strong> cannot be used in combination with <strong>-R</strong>.</p>
  
  <pre><code>   -R, --regions-file FILE
  </code></pre>
    
    <p>Regions can be specified either on command line or in a VCF, BED, or tab-delimited file (the default). The columns of the
               tab-delimited file are:
               CHROM, POS, and, optionally, POS_TO, where positions are 1-based and inclusive. Uncompressed files are stored in memory, while
    bgzip-compressed and
               tabix-indexed region files are streamed. Note that sequence names must match exactly, "chr20" is not the same as "20". Also note
    that chromosome
               ordering in FILE will be respected, the VCF will be processed in the order in which chromosomes first appear in FILE.
    However, within chromosomes, the
               VCF will always be processed in ascending genomic coordinate order no matter what order they appear in FILE. Note that
    overlapping regions in FILE can
               result in duplicated out of order positions in the output. This option requires indexed VCF/BCF files. Note that <strong>-R</strong> cannot be
    used in combination with <strong>-r</strong>.</p>
  </blockquote>
  
  <p>Here's an example:</p>
  
  <pre><code>$ bcftools filter vcf_nocomp_merge_geno98.vcf.gz -r 4 | head -n 38 | colrm 100 1000000
  ##fileformat=VCFv4.2
  ##FILTER=&lt;ID=PASS,Description="All filters passed"&gt;
  ##fileDate=20171222
  ##source=PLINKv1.90
  ##contig=&lt;ID=1,length=249167691&gt;
  ##contig=&lt;ID=2,length=242695901&gt;
  ##contig=&lt;ID=3,length=197800245&gt;
  ##contig=&lt;ID=4,length=190915651&gt;
  ##contig=&lt;ID=5,length=180666277&gt;
  ##contig=&lt;ID=6,length=170877445&gt;
  ##contig=&lt;ID=7,length=159086805&gt;
  ##contig=&lt;ID=8,length=146293415&gt;
  ##contig=&lt;ID=9,length=141018424&gt;
  ##contig=&lt;ID=10,length=135434552&gt;
  ##contig=&lt;ID=11,length=134938471&gt;
  ##contig=&lt;ID=12,length=133763353&gt;
  ##contig=&lt;ID=13,length=115045730&gt;
  ##contig=&lt;ID=14,length=107285438&gt;
  ##contig=&lt;ID=15,length=102369712&gt;
  ##contig=&lt;ID=16,length=90141356&gt;
  ##contig=&lt;ID=17,length=81006630&gt;
  ##contig=&lt;ID=18,length=78014583&gt;
  ##contig=&lt;ID=19,length=59071322&gt;
  ##contig=&lt;ID=20,length=62906515&gt;
  ##contig=&lt;ID=21,length=48077813&gt;
  ##contig=&lt;ID=22,length=51156934&gt;
  ##contig=&lt;ID=23,length=154847490&gt;
  ##INFO=&lt;ID=PR,Number=0,Type=Flag,Description="Provisional reference allele, may not be based on rea
  ##FORMAT=&lt;ID=GT,Number=1,Type=String,Description="Genotype"&gt;
  ##bcftools_filterVersion=1.2-187-g1a55e45+htslib-1.2.1-256-ga356746
  ##bcftools_filterCommand=filter -r 4 vcf_nocomp_merge_geno98.vcf.gz
  #CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  1_125   3_5555  4_7565  5_9
  4   71566   rs13125929  T   C   .   PASS    PR  GT  0/1 0/1 0/0
  4   87113   rs2006748   T   C   .   PASS    PR  GT  0/1 0/0 0/0
  4   110646  rs11727494  C   T   .   PASS    PR  GT  0/1 0/0 1/1
  4   142550  rs11735742  T   C   .   PASS    PR  GT  0/1 0/0 0/0
  4   200631  rs6826124   A   G   .   PASS    PR  GT  0/1 1/1 0/1
  4   221623  rs7695945   G   T   .   PASS    PR  GT  0/0 0/0 0/0
  </code></pre>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2018-01-25T20:23:28.003Z
Answers.LastActivityDate: 2018-01-25T20:23:28.003Z
Answers.UserDisplayName: gringer
Body: <p>I would like to subset a VCF which only has chromosome 2. </p>
  
  <p>The problem with using various <code>grep</code> commands, e.g. </p>
  
  <pre><code>grep -w '^#\|^2' my.vcf &gt; my_new.vcf
  </code></pre>
  
  <p>or if there's a 'chr' prefix</p>
  
  <pre><code>grep -w '^#\|^chr2' my.vcf &gt; my_new.vcf
  </code></pre>
  
  <p>is that this will remove the header. </p>
  
  <p>I've been trying to use <code>vcftools</code> with the following command:</p>
  
  <pre><code>vcftools --vcf input.vcf  --chr 2 --out subset
  </code></pre>
  
  <p>which does not output a VCF <code>subset.vcf</code> as expected, but rather <code>subset.log</code>:</p>
  
  <pre><code>VCFtools - 0.1.15
  (C) Adam Auton and Anthony Marcketta 2009
  
  Parameters as interpreted:
      --vcf input.vcf
      --chr 2
      --out subset
  
  After filtering, kept 2 out of 2 Individuals
  After filtering, kept 80 out of a possible 720 Sites
  Run Time = 0.00 seconds
  </code></pre>
  
  <p>Given that the run time is <code>Run Time = 0.00 seconds</code>, I guess there's an error here. </p>
  
  <p>How can I subset a VCF and keep the header? Is there an option with bcftools instead maybe?</p>
  
Tags: <vcf><vcftools>
Title: How to subset a VCF by chromosome and keep the header?
ThreadText: # How to subset a VCF by chromosome and keep the header?
  
  (Asked by: username_0 on 2018-01-25)
  
  I would like to subset a VCF which only has chromosome 2.
  
  The problem with using various `grep` commands, e.g.
  
  ```
  grep -w '^#\|^2' my.vcf > my_new.vcf
  
  ```
  
  or if there's a 'chr' prefix
  
  ```
  grep -w '^#\|^chr2' my.vcf > my_new.vcf
  
  ```
  
  is that this will remove the header.
  
  I've been trying to use `vcftools` with the following command:
  
  ```
  vcftools --vcf input.vcf  --chr 2 --out subset
  
  ```
  
  which does not output a VCF `subset.vcf` as expected, but rather `subset.log`:
  
  ```
  VCFtools - 0.1.15
  (C) Adam Auton and Anthony Marcketta 2009
  
  Parameters as interpreted:
      --vcf input.vcf
      --chr 2
      --out subset
  
  After filtering, kept 2 out of 2 Individuals
  After filtering, kept 80 out of a possible 720 Sites
  Run Time = 0.00 seconds
  
  ```
  
  Given that the run time is `Run Time = 0.00 seconds`, I guess there's an error here.
  
  How can I subset a VCF and keep the header? Is there an option with bcftools instead maybe?
  
  **Tags:** vcf, vcftools
  
  
  ---
  
  ## 3 answers
  
  ### Answer by username_2 on 2018-01-25. Score: 11 | Accepted answer
  
  This seems to work for me with `bcftools filter` and the `-r` or `-R` argument.
  
  > ```
  >    -r, --regions chr|chr:pos|chr:from-to|chr:from-[,...]
  >
  > ```
  >
  > Comma-separated list of regions, see also **-R**, **--regions-file**. Note that **-r** cannot be used in combination with **-R**.
  >
  > ```
  >    -R, --regions-file FILE
  >
  > ```
  >
  > Regions can be specified either on command line or in a VCF, BED, or tab-delimited file (the default). The columns of the
  > tab-delimited file are:
  > CHROM, POS, and, optionally, POS\_TO, where positions are 1-based and inclusive. Uncompressed files are stored in memory, while
  > bgzip-compressed and
  > tabix-indexed region files are streamed. Note that sequence names must match exactly, "chr20" is not the same as "20". Also note
  > that chromosome
  > ordering in FILE will be respected, the VCF will be processed in the order in which chromosomes first appear in FILE.
  > However, within chromosomes, the
  > VCF will always be processed in ascending genomic coordinate order no matter what order they appear in FILE. Note that
  > overlapping regions in FILE can
  > result in duplicated out of order positions in the output. This option requires indexed VCF/BCF files. Note that **-R** cannot be
  > used in combination with **-r**.
  
  Here's an example:
  
  ```
  $ bcftools filter vcf_nocomp_merge_geno98.vcf.gz -r 4 | head -n 38 | colrm 100 1000000
  ##fileformat=VCFv4.2
  ##FILTER=<ID=PASS,Description="All filters passed">
  ##fileDate=20171222
  ##source=PLINKv1.90
  ##contig=<ID=1,length=249167691>
  ##contig=<ID=2,length=242695901>
  ##contig=<ID=3,length=197800245>
  ##contig=<ID=4,length=190915651>
  ##contig=<ID=5,length=180666277>
  ##contig=<ID=6,length=170877445>
  ##contig=<ID=7,length=159086805>
  ##contig=<ID=8,length=146293415>
  ##contig=<ID=9,length=141018424>
  ##contig=<ID=10,length=135434552>
  ##contig=<ID=11,length=134938471>
  ##contig=<ID=12,length=133763353>
  ##contig=<ID=13,length=115045730>
  ##contig=<ID=14,length=107285438>
  ##contig=<ID=15,length=102369712>
  ##contig=<ID=16,length=90141356>
  ##contig=<ID=17,length=81006630>
  ##contig=<ID=18,length=78014583>
  ##contig=<ID=19,length=59071322>
  ##contig=<ID=20,length=62906515>
  ##contig=<ID=21,length=48077813>
  ##contig=<ID=22,length=51156934>
  ##contig=<ID=23,length=154847490>
  ##INFO=<ID=PR,Number=0,Type=Flag,Description="Provisional reference allele, may not be based on rea
  ##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
  ##bcftools_filterVersion=1.2-187-g1a55e45+htslib-1.2.1-256-ga356746
  ##bcftools_filterCommand=filter -r 4 vcf_nocomp_merge_geno98.vcf.gz
  #CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  1_125   3_5555  4_7565  5_9
  4   71566   rs13125929  T   C   .   PASS    PR  GT  0/1 0/1 0/0
  4   87113   rs2006748   T   C   .   PASS    PR  GT  0/1 0/0 0/0
  4   110646  rs11727494  C   T   .   PASS    PR  GT  0/1 0/0 1/1
  4   142550  rs11735742  T   C   .   PASS    PR  GT  0/1 0/0 0/0
  4   200631  rs6826124   A   G   .   PASS    PR  GT  0/1 1/1 0/1
  4   221623  rs7695945   G   T   .   PASS    PR  GT  0/0 0/0 0/0
  
  ```
  
  
  ---
  
  ### Answer by username_1 on 2018-01-25. Score: 4
  
  You can use `grep` for this. The problem is the `-w` flag which means "match entire words", but VCF header lines start with `##`, not `#`. Here are some alternative `grep` commands:
  
  ```
  grep -E '^(#|2[[:space:]])' my.vcf 
  
  ```
  
  Or, if you have GNU `grep` (default on Linux):
  
  ```
  grep -P '^(#|2\t)' my.vcf 
  
  ```
  
  To match either `2` or `chr2`:
  
  ```
  grep -E '^(#|(chr2|2)[[:space:]])' my.vcf 
  
  ```
  
  or
  
  ```
  grep -P '^(#|(chr){0,1}2\t)' my.vcf 
  
  ```
  
  Alternatively, you can use `awk`:
  
  ```
  awk '/^#/ || $1=="2" || $2=="chr2"' my.vcf 
  
  ```
  
  Finally, if you want to use `vcftools`, you need to also use the `--recode` option:
  
  > --recode
  >
  > --recode-bcf
  >
  > These options are used to generate a new file in either VCF or BCF from the input VCF or BCF file after applying the
  > filtering options specified by the user. The output file has the
  > suffix ".recode.vcf" or ".recode.bcf". By default, the INFO
  > fields are removed from the output file, as the INFO values may be
  > invalidated by the recoding
  
  (e.g. the total depth may need to be recalculated if individuals are removed). This behavior may be overriden by
  the following options. By default, BCF files are written out as BGZF compressed files.
  
  But that will remove the INFO fields as explained above, so just use normal text-parsing tools like the ones I suggested.
  
  > **Comments:**
  > **username_2**: grep (at least the version on Linux) can use multiple `-e` commands as expressions to search for; lines matching any expression will be displayed.  
  > **username_1**: @username_2 yes, that's [POSIX](http://pubs.opengroup.org/onlinepubs/9699919799/utilities/grep.html). Personally, I prefer ERE or PCRE syntax over multiple `-e`, but both work.  
  > **username_0**: These are both great answers. I will probably mark @username_2 as the accepted as there are less opportunities for readers to make mistakes, but personally, I'm using these grep commands. Thank you!  
  > **username_1**: @username_0 you're very welcome, and please feel very free to mark whichever you prefer. I personally prefer rolling my own parser scripts, but with formats as complex as vcf, it is indeed usually better to use dedicated tools. I would go for username_2's answer as well in your place :)  
  
  ---
  
  ### Answer by username_3 on 2020-04-23. Score: 3
  
  You can do this more easily with the GATK tool SelectVariants, which allows you to subset variants based on various criteria. Here's the doc page:
  
  <https://gatk.broadinstitute.org/hc/en-us/articles/360041849831-SelectVariants>
  
  For this case (subsetting by chromosome) you would simply do:
  
  ```
  gatk SelectVariants \
       -R reference.fasta \
       -V input.vcf \
       -L chr2 \
       -O output.chr2.vcf
  
  ```
  
  If you wanted to do this for several chromosomes, or regions within chromosomes, you would give the tool a file of intervals with `-L` instead of the chromosome name. Suitable interval file formats are described here:
  
  <https://gatk.broadinstitute.org/hc/en-us/articles/360035531852-Intervals-and-interval-lists>
  
  > **Comments:**
  > **username_4**: In the folder of the command ( ./gatk ), the Fasta index file and the Fasta dict file of the reference.fasta should also be included.  
  
  ---
  
  