text: user148: I'm interested in obtaining coding sequences of my favourite gene in all individuals from the 1000Genomes (and similar projects). I use GATK to get the right subset of variants, vcf-consensus to map these variants onto the reference genome and finally samtools to extract the individual exons. This works fine if the variants are SNPs but if there are any indels, this changes the coordinates of exons and I end up getting the wrong region. Is there any generic way of remapping genomic coordinates to account for the changes created by indels?
  
  user1310: If the goal is to extract the consensus sequence for given regions, this is how it works in 2019.
  bgzip and index your vcf file.
  $ bgzip -c input.vcf > input.vcf.gz
  $ tabix input.vcf.gz
  
  Create a regions.txt which contains one region per line in the format chr:from-to. If you have already a bed file, you can use this little awkscript to create it:
  $ awk '{print $1":"$2+1"-"$3}' input.bed > regions.txt
  
  Now run this combination of samtools and bcftools:
  $ samtools faidx -r regions.txt genome.fa|bcftools consensus input.vcf.gz -o consensus.fa
  
  
  user425: I think that you need a LiftOver Chain file to transform your coordinates. You can obtain such a file using bcftools consensus with the -c parameter:
  -c, --chain <file>         write a chain file for liftover
  
  Then you can use it to transform coordinates in various genomic formats using CrossMap.
  
original_text: 148: <p>I'm interested in obtaining coding sequences of my favourite gene in all individuals from the 1000Genomes (and similar projects). I use GATK to get the right subset of variants, vcf-consensus to map these variants onto the reference genome and finally samtools to extract the individual exons. This works fine if the variants are SNPs but if there are any indels, this changes the coordinates of exons and I end up getting the wrong region. Is there any generic way of remapping genomic coordinates to account for the changes created by indels?</p>
   1310: <p>If the goal is to extract the consensus sequence for given regions, this is how it works in 2019.</p>
  
  <p><code>bgzip</code> and index your <code>vcf</code> file.</p>
  
  <pre><code><span class="math-container">$ bgzip -c input.vcf &gt; input.vcf.gz
  $</span> tabix input.vcf.gz
  </code></pre>
  
  <p>Create a <code>regions.txt</code> which contains one region per line in the format <code>chr:from-to</code>. If you have already a <code>bed</code> file, you can use this little <code>awk</code>script to create it:</p>
  
  <pre><code><span class="math-container">$ awk '{print $1":"$2+1"-"$</span>3}' input.bed &gt; regions.txt
  </code></pre>
  
  <p>Now run this combination of <code>samtools</code> and <code>bcftools</code>:</p>
  
  <pre><code>$ samtools faidx -r regions.txt genome.fa|bcftools consensus input.vcf.gz -o consensus.fa
  </code></pre>
   425: <p>I think that you need a LiftOver Chain file to transform your coordinates. You can obtain such a file using <code>bcftools consensus</code> with the <code>-c</code> parameter:</p>
  
  <pre><code>-c, --chain &lt;file&gt;         write a chain file for liftover
  </code></pre>
  
  <p>Then you can use it to transform coordinates in various genomic formats using <a href="http://crossmap.sourceforge.net/" rel="nofollow noreferrer">CrossMap</a>.</p>
  