Hi, I try to use samtools in a pipeline for variant calling from RNA-seq reads. To improve the quality of the calls, I preprocess the mapped reads with a tool [1] that splits reads containing a skipped region (N in cigar string) and removes parts of the alignment which are very short, e.g., from 93M1394N7M the 7M is removed by changing the alignment to 93M7H. When calling variants on that edited file, samtools mpileup misses a large number of SNPs which should at least be reported. I narrowed it down to the calculation of the BAQ so I read your papers about it but still do not understand why the BAQ for the bases is lowered so much. Please find attached SAM files showing the original reads and the “splitted reads”. I added a visualization from the IGV in which I marked the questionable SNV by a red arrow (pos chr12:85277730). Note that I have chosen this example because of the low number of reads and that examples with a higher coverage exist.
The pileups using different parameters display a low BAQ for the mutation after
the splitting but not in the original files:
$ samtools mpileup -f genome.fa -r chr12:85277730-85277730 original.bam
chr12 85277730 C 15 ,......,...,TTT
FHEJJIIDJHGCDFF
$ samtools mpileup -f genome.fa -r chr12:85277730-85277730 splitted.bam
chr12 85277730 C 12 ,......,...,
FHEJJIIDJHGC
$ samtools mpileup -Q 0 -f genome.fa -r chr12:85277730-85277730 splitted.bam
chr12 85277730 C 15 ,......,...,TTT
FHEJJIIDJHGC!!!
$ samtools mpileup -B -f genome.fa -r chr12:85277730-85277730 splitted.bam
chr12 85277730 C 15 ,......,...,TTT
FHEJJIIDJHGCDFF
The behavior is reproducible with samtools 0.1.19 and 1.2
To conclude, I do not understand this severe drop of quality just because the
last 7 bases of the reads are trimmed. Can you probably explain this or do you
recommend to not using the BAQ?
Thanks for your help in advance
Manuel
[1] https://github.com/mjafin/splitNreads (results equal GATKs SplitNCigarReads)
[cid:[email protected]]
Dr. Manuel Landesfeind
Bioinformatician
Oncotest GmbH, Am Flughafen 12-14, 79108 Freiburg, Germany
phone: +49 (0) 761 51559 591 ▪ fax: +49 (0) 761 51559 55
email: [email protected]<mailto:[email protected]>
website: www.oncotest.com<http://www.oncotest.com/>
This message and any attachments are confidential and may be privileged or
otherwise protected from disclosure. lf you are not the intended recipient, you
must not copy this message or attachment or disclose the contents to any other
person. lf you have received this transmission in error, please notify the
sender immediately and delete the message and any attachment from your system.
original.sam
Description: original.sam
splitted.sam
Description: splitted.sam
------------------------------------------------------------------------------ Don't Limit Your Business. Reach for the Cloud. GigeNET's Cloud Solutions provide you with the tools and support that you need to offload your IT needs and focus on growing your business. Configured For All Businesses. Start Your Cloud Today. https://www.gigenetcloud.com/
_______________________________________________ Samtools-help mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/samtools-help
