Question: How Do You Get The Quality Score And Coverage For Every Single Position Of A Reference Assembly
gravatar for Joseph Hughes
9.1 years ago by
Joseph Hughes2.9k
Scotland, UK
Joseph Hughes2.9k wrote:


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 | varFilter -D100 > var.flt.vcf

but this only provides the sites for which there are SNPs. Any advice greatly appreciated.


ADD COMMENTlink written 9.1 years ago by Joseph Hughes2.9k
gravatar for Pierre Lindenbaum
9.1 years ago by
France/Nantes/Institut du Thorax - INSERM UMR1087
Pierre Lindenbaum134k wrote:

Don't use use bcftools with the "calling options" ( -v , -c -g ). You should get an output for all the positions

ADD COMMENTlink written 9.1 years ago by Pierre Lindenbaum134k

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

ADD REPLYlink written 9.1 years ago by Joseph Hughes2.9k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1673 users visited in the last hour