I am fairly new to variant calling and somehow I cannot find an answer to my question. Maybe I use the wrong the search terms...
I have a VCF file from freebayes and I want to know the average base quality of the bases that support that certain variant.
So for example I have this:
POS ID REF ALT
200 . T C
From all the C bases that are mapped on position 200 I want to calculate the average base quality score that you normally see in the FASTQ file (And QUAL column of SAM file). At this moment I am not sure if a tool to do this already exist or that I need create something myself.
EDIT:
Started with making my own script using mpileup but ran into another problem. In short, the length of the base column is not the same as the length of the quality score column. After some filtering I got a position with 101 bases and 100 quality scores. In the bases I noticed this (example):
.^.^A^,^,^,^,g^,^,^,^,^W,
One of my filter step is:
re.sub(r'\^[^actgnNACTG\.\,]', '', bases)
bases = bases.replace("^", "")
But with this filter step I remove one base to many. How can I know the character after the ^ is a mismatch/match or a ASCII score.
0 answers
No answers yet.
Log in to answer this question.
-
Get average quality score per position in an aligned bam file
written by cedric.blais •
-
Calling zero mapping quality variant
written by sonsunjirachote •
-
What is the relationship between the per-base sequence quality of a FASTQ file and the QUAL of a VC…
written by ManuelDB
-
trimming prinseq trimmomatic
written by Rob •
-
Fastqc of the fastq file
written by sgadila •
-
Filtering low quality score bases from BAM file
written by tayyabaalvi82 •
-
what is the conventional apporach to calculate depth of coverage?
written by serpalma.v •
-
samtools pileup format: insertions and deletions' quality
written by jing.mengrabbit •
-
Changing samtools mpileup I16 calculations?
written by emilywynn6 •
-
Need an advice on quality control of illumina reads
written by sumudu_rangika •
Hello gb ,
not an answer to your question, but why do you want to calculate the average base quality?
fin swimmer
I know there is a
--min-base-qualityoption but I want to have that average to be able to filter afterwards. I also have variant calling output from clc (not generated by me) that does contain this average value and I want (need) to compare it with freebayes.EDIT2:
Maybe some other background information, as I am reading and learning more about variant calling. The QUAL value seems to be important, but we sequence very deep and also want to see variants with a frequency of 2%. As far as I know now the depth heavily influence that score. So I need other values to determine the quality of that variant position.
Instead of second-guessing FreeBayes, I'd look for a variant caller that specializes in finding low % variants, and just trust that.