Get genotype from reference fasta

Question

[win]

Hi all,

In NGS analysis we use a reference fasta and I wanted to know if it would be possible to get genotype at a certain genomic location from that file?

Answer 1

[Devon Ryan]

samtools faidx foo.fa chromosome:pos1-pos2

[Matt Shirley]

And if you need to do it in a script or use Windows

pip install --user pyfaidx 
faidx foo.fa chr:pos1-pos2

[win]

This command produced a binary sai file. So I am still unsure how I could get genotype at pos1-pos2?

[Devon Ryan]

It prints the sequence to the screen. It just has to build an index the first time you use it.

[win]

I had to run samtools faidx file.fa and then your command. Now I can see the output in FASTA format.

[Matt Shirley]

Yes, an unfortunate behavior of samtools is that you must invoke faidx once to generate the index, then again to fetch sequence. This is something I did not emulate with pyfaidx.

[Devon Ryan]

This must have already gotten fixed, since the version of samtools I have (granted, I track the github repo) will create the index (if needed) and then load it in a single go.

[Matt Shirley]

You might be right, or I might be remembering incorrectly :)