multi-mapped reads counting methods with RNA-seq

Question

[guillaume.rbt]

Hi all,

I'm currently doing RNA-seq data analysis to assess differential expression between genes.

If I understand well, standard methods for counting reads does not take into account reads mapped to multiple loci, which doesn't allow to evaluate expression of duplicated genes.

I saw that there is a new tool that deals with multiple alignments (mmquant).

I was wondering if I still use a "classical" counting tool for uniquely mapped reads, like featurecounts, and then add the results of mmquant for the multiple alignments.

Does that make any sense, or is it statistically stupid and will mess with the normalization later?

Thanks,

Guillaume

[Devon Ryan]

Don't mix and match methods for different genes/types of reads, that's a recipe for messing up the downstream analysis.

[guillaume.rbt]

Thank for the tip, I had a bad feeling about this idea.

Answer 1

[h.mon]

Use just one tool, there is no reason to complicate matters. mmquant should give identical counts to uniquely mapped genes as featureCounts or HTSeq - the difference is only for multi-mappers. As the paper states:

mmquant is a drop-in replacement of the widely used tools htseq-count and featureCounts that handles multi-mapping reads in an unabiased way.

[Devon Ryan]

The method used by mmquant is rather odd. It's incredibly unclear what effect including counts for non-existent fusion genes (due to their sharing reads) would have downstream when one is performing differential expression. My primary concern is that this would both lead to an inflated number of tests (often with lower power). Any 3' or 5' bias in a sample could have dramatic consequences on this.

In general, I would personally think it's better to use RSEM/salmon/kallisto followed by tximport in R to get more accurate handling of multimappers on the gene level.

[guillaume.rbt]

Does those tools handle multi-mapping reads? Which may be the best for this task?

[h.mon]

All of them handle multi-mappers using an expectation-maximization algorithm. Salmon and kallisto are the fastest. However, they are not drop-ins replacement to featureCounts or HTSeq, as they do not use as input a bam file of reads aligned to a reference genome. Instead, Salmon and RSEM can count from a bam aligned to a reference transcriptome; and Salmon and kallisto can count reads directly to a reference transcriptome, without the need to align - which is the fastest option.