Answers.Body: <p>First of all, bcftools might be a bit of overkill if you want one specific variant. All you need is GNU <code>grep</code> (the default on Linux):</p>
  <pre><code>zgrep -P '^(#|1\t630128\t)' file.vcf.gz 
  </code></pre>
  <p>Or, with any <code>grep</code>:</p>
  <pre><code>zgrep -E '^(#|1'$'\t''630128'$'\t'')' file.vcf.gz 
  </code></pre>
  <p>That will print all lines that start with a <code>#</code> (the headers) or that start with <code>1</code>, followed by a tab, then <code>630128</code> and another tab.</p>
  <p>Alternatively, you can use <code>awk</code>:</p>
  <pre><code>zcat file.vcf.gz |
      awk -F'\t' '/^#/ || ( $1==1 &amp;&amp; $2==630128)'
  </code></pre>
  <hr />
  <p>All that said, yes, you can also do it with <code>bcftools</code> and in most cases that is the right tool for the job so you might as well use it here too to get used to it. However, you need the right command and that's <code>view</code>, not <code>filter</code>:</p>
  <pre><code>bcftools view -r 1:630128 /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  </code></pre>
  <p>Note that the file needs to be indexed, so if you have not done so already, run</p>
  <pre><code>tabix file.vcf.gz
  </code></pre>
  <hr />
  <p>Also, the reason your unquoted command failed as it did is because the <code>&amp;</code> has special meaning to the shell, it means &quot;run this command in the background&quot;. Therefore, when you ran <code>bcftools filter -e CHROM=1&amp;POS=63018 file.vcf.gz</code> that actually ran two separate commands:</p>
  <pre><code>bcftools filter -e CHROM=1&amp;
  </code></pre>
  <p>And then:</p>
  <pre><code>POS=63018 file.vcf.gz
  </code></pre>
  <p>Which means &quot;set the variable <code>POS</code> to the value <code>63018</code> and execute the command <code>file.vcf.gz</code>. Since the file isn't executable, you got:</p>
  <pre><code>sh: 1: file.vcf.gz: Permission denied
  </code></pre>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-11-07T15:50:59.613Z
Answers.Text: First of all thanks for the quite detailed answer! This gives me a good introduction into bcftools environment. So when I use view to filter SNP's the reason for  `-filter` is then to apply some standard filter to the file? 
  Second question: If I want to remove this single variant from the file is the correct way to do this `bcftools view -r 1:630128 -e file.vcf.gz` ?
