2
$\begingroup$

I am working on a sample which was sequenced using an amplicon kit and I would like to ensure that no variants are called in the PCR primer sequences. I know that the primer sequences are part of the genome and a variant there is no less valid than a variant anywhere else, but for reasons I cannot control, I need to ensure that no variants can be called in the primer regions.

My first thought was to use a soft-clipping tool like BAMClipper. However, the variant caller I am using (Sentieon's implementation of GATK) preforms a local realignment when calling and if it finds that soft-clipped regions align perfectly against the genome, it will remove the soft clips. This makes perfect sense for cases where the softclips were added by the aligner, but sadly breaks my workflow since I need the softclipped primers to remain clipped.

I would therefore like to hard clip, to completely remove the primer sequences from the bam file before variant calling. I have found one tool that seems to be able to do this, TrimPrimers by Fulcrum Genomics, but this fails on my input bam file with an error I can't really make much of.

Is there any other tool that can take a list of regions (e.g. a bedpe file) and a BAM and hard-clip the regions from the BAM?

Improve this question
$\endgroup$

2 Answers 2

Reset to default
2
$\begingroup$

Another option would be to adjust the quality of the primer bases to be below the threshold that your variant caller uses. In my case I have "masked" the primer bases by dropping the quality to 0. This prevents the variants from being called based that are synthetic but if you have an overlapping amplicon it will call bases from those reads.

Improve this answer
$\endgroup$
2
  • $\begingroup$ That sounds promising. How would I go about that given a bedpe file of primer locations and a bam file? Is there a tool that can take those two files as input and set the MAPQ to 0? Or do you mean changing the quality score in the original fastq file? $\endgroup$
    – terdon ♦
    Commented Jul 23, 2019 at 8:22
  • $\begingroup$ The tools we use are custom built but in general we detect the primer for each amplicon and then edit the synthetic (primer) baseQs not the MapQ as 0 for those reads. You could do it either in the FASTQs prior to alignment or on the BAMs before variant calling. If you don't have overlapping amplicons this would probably be overkill and filtering the result VCF or redefine the regions of interested to the variant caller might be easier $\endgroup$
    – Bioathlete
    Commented Jul 24, 2019 at 14:28
1
$\begingroup$

I know that the primer sequences are part of the genome and a variant there is no less valid than a variant anywhere else,

You won't see variants under PCR primers. You get the primer sequence itself in the vast majority of reads.

I wouldn't clip the bam file, I'd clip either the fastq, or clip within bcl2fastq while making the fastqs.

Improve this answer
$\endgroup$
7
  • $\begingroup$ Thanks, but I do actually get variants in the primer. That's the problem :( And I really don't want to clip the fastq. Both because the primers can help the alignment process and because that's harder to do because I would require fasta sequences of the primers and not regions which is what I have and what assay manufacturers provide. Finally, I never have the bcl files, I will always start from FASTQ (people send me their fastq files). $\endgroup$
    – terdon ♦
    Commented Jul 22, 2019 at 19:56
  • $\begingroup$ Variants under the primers can't be real. Just look at how PCR works. You will only see a variant under the forward primer if amplification happens from the reverse primer sitting on the orignal template. Once you get amplification happening off of templates made by extending off of the forward primer, the variant is gone. The orignal template will be outnumbered a billion to 1. $\endgroup$
    – swbarnes2
    Commented Jul 22, 2019 at 20:13
  • $\begingroup$ @swbarnes2 primers may be inaccurately synthetized and lead to false positives, theoretically, but I would be more worried about false negatives. A true variant under the primer will be lost, for which you should not call a 0/0 but a ./. $\endgroup$
    – Wouter De Coster
    Commented Jul 23, 2019 at 5:50
  • $\begingroup$ @swbarnes2 that's a very valid point. But my main issue is that I have a client who doesn't want to see any variants in the primers, ever. And I am calling variants in th primers, for whatever reason. It's possible that the primer itself has a mismatch, I guess, and that's what I'm seeing. But I need to completely hide these regions from the variant caller. $\endgroup$
    – terdon ♦
    Commented Jul 23, 2019 at 8:09
  • $\begingroup$ Wouldn't it be easier to filter the vcf? $\endgroup$
    – swbarnes2
    Commented Jul 23, 2019 at 15:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.