Hi
I'm trying to make GATK UnifiedGenotyper to detect indels inserted by wgsim short reads simulator.
Broad employees just informed me that indels are not called because wgsim creates reads of Q17 and everything below 20 is thrown away.
I checked the fastq files, and indeed, all qualities are "2" which is if I understood well the ASCII translation corresponds to 17.
Is there a way to make wgsim output reads with a configurable base quality?
I thought that a trick could be to replace these "2222...2" strings by somthing like "????...?" Or is there a more elegant way to proceed?
Thanks.
3 answers
You could use the following awk script to change your 2s to Is:
awk ' /222222222/ { gsub("2", "I"); print $0; next } { print } '
Wait a minute: I actually used the wgsim available at https://github.com/lh3/wgsim and this produces only "2" qualities ?!
my mistake - I have a file like that I thought I generated it with wgsim - you are correct it does put in 2s. I will correct the answer
I forked wgsim and made the quality score configurable using the -Q command line argument. I've made I the default value for the quality score (Phred quality score of 40). Feel free to use my fork at https://github.com/hammer/wgsim to alleviate the need for the post-processing step with awk.
Previously the quality score was scaled with the error rate, so that as the error rate went up, the quality score went down. The precise formula used is (-10.0 * log(ERR_RATE) / log(10.0) + 0.499) + 33 (see https://github.com/lh3/wgsim/blob/master/wgsim.c#L248). Because ERR_RATE is set to 0.02 by default, the quality score comes out to (-10.0 * log(0.02) / log(10.0) + 0.499) + 33 = 50.4887..., which is encoded using the character 2, which has the ASCII character 50.
To make what Jeff said easier, this table gives you error probabilities for each quality score. If you set the -e flag in wgsim to the value in the P_error column you will get the respective base quality. If you want a given error rate with a different quality then you are better off replacing the base qualities like Istvan said.
Log in to answer this question.
