This is a test version of Biostars. For the public version, visit https://www.biostars.org.
bam_stat.py and featureCounts stats don't seem to make sense (to me)

Hello,

Tl;dr - can low quality (non-unique) alignments (according to bam_stat.py) be considered multimappers by featureCounts?

I have an RNAseq bam file (subsetted for troubleshooting) that when I run featureCounts, it tells me I have multimapping reads. However, this bam file has been filtered for primary reads only (-F 256), and my bam_stat output doesn't indicate any multimapped reads:

bam_stat.py -i UWNsub100k_256_noribos.bam

Total records:                          103011

QC failed:                              0
Optical/PCR duplicate:                  0
Non primary hits                        0
Unmapped reads:                         54891
mapq < mapq_cut (non-unique):           11891

mapq >= mapq_cut (unique):              36229
Read-1:                                 18198
Read-2:                                 18031
Reads map to '+':                       18009
Reads map to '-':                       18220
Non-splice reads:                       36017
Splice reads:                           212
Reads mapped in proper pairs:           32780
Proper-paired reads map to different chrom:0

From the above, it looks like I have 11891 non-unique reads.

When I run featureCounts, unmapped reads and assigned reads appear to agree with the bam_stat output, but featureCounts seems to by considering my non-unique reads as multimappers. Should it not be categorizing them as Unassigned_MappingQuality? And shouldn't multimapping reads be filtered out with -F 256?

featureCounts -T 4 -p -a ~/path/sacCer3.ncbiRefSeq.gtf -t exon -g gene_id -o ~/path/UWNsub100temp.txt UWNsub100k_256_noribos.bam

Assigned            14434
Unassigned_Unmapped     27144
Unassigned_Read_Type        0
Unassigned_Singleton        0
Unassigned_MappingQuality   0
Unassigned_Chimera      0
Unassigned_FragmentLength   0
Unassigned_Duplicate        0
Unassigned_MultiMapping     6034
Unassigned_Secondary        0
Unassigned_NonSplit     0
Unassigned_NoFeatures       3657
Unassigned_Overlapping_Length   0
Unassigned_Ambiguity        281

Any insights are much appreciated - thank you!

rna-seq alignment

0 answers

No answers yet.

Log in to answer this question.