Answers.UserDisplayName: mugdi
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-11-07T16:08:35.750Z
Answers.Text: I think `filter` is more to get things like "show me variants where the AB is X and the GT is Y" and things like this, while in this case you are just focusing on a specific region, albeit a region of length 1, so `view` is the command for this. Reversing the match is a completely different question. I don't know if bcftools can do that, but `zcat file.vcf.gz | grep -vP '^1\t630128\t'` will get rid of all variants not matching that specific position.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-11-07T14:56:26.803Z
Answers.LastActivityDate: 2022-11-07T14:56:26.803Z
Answers.UserDisplayName: terdon
Body: <p>I want to filter a SNP, specifically <code>CHROM1:POS:630128</code>. I tried to use <code>bcftools</code> for this.</p>
  <p>I came up with the following variant of a bcftools call:</p>
  <pre><code>bcftools filter -e &quot;CHROM=1&amp;POS=63018&quot; /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  </code></pre>
  <p>which leads to:</p>
  <pre><code>Error: cannot use arithmetic operators to compare strings and numbers
  </code></pre>
  <p>and the call with unquoted parameters leads to</p>
  <pre><code>$ bcftools filter -e CHROM=1&amp;POS=63018 /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  
  sh: 1: /home/oem/Documents/GitHub/Genomics_prac_guide/trio/file/genome_file.vcf.gz: Permission denied
  Failed to read from standard input: unknown file type
  </code></pre>
  <ul>
  <li>Q1: Is bfctools filter the correct tool here to use?</li>
  <li>Q2: Regarding <a href="https://github.com/samtools/bcftools/issues/1028#issuecomment-497651048" rel="nofollow noreferrer">this comment</a> the first error might indicate that I have to annotate the CHROM &amp; POS Tags in the header file. Is this correct?</li>
  </ul>
  <p>An example of the file:</p>
  <pre><code>##fileformat=VCFv4.2
  ##fileDate=20221104
  ##source=PLINKv1.90
  ##contig=&lt;ID=1,length=248928329&gt;
  ##contig=&lt;ID=2,length=242157857&gt;
  ##contig=&lt;ID=3,length=198157342&gt;
  ##contig=&lt;ID=4,length=190106194&gt;
  ##contig=&lt;ID=5,length=181288140&gt;
  ##contig=&lt;ID=6,length=170595721&gt;
  ##contig=&lt;ID=7,length=159331484&gt;
  ##contig=&lt;ID=8,length=145075470&gt;
  ##contig=&lt;ID=9,length=138211490&gt;
  ##contig=&lt;ID=10,length=133659511&gt;
  ##contig=&lt;ID=11,length=135075227&gt;
  ##contig=&lt;ID=12,length=133261768&gt;
  ##contig=&lt;ID=13,length=114341522&gt;
  ##contig=&lt;ID=14,length=106873816&gt;
  ##contig=&lt;ID=15,length=101922277&gt;
  ##contig=&lt;ID=16,length=90172890&gt;
  ##contig=&lt;ID=17,length=83203771&gt;
  ##contig=&lt;ID=18,length=80252738&gt;
  ##contig=&lt;ID=19,length=58586567&gt;
  ##contig=&lt;ID=20,length=64323573&gt;
  ##contig=&lt;ID=21,length=46679699&gt;
  ##contig=&lt;ID=22,length=50776369&gt;
  ##INFO=&lt;ID=PR,Number=0,Type=Flag,Description=&quot;Provisional reference allele, may not be based on real reference genome&quot;&gt;
  ##FORMAT=&lt;ID=GT,Number=1,Type=String,Description=&quot;Genotype&quot;&gt;
  #CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  FAM001_ID001
  1   69869   rs548049170 T   .   .   .   PR  GT  0/0
  1   74792   rs13328684  N   .   .   .   PR  GT  ./.
  1   630128  rs9283150   A   .   .   .   PR  GT  0/0
  1   791532  i713426 N   .   .   .   PR  GT  ./.
  1   792461  rs116587930 G   .   .   .   PR  GT  0/0
  1   817341  rs3131972   G   .   .   .   PR  GT  0/0
  1   818725  rs12184325  C   .   .   .   PR  GT  0/0
  1   820888  rs12567639  A   .   .   .   PR  GT  0/0
  1   823656  rs114525117 G   .   .   .   PR  GT  0/0
  1   841166  rs12124819  G   A   .   .   PR  GT  0/1
  1   858952  rs12127425  G   .   .   .   PR  GT  0/0
  1   866156  rs79373928  T   .   .   .   PR  GT  0/0
  1   880041  rs72888853  N   .   .   .   PR  GT  ./.
  </code></pre>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-11-07T14:39:17.783Z
