Tools to detect SNP supporting read count in each bam file

Question

[sckinta]

Hi,

I have 45 pool-seq samples and aligned to reference (bam file). Those 45 pool-seq belong to 3 groups, say group 1, group 2 and group 3, thus each group have 15 samples. I combined 45 pool-seq samples and change the RGSM to Group number, ran freebayes to detect SNP. The vcf gives the result of genotype for each group, but not for each sample.

Now with given SNP information (chr and pos), I have to detect read count number supporting for each SNP in each pool-seq samples instead of group.

e.g.,

Given variant at chr1:2000, major allele is A and minor is G. I want to know how many A and how many G in pool-seq Sample_1 respectively at that position.

Any tools can do that?

Thank you

Answer 1

[ernfrid]

I have not run freebayes much myself, but looking at its code, it seems that there should be RO and AO tags in the FORMAT columns for each pool that answer your question.

Re-running freebayes is likely to be your most comprehensive bet, especially if you become interested in indels. That said, you can get counts by either parsing samtools mpileup output or using bam-readcount (https://github.com/genome/bam-readcount). Fair disclosure, I am the author of the latter tool.

[sckinta]

Yes, it has RO and AO for each sample group (each containing 15 individual pool-seq samples), but not for individual pool-seq samples/libraries. I can rerun freebayes on reheaded each library RGSM to distinguished names (instead of sample group name), but freebayes itself takes really really long time to finish.

[ernfrid]

I understand now. Answer edited with two possible methods.

[sckinta]

Thank you! I will try it. I wrote a script too, but it is super slow.