Hi everyone!
I am running an RNA-Seq experiment and was wondering if anyone could help me with a particular query i had.
So my RNA-Seq sample had 29,775,906 reads.
I ran this sample through tophat (no coverage search, unstranded, w/bowtie2). I checked the logs folder in the output and got the following stats :-
[2013-04-04 10:57:02] Preparing reads
left reads: min. length=50, max. length=50, 29767088 kept reads (8818 discarded)
Now, 29767088 + 8818 = 29,775,906 . So until now it makes sense.
From the bowtie files i got the following information :-
In,bowtie.left_kept_reads.log
29767088 (100.00%) were unpaired; of these:
4316521 (14.50%) aligned 0 times
17118264 (57.51%) aligned exactly 1 time
8332303 (27.99%) aligned >1 times
In,bowtie.left_kept_reads_seg1.log
4316521 (100.00%) were unpaired; of these:
3198237 (74.09%) aligned 0 times
1037571 (24.04%) aligned exactly 1 time
80713 (1.87%) aligned >1 times
In bowtie.left_kept_reads_seg2.log
4316521 (100.00%) were unpaired; of these:
3177510 (73.61%) aligned 0 times
1069186 (24.77%) aligned exactly 1 time
69825 (1.62%) aligned >1 times
So at this point , correct me if i am wrong, we can conclude :-
total reads aligned 1 time = 17118264 + 1069186 + 1037571 = 19225021 total reads aligned >1 time = 8332303 + 80713 + 69825 = 8482841
But then i converted the accepted_hits.bam file in a .sam file. On this .sam file I did the following two operations :-
awk '{print $1}' 11129_nc_u.sam | sort | uniq -u | wc -l
gave me 21638486 (i.e total number of uniques)
and
awk '{print $1}' 11129_nc_u.sam | sort | uniq -d | wc -l
gave me 5519946 (i.e. total duplicates.)
Now why is the number of uniques from the .sam folder (21638486 ) greater than 19225021 (i.e. obtained from the logs folder)? In the same breath, why is the number of duplicates from .sam folder (5519946) less than the number in the logs folder (8482841) ?
Can someone please clarify this point for me?
Thanks in Advance :)
