Samtools Mpileup Defined Regions

Question

[Vova Naumov]

Hi! I need to get genotypes in a list of positions (for example BED file) having BAM file. Trying to use tool samtools mpileup. I know I can do it by creating consensus fastq using

samtools mpileup -uf hg19.fa file.bam | \
bcftools view -cg - | \
perl vcfutils.pl vcf2fq > consensus.fastq

and then getting genotypes from fastq file, but it's ineffective. Do you know how to get differences from reference, genotype calls at interesting positions simultaneously?

[Zev.Kronenberg]

Please take a closer look at mpileup in samtools. The documentation clearly states how to call variants and use a bed file.

[Vova Naumov]

Thank's a lot but an example would be very helpful

Answer 1

[swbarnes2]

I'd use BEDTools intersectBed to get the intersection between your all-points vcf and your bed file coordinates. That's way simpler than making a consensus fastq.