Supplementary aligments in VAF

Question

This question has also been asked on Biostars

I have a doubt, are supplementary alignment usually considered when the variant allele frequency is calculated? Thanks a lot.

In my case I have some Illumina short data edited with Crisp-Cas9 and I want to calculate for a specific position all the frequency of the 4 nucleotides. I would like to understand if the supplementary reads aligned in that position must be taken into consideration.

I've used pysamstats to calculate and it consider them, but my doubt is more if it right to consider them, if it is normally done or not. We are doing it on human and we expect only a change of a single nucleotide variant in the region if the editing works.

I'm interested in what is done generally, if there is a "rule" to consider or not the supplementary alignment when the AF is calculated. In my case I have some Illumina short data edited with Crisp-Cas9 and I want to calculate for a specific position all the frequency of the 4 nucleotides. I was trying to understand if the supplementary reads aligned in that position must be taken into consideration.

GATK HaplotypeCaller and Mutect2 take them into consideration? I'm not looking for a specific variant caller, but if I want to take one variant caller as "reference", I would probably look at them.

Thanks a lot.

Answer 1

It depends on the tool and the use case but generally they should be.

It's not as important with short reads as only a small portion of reads have supplementary mapping (unless you're looking at breast cancer post chromothripsis, for example). But imagine the case where you're working with long reads (as I do). I have samples with a mean read length of 60kb. The human reference is structurally different to any given sample meaning that the chance of any given long read having supplementary mapping is high. If supplementary mapping didn't count towards allele counts then many would be artificially low. This is my use case, others are different.

I would count them unless you have some reason to suspect that they are are artefact of the editing and or the mapping. Sorry to never give a straight answer but that's bioinformatics! From what you described I would guess 99% are real but there's normally some noise and if that's of concern your code should include some QC. I.e, consider if multiple reads support the sup mapping and their quality.