Obtain contex sequence from a vcf file

Question

[shinken123]

Hi

I have a vcf file that looks like this:

4   169239903   M1  G   A   .   PASS    .   GT  0/0 1/1 1/1 1/1 1/1 1/1 1/1 1/1 ./. 1/1 1/1 1/1 1/1 1/1 ./. 1/1 1/0 1/1 1/1 ./.
4   169832919   M2  A   C   .   PASS    .   GT  0/0 ./. 1/1 1/1 ./. 1/1 1/1 1/1 1/1 1/1 1/1 1/1 1/1 ./. 1/1 ./. ./. 1/1 1/1 ./.
4   169989920   M3  C   T   .   PASS    .   GT  0/0 ./. ./. 1/1 1/1 1/1 1/1 1/1 1/1 1/1 1/1 1/1 1/1 1/1 ./. 1/1 ./. 1/1 ./. ./.

Now I will like to obtain a contex sequence with around 15 nucleotides upstream the SNP and 15 nucleotides downstream the SNP from the reference genome to have something like this:

>4_169239903
ATCCCGGTTCCCTCTgTTTCCTTTCCTTGCT
>4_169832919
ATCGTAGTTAGCTAAaTTTCGTACTGACTGA
>4_169989920
CGTGCTGCATGATCCcTTGGCAGTACGTACC

Do you known how to obtain it?

Thank you

Answer 1

[Pierre Lindenbaum]

• convert your vcf to bed vcf to bed
• extends the bed with bedtools slop http://bedtools.readthedocs.io/en/latest/content/tools/slop.html
• extract the fasta sequence with bedtools getfasta http://bedtools.readthedocs.io/en/latest/content/tools/getfasta.html

[shinken123]

Thank you very much Pierre,

This is exactly what I want.