<p>Just to highlight that all the steps can be done within bcftools capabilities, and since I can't just comment on @blmoore 's answer: </p>

<pre><code>bcftools view --types indels &lt;vcf&gt; |
  bcftools norm -m - |
  bcftools filter --include 'strlen(REF)&lt;strlen(ALT)' |
  bcftools view -H
</code></pre>

<p>more bult-in functions for 'bcftools filter' <a href="https://samtools.github.io/bcftools/bcftools.html" rel="noreferrer">here</a></p>

CC BY-SA 4.0
2018-09-22T22:00:55.567Z
Very nice, I had only ever used bcftools view to uncompress and then parse the output into either the header, or the data records.
mdperry
CC BY-SA 4.0
2018-06-29T17:03:56.600Z
2018-07-02T13:49:52.417Z
2018-07-02T13:49:52.417Z
JRodrigoF
<p>One method is to decompose multi-allelic records so that they're represented as one-allele, one-record using <a href="https://github.com/atks/vt" rel="nofollow noreferrer">vt</a> or similar:</p>

<pre><code>bcftools view -v indels &lt;vcf&gt; |
  vt decompose - |
  bcftools view -H |
  awk '{if(length($5)&gt;length($4)) print}'
</code></pre>

<p>Some decomposed alleles will come with excess reference padding. To left-shift and trim these, add a <code>normalize</code> step (NB matching these back to your input VCF becomes non-trivial):</p>

<pre><code>bcftools view -v indels &lt;vcf&gt; |
  vt decompose - |
  vt normalize -r &lt;reference.fasta - |
  awk '{if(length($4)==1) print}'
</code></pre>

<p><strong>edit:</strong> As gringer suggests, this can also be done without vt:</p>

<pre><code>bcftools view -Ou -v indels &lt;vcf&gt; |
  bcftools norm -Ou -Nm - |
  bcftools view -H |
  awk '{if(length($5)&gt;length($4)) print}'
</code></pre>

<p>To also include complex alleles, use <code>view -V snps</code> (!snvs) instead of <code>-v indels</code></p>

CC BY-SA 3.0
2017-06-16T12:24:06.930Z
`bcftools norm` will also split multi-allelic records with the `-m -` option
gringer
CC BY-SA 3.0
2017-06-16T11:24:00.487Z
2017-06-16T14:21:01.260Z
2017-06-16T14:21:01.260Z
blmoore
<p>Using <a href="http://lindenb.github.io/jvarkit/VCFFilterJS.html" rel="nofollow noreferrer">vcffilterjs</a></p>

<ul>
<li>get the length of the REF;</li>
<li>loop over the ALT, ignore the symbolic</li>
<li>accept the variant if it's an insertion , eq:  len(ALT)>len(REF) </li>
</ul>

<p>.</p>

<pre><code>java -jar dist/vcffilterjs.jar -e 'function accept(vc){var a=vc.getAlleles();var lenRef=a.get(0).length();for(i=1;i&lt;a.size();++i) {var alt=a.get(i);if(alt.isSymbolic()) continue;var lenAlt=alt.length(); if(lenRef&lt;lenAlt) return true; } return false; }accept(variant);' input.vcf
</code></pre>

<p>.</p>

CC BY-SA 3.0
2017-06-16T11:24:31.330Z
2017-06-16T11:24:31.330Z
Pierre
<p>If you're doing set operations, you could use <a href="http://bedops.readthedocs.io/en/latest/content/reference/file-management/conversion/vcf2bed.html" rel="nofollow noreferrer"><code>vcf2bed</code></a>:</p>

<pre><code>$ vcf2bed --insertions &lt; in.vcf &gt; out.bed
</code></pre>

<p>Based on the VCF v4.2 specification, <code>--snvs</code>, <code>--insertions</code>, and <code>--deletions</code> are options available to filter input. In each case, the length of the reference and alternate alleles is used to determine which type of variant is being handled.</p>

CC BY-SA 3.0
2017-06-16T23:03:39.673Z
2017-06-16T23:03:39.673Z
Alex Reynolds
<p>A one-liner:</p>

<pre class="lang-bash prettyprint-override"><code>zcat my.vcf.gz |
  perl -ane '$x=0;for $y (split(",",$F[4])){$x=1 if length($y)&gt;length($F[3])}print if /^#/||$x'
</code></pre>

<p>or equivalently</p>

<pre class="lang-bash prettyprint-override"><code>zcat my.vcf.gz |
  perl -ane '$x=0;map{$x=1 if length&gt;length($F[3])}split(",",$F[4]);print if /^#/||$x'
</code></pre>

<p>For simple VCF operations, I generally recommend to write a script. This may be faster than those using heavy libraries. With a script, you only parse fields you care about; most libraries unnecessarily parse every field.</p>

<hr>

<p>On a related note, I recommend <strong>not</strong> to decompose multi-allelic sites unless necessary. Decomposing is tricky, makes VCF harder to parse and to understand and may be a potential source of errors. Here is an example:</p>

