Hi,
I am trying to extract the coverage and the average quality score for each position of a reference assembly in bam/sam format. I have managed to get the coverage using BEDtools
genomeCoverageBed -ibam mybamfile.bam -g my_genome -d > my_coverage.txt
but am at a loss on how to get some measure of the quality of the base calls at each position. I was thinking that I could use the bcftools to get a variant call formatted file
samtools mpileup -uf ref.fa mybamfile.bam | bcftools view -bvcg - > var.raw.bcf
bcftools view var.raw.bcf | vcfutils.pl varFilter -D100 > var.flt.vcf
but this only provides the sites for which there are SNPs. Any advice greatly appreciated.
Joseph
•
written
by
Joseph Hughes
Don't use use bcftools with the "calling options" ( -v , -c -g ). You should get an output for all the positions
Hi Pierre, I have tried that, there are more positions but not all of them. I am going to try setting -d10000000 in the mpileup command to see if that changes anything. Thanks
Log in to answer this question.
