Hello, I’m trying to figure out why the number of reads in the STAR log and the RSEM log don’t seem to tally. If my understanding is right, based on the logs below, I’m seeing the following results from my nf-core/rnaseq run 1. Total Reads a. STAR: 14227877 b. RSEM: 10444921 2. Uniquely mapped reads a. STAR: 9958596 b. RSEM: 8618444 3. Multi-mapped reads a. STAR: 3183037 b. RSEM: 740233 Below are the raw outputs for reference: According to RSEM documentation,

N0 N1 N2 N_tot   
#  N0, number of unalignable reads; N1, number of alignable reads; N2, number of filtered reads due to too many alignments; N_tot = N0 + N1 + N2     
    
nUnique nMulti nUncertain   
# nUnique, number of reads aligned uniquely to a gene; nMulti, number of reads aligned to multiple genes; nUnique + nMulti = N1;                             
# nUncertain, number of reads aligned to multiple locations in the given reference sequences, which include isoform-level multi-mapping reads

From the RSEM log:

> head sample1.stat/sample1.cnt
1086244 9358677 0 10444921
8618444 740233 6960702
39948491 3
0	1086244
1	2397975
2	1357491
3	1225020
4	1145903
5	777761
6	626725

From the STAR log:

❯ cat log/sample1.log
                                 Started job on |	Sep 03 02:24:13
                             Started mapping on |	Sep 03 02:24:41
                                    Finished on |	Sep 03 02:37:22
       Mapping speed, Million of reads per hour |	67.31

                          Number of input reads |	14227877
                      Average input read length |	239
                                    UNIQUE READS:
                   Uniquely mapped reads number |	9958596
                        Uniquely mapped reads % |	69.99%
                          Average mapped length |	234.02
                       Number of splices: Total |	10169234
            Number of splices: Annotated (sjdb) |	9979952
                       Number of splices: GT/AG |	10051968
                       Number of splices: GC/AG |	68370
                       Number of splices: AT/AC |	7487
               Number of splices: Non-canonical |	41409
                      Mismatch rate per base, % |	0.37%
                         Deletion rate per base |	0.03%
                        Deletion average length |	2.07
                        Insertion rate per base |	0.00%
                       Insertion average length |	2.46
                             MULTI-MAPPING READS:
        Number of reads mapped to multiple loci |	3183037
             % of reads mapped to multiple loci |	22.37%
        Number of reads mapped to too many loci |	28976
             % of reads mapped to too many loci |	0.20%
                                  UNMAPPED READS:
  Number of reads unmapped: too many mismatches |	214036
       % of reads unmapped: too many mismatches |	1.50%
            Number of reads unmapped: too short |	779433
                 % of reads unmapped: too short |	5.48%
                Number of reads unmapped: other |	63799
                     % of reads unmapped: other |	0.45%
                                  CHIMERIC READS:
                       Number of chimeric reads |	0
                            % of chimeric reads |	0.00%

Would greatly appreciate any insights into what might be causing these differences. Thanks!

Excellent question, unfortunately I don't know, either, and I am not very familiar with RSEM. I looked into this a bit, and here are my observations: • The "unmapped reads" agree between the two (if for STAR you include the

reads mapped to too many loci

). • The RSEM reads are in total 3,782,956 fewer than STAR. • All of these are differences in "alignable" reads, i.e. if you add up the differences in "unique" and "multimapped" (

nMulti

in RSEM), you'll get to this number as well. • RSEM has additionally 6,960,702 in

nUncertain, number of reads aligned to multiple locations in the given reference sequences, which include isoform-level multi-mapping reads

. I think this is the key. I am not convinced these are "reads" (as in entries in the FASTQ), because if you add them to

nUnique

and

nMulti

, the total becomes much too high (even higher than total STAR reads). Maybe these are actually alignments (where one read can have multiple)? In any case, my best guess is this is how RSEM accounts for the reads missing from whatever STAR counts. I hope someone else knows more, ohterwise, if you really, really want to know, the RSEM source code should have the answer ...