<pre><code>#CHROM POS ID REF  ALT               QUAL FILTER  INFO  FORMAT  S1  S2  S3  S4
11     101 .  GCGT G,GCGA,GTGA,CCGT  199  PASS    .     GT      0/1 1/2 2/3 2/4
</code></pre>

<p>vt decompose+normalize produces the following VCF:</p>

<pre><code>#CHROM POS ID REF  ALT QUAL FILTER INFO FORMAT  S1   S2   S3   S4
11     101 .  GCGT G   199  PASS   .    GT      0/1  1/.  ./.  ./.
11     101 .  G    C   199  PASS   .    GT      0/.  ./.  ./.  ./1
11     102 .  CGT  TGA 199  PASS   .    GT      0/.  ./.  ./1  ./.
11     104 .  T    A   199  PASS   .    GT      0/.  ./1  1/.  1/.
</code></pre>

<p>In theory, you can reconstruct the original VCF from this output. However, it is very challenging for a program to do that. When you compute allele frequency line-by-line, this VCF will give you wrong results. <code>bcftools norm -m-</code> replaces "." with "0". You can get a correct ALT allele frequency from the bcftools output, but a wrong REF allele frequency. Furthermore, vt is also imperfect in that "CGT=>TGA" is not decomposed.</p>

<p>My preferred output is:</p>

<pre><code>#CHROM POS    ID     REF    ALT    QUAL   FILTER INFO   FORMAT S1     S2     S3     S4
11     101    .      GCGT   G,&lt;M&gt;  0      .      .      GT     0/1    1/2    2/2    2/2
11     101    .      G      C,&lt;M&gt;  0      .      .      GT     0/2    2/0    0/0    0/1
11     102    .      C      T,&lt;M&gt;  0      .      .      GT     0/2    2/0    0/1    0/0
11     104    .      T      A,&lt;M&gt;  0      .      .      GT     0/2    2/1    1/1    1/0
</code></pre>

<p>Here we use a symbolic allele <code>&lt;M&gt;</code> to represent "another ALT allele". You can calculate the allele frequency by looking at one line, and won't confuse other ALT alleles with REF. bgt can produce such a VCF indirectly. However, it discards all INFO, so is not a practical solution, either.</p>

<p>In summary, it is very difficult to decompose multi-allelic sites. When you get decomposing wrong, your downstream analyses may be inaccurate. Decomposition should be used with caution.</p>

CC BY-SA 3.0
2017-06-16T14:45:22.267Z
2017-06-16T14:45:22.267Z
user172818
<p>I'm looking to subset a standard VCF file to generate one which only includes insertions (i.e. not indels). </p>

<p>I can get part of the way there with:</p>

<pre><code>bcftools view -v indels &lt;vcf&gt; | awk '{if(length($4) == 1) print}'
</code></pre>

<p>However this wouldn't catch an insertion that was part of a multi-allelic record with an insertion and deletion, where the reference length would be greater than 1 bp. Potentially then one way to go is a chain of decomposition, normalisation and then the same reference length filtering — surely there's a better way in one of the many VCF manipulation utils?</p>

CC BY-SA 3.0
2017-06-16T10:59:34.627Z
What does your VCF look like? Many tools produce VCFs that include the type of variant.
terdon
CC BY-SA 3.0
2017-06-16T11:09:59.950Z
Safe to assume there's not a type annotation you can just grep out, though running an annotation tool first could be one answer.
blmoore
CC BY-SA 3.0
2017-06-16T11:14:02.000Z
Oh well, things are never that easy. Still, could you please [edit] your question and give an example of the kind of multi-allelic record you're thinking of? Ideally, give us a minimal example of your vcf file and include a couple of entries you want to keep and a couple you want to skip and show the desired output. That way we can use your example to test our solutions and we can be sure we're giving you what you need.
terdon
<vcf>
How can I extract only insertions from a VCF file?
# How can I extract only insertions from a VCF file?

(Asked by: username_0 on 2017-06-16)

I'm looking to subset a standard VCF file to generate one which only includes insertions (i.e. not indels). 
I can get part of the way there with:
bcftools view -v indels  | awk '{if(length($4) == 1) print}'

However this wouldn't catch an insertion that was part of a multi-allelic record with an insertion and deletion, where the reference length would be greater than 1 bp. Potentially then one way to go is a chain of decomposition, normalisation and then the same reference length filtering — surely there's a better way in one of the many VCF manipulation utils?

**Tags:** vcf
> **Comments:**
> **username_1**: What does your VCF look like? Many tools produce VCFs that include the type of variant.  
> **username_0**: Safe to assume there's not a type annotation you can just grep out, though running an annotation tool first could be one answer.  
> **username_1**: Oh well, things are never that easy. Still, could you please [edit] your question and give an example of the kind of multi-allelic record you're thinking of? Ideally, give us a minimal example of your vcf file and include a couple of entries you want to keep and a couple you want to skip and show the desired output. That way we can use your example to test our solutions and we can be sure we're giving you what you need.  

---

## 5 answers

