No need to sort in the second command.
• 0 views
•
link
written
by
GouthamAtla
Why the file unmapped.bam, from tophat has all reads with mapq =/> 30?
Used the samstat tool (http://samstat.sourceforge.net/) for seeing the mapq scores of my alignment;
First, I used with the accepted_hits.bam, and found some reads classified as unmapped. And, I thought strange. Shouldn't the accepted hits contain only mapped reads?
Second, When I use the samstat against unmapped reads.bam, all reads are 30 + mapq scored.
Can someone explain this to me?
• 4,463 views
•
link
updated
by
Ram
•
written
by
tiago211287
•
written
by
tiago211287
1 answer
As the following command gave you 756288 255, this indicates all the reads in unmapped.bam have mapping quality of 255, which indicates that the mapping quality could not be assigned to them, hence they are all can be considered as unmapped. If you run the same command on accepted_hits.bam, you should see MAPQ of 0,1,3,50.
samtools view unmapped.bam | cut -f5 | sort | uniq -c
You can also run the following command and see if all the reads have reference name as *. This also indicates the reads mapped to no chromosomes.
samtools view unmapped.bam | cut -f3 | sort | uniq -c
• 0 views
•
link
updated
by
Ram
•
written
by
GouthamAtla
•
written
by
GouthamAtla
Can you tell me how can be possible that, all MAPQ in all bases of all reads be the same? like this plot on the samstat are saying?
Saw here that "The counts and proportions should be almost invariant across read positions"
• 0 views
•
link
updated
by
Ram
•
written
by
tiago211287
•
written
by
tiago211287
More one thing, shouldn't the accepted_hits.bam have 0 unmapped reads? Must have something wrong with my alignment. as you can see here:
There are 0.5% unmapped reads on the accepted_hits.bam: http://s28.postimg.org/wunusnl71/Untitled.png)
• 0 views
•
link
updated
by
Ram
•
written
by
tiago211287
•
written
by
tiago211287
It should have 0 unmapped reads. If you run the command mentioned:
samtools view accepted_hits.bam | cut -f5 | sort | uniq -c
you should not see any reads with MAPQ 255. I never used any software for these calculations. I use only samtools with various flags.
• 0 views
•
link
written
by
GouthamAtla
MAPQ is not for all bases. its for read. All unmapped reads in tophat output will have MAPQ of 255, which indicates that the MAPQ can not be calculated for them.
• 0 views
•
link
written
by
GouthamAtla
Log in to answer this question.
More posts like this
-
Unexpected results when alignment FASTQ files with different read alignments
written by ManuelDB
-
how to identify where multi-mappers, or "non primary hits" are mapping to
written by bertb •
-
Filtering BAM file to include PE reads that did not map, even though their mate did
written by nute11a •
-
Mismatch number of reads after bedtools bam2fastq on RNA-seq data
written by aditi.qamra
-
Extracting mapped paired end FastQ reads from BAM file using Samtools - trying to understand -f 0x0…
written by acxcrew007 •
-
Same alignment but with different output on tophat2, any ideas?
written by tiago211287
-
Is correct finding unique reads with the flag -F 0x104?
written by miquinhap •
-
What Is A Good Strategy For Finding Common Bacterial Sequences In A Potentially Contaminated Human …
written by Dan D
-
Paired-End Sam Flag Mismatch During Filtering
written by Andreas
-
Trying To Calculate Percent Of Mapped And Unmapped Vs.Total Reads
written by erichpeterson •
Can you post the output of
I found strange that my accepted_hits has unmapped reads. Also, The unmapped.bam has 30+ mapq.