Comments.Text: Can you give us an example of your file?
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-11-07T14:47:19.367Z
Comments.Text: Done! See edited post.
Comments.UserDisplayName: mugdi
Tags: <vcf><filtering><bcftools>
Title: Filter with bcftools
ThreadText: # Filter with bcftools
  
  (Asked by: username_0 on 2022-11-07)
  
  I want to filter a SNP, specifically `CHROM1:POS:630128`. I tried to use `bcftools` for this.
  
  I came up with the following variant of a bcftools call:
  
  ```
  bcftools filter -e "CHROM=1&POS=63018" /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  
  ```
  
  which leads to:
  
  ```
  Error: cannot use arithmetic operators to compare strings and numbers
  
  ```
  
  and the call with unquoted parameters leads to
  
  ```
  $ bcftools filter -e CHROM=1&POS=63018 /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  
  sh: 1: /home/oem/Documents/GitHub/Genomics_prac_guide/trio/file/genome_file.vcf.gz: Permission denied
  Failed to read from standard input: unknown file type
  
  ```
  
  * Q1: Is bfctools filter the correct tool here to use?
  * Q2: Regarding [this comment](https://github.com/samtools/bcftools/issues/1028#issuecomment-497651048) the first error might indicate that I have to annotate the CHROM & POS Tags in the header file. Is this correct?
  
  An example of the file:
  
  ```
  ##fileformat=VCFv4.2
  ##fileDate=20221104
  ##source=PLINKv1.90
  ##contig=<ID=1,length=248928329>
  ##contig=<ID=2,length=242157857>
  ##contig=<ID=3,length=198157342>
  ##contig=<ID=4,length=190106194>
  ##contig=<ID=5,length=181288140>
  ##contig=<ID=6,length=170595721>
  ##contig=<ID=7,length=159331484>
  ##contig=<ID=8,length=145075470>
  ##contig=<ID=9,length=138211490>
  ##contig=<ID=10,length=133659511>
  ##contig=<ID=11,length=135075227>
  ##contig=<ID=12,length=133261768>
  ##contig=<ID=13,length=114341522>
  ##contig=<ID=14,length=106873816>
  ##contig=<ID=15,length=101922277>
  ##contig=<ID=16,length=90172890>
  ##contig=<ID=17,length=83203771>
  ##contig=<ID=18,length=80252738>
  ##contig=<ID=19,length=58586567>
  ##contig=<ID=20,length=64323573>
  ##contig=<ID=21,length=46679699>
  ##contig=<ID=22,length=50776369>
  ##INFO=<ID=PR,Number=0,Type=Flag,Description="Provisional reference allele, may not be based on real reference genome">
  ##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
  #CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  FAM001_ID001
  1   69869   rs548049170 T   .   .   .   PR  GT  0/0
  1   74792   rs13328684  N   .   .   .   PR  GT  ./.
  1   630128  rs9283150   A   .   .   .   PR  GT  0/0
  1   791532  i713426 N   .   .   .   PR  GT  ./.
  1   792461  rs116587930 G   .   .   .   PR  GT  0/0
  1   817341  rs3131972   G   .   .   .   PR  GT  0/0
  1   818725  rs12184325  C   .   .   .   PR  GT  0/0
  1   820888  rs12567639  A   .   .   .   PR  GT  0/0
  1   823656  rs114525117 G   .   .   .   PR  GT  0/0
  1   841166  rs12124819  G   A   .   .   PR  GT  0/1
  1   858952  rs12127425  G   .   .   .   PR  GT  0/0
  1   866156  rs79373928  T   .   .   .   PR  GT  0/0
  1   880041  rs72888853  N   .   .   .   PR  GT  ./.
  
  ```
  
  **Tags:** vcf, filtering, bcftools
  
  > **Comments:**
  > **username_1**: Can you give us an example of your file?  
  > **username_0**: Done! See edited post.  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_1 on 2022-11-07. Score: 3 | Accepted answer
  
  First of all, bcftools might be a bit of overkill if you want one specific variant. All you need is GNU `grep` (the default on Linux):
  
  ```
  zgrep -P '^(#|1\t630128\t)' file.vcf.gz 
  
  ```
  
  Or, with any `grep`:
  
  ```
  zgrep -E '^(#|1'$'\t''630128'$'\t'')' file.vcf.gz 
  
  ```
  
  That will print all lines that start with a `#` (the headers) or that start with `1`, followed by a tab, then `630128` and another tab.
  
  Alternatively, you can use `awk`:
  
  ```
  zcat file.vcf.gz |
      awk -F'\t' '/^#/ || ( $1==1 && $2==630128)'
  
  ```
  
  ---
  
  All that said, yes, you can also do it with `bcftools` and in most cases that is the right tool for the job so you might as well use it here too to get used to it. However, you need the right command and that's `view`, not `filter`:
  
  ```
  bcftools view -r 1:630128 /home/oem/Documents/GitHub/Genomics_prac_guide/trio/files/genome_file.vcf.gz 
  
  ```
  
  Note that the file needs to be indexed, so if you have not done so already, run
  
  ```
  tabix file.vcf.gz
  
  ```
  
  ---
  
  Also, the reason your unquoted command failed as it did is because the `&` has special meaning to the shell, it means "run this command in the background". Therefore, when you ran `bcftools filter -e CHROM=1&POS=63018 file.vcf.gz` that actually ran two separate commands:
  
  ```
  bcftools filter -e CHROM=1&
  
  ```
  
  And then:
  
  ```
  POS=63018 file.vcf.gz
  
  ```
  
  Which means "set the variable `POS` to the value `63018` and execute the command `file.vcf.gz`. Since the file isn't executable, you got:
  
  ```
  sh: 1: file.vcf.gz: Permission denied
  
  ```
  
  > **Comments:**
  > **username_0**: First of all thanks for the quite detailed answer! This gives me a good introduction into bcftools environment. So when I use view to filter SNP's the reason for  `-filter` is then to apply some standard filter to the file?  Second question: If I want to remove this single variant from the file is the correct way to do this `bcftools view -r 1:630128 -e file.vcf.gz` ?  
  > **username_1**: I think `filter` is more to get things like "show me variants where the AB is X and the GT is Y" and things like this, while in this case you are just focusing on a specific region, albeit a region of length 1, so `view` is the command for this. Reversing the match is a completely different question. I don't know if bcftools can do that, but `zcat file.vcf.gz | grep -vP '^1\t630128\t'` will get rid of all variants not matching that specific position.  
  
  ---
  
  