### Answer by username_7 on 2017-06-16. Score: 4 | Accepted answer

A one-liner:
zcat my.vcf.gz |
  perl -ane '$x=0;for $y (split(",",$F[4])){$x=1 if length($y)>length($F[3])}print if /^#/||$x'

or equivalently
zcat my.vcf.gz |
  perl -ane '$x=0;map{$x=1 if length>length($F[3])}split(",",$F[4]);print if /^#/||$x'

For simple VCF operations, I generally recommend to write a script. This may be faster than those using heavy libraries. With a script, you only parse fields you care about; most libraries unnecessarily parse every field.

On a related note, I recommend not to decompose multi-allelic sites unless necessary. Decomposing is tricky, makes VCF harder to parse and to understand and may be a potential source of errors. Here is an example:
#CHROM POS ID REF  ALT               QUAL FILTER  INFO  FORMAT  S1  S2  S3  S4
11     101 .  GCGT G,GCGA,GTGA,CCGT  199  PASS    .     GT      0/1 1/2 2/3 2/4

vt decompose+normalize produces the following VCF:
#CHROM POS ID REF  ALT QUAL FILTER INFO FORMAT  S1   S2   S3   S4
11     101 .  GCGT G   199  PASS   .    GT      0/1  1/.  ./.  ./.
11     101 .  G    C   199  PASS   .    GT      0/.  ./.  ./.  ./1
11     102 .  CGT  TGA 199  PASS   .    GT      0/.  ./.  ./1  ./.
11     104 .  T    A   199  PASS   .    GT      0/.  ./1  1/.  1/.

In theory, you can reconstruct the original VCF from this output. However, it is very challenging for a program to do that. When you compute allele frequency line-by-line, this VCF will give you wrong results. bcftools norm -m- replaces "." with "0". You can get a correct ALT allele frequency from the bcftools output, but a wrong REF allele frequency. Furthermore, vt is also imperfect in that "CGT=>TGA" is not decomposed.
My preferred output is:
#CHROM POS    ID     REF    ALT    QUAL   FILTER INFO   FORMAT S1     S2     S3     S4
11     101    .      GCGT   G,  0      .      .      GT     0/1    1/2    2/2    2/2
11     101    .      G      C,  0      .      .      GT     0/2    2/0    0/0    0/1
11     102    .      C      T,  0      .      .      GT     0/2    2/0    0/1    0/0
11     104    .      T      A,  0      .      .      GT     0/2    2/1    1/1    1/0

Here we use a symbolic allele  to represent "another ALT allele". You can calculate the allele frequency by looking at one line, and won't confuse other ALT alleles with REF. bgt can produce such a VCF indirectly. However, it discards all INFO, so is not a practical solution, either.
In summary, it is very difficult to decompose multi-allelic sites. When you get decomposing wrong, your downstream analyses may be inaccurate. Decomposition should be used with caution.


---

### Answer by username_2 on 2018-06-29. Score: 5

Just to highlight that all the steps can be done within bcftools capabilities, and since I can't just comment on @blmoore 's answer: 
bcftools view --types indels  |
  bcftools norm -m - |
  bcftools filter --include 'strlen(REF)
more bult-in functions for 'bcftools filter' here

> **Comments:**
> **username_3**: Very nice, I had only ever used bcftools view to uncompress and then parse the output into either the header, or the data records.  

---

### Answer by username_0 on 2017-06-16. Score: 3

One method is to decompose multi-allelic records so that they're represented as one-allele, one-record using vt or similar:
bcftools view -v indels  |
  vt decompose - |
  bcftools view -H |
  awk '{if(length($5)>length($4)) print}'

Some decomposed alleles will come with excess reference padding. To left-shift and trim these, add a normalize step (NB matching these back to your input VCF becomes non-trivial):
bcftools view -v indels  |
  vt decompose - |
  vt normalize -r 
edit: As gringer suggests, this can also be done without vt:
bcftools view -Ou -v indels  |
  bcftools norm -Ou -Nm - |
  bcftools view -H |
  awk '{if(length($5)>length($4)) print}'

To also include complex alleles, use view -V snps (!snvs) instead of -v indels

> **Comments:**
> **username_4**: `bcftools norm` will also split multi-allelic records with the `-m -` option  

---

### Answer by username_5 on 2017-06-16. Score: 2

Using vcffilterjs

get the length of the REF;
loop over the ALT, ignore the symbolic
accept the variant if it's an insertion , eq:  len(ALT)>len(REF) 

.
java -jar dist/vcffilterjs.jar -e 'function accept(vc){var a=vc.getAlleles();var lenRef=a.get(0).length();for(i=1;i
.


---

### Answer by username_6 on 2017-06-16. Score: 0

If you're doing set operations, you could use vcf2bed:
$ vcf2bed --insertions < in.vcf > out.bed

Based on the VCF v4.2 specification, --snvs, --insertions, and --deletions are options available to filter input. In each case, the length of the reference and alternate alleles is used to determine which type of variant is being handled.


---

