Why is eXpress producing strangely low counts?

Question

[bagi.m]

I have aligned our RNA-seq to reference transcriptome (there is no genome) using bwa mem. The reference transcriptome is from the same species but from different population, so some mismatches and indels are expected.

Then I used eXpress to count reads.

The results are strangely low. No reference contig has more than 60 tot_counts. However, when I look at the alignment in tablet or IGV some contigs have up to 500000 reads aligned to them. Moreover, uniq_counts is equal to tot_counts for all contigs. But cca 0.1% of reads in my bam file are multimapped.

Any suggestions what could be wrong?

Answer 1

[bagi.m]

Ok, this one is embarrassing. I completely misunderstood the eXpress documentation regarding sorting the input .bam file and the fact that you are _not_ supposed to sort it.

Answer 2

[Istvan Albert]

Use kallisto instead. The creators of eXpress recommend as a better alternative.

As for what is going on - well I will say that it is difficult to track down not know what the data looks like - and since there are much better alternatives it is not quite worth doing so

For organisms that have significant splicing going on a lot many more of your reads should be multi-mapping in that case something is off. If your organism does not have splicing then there is no reason to use eXpress.

[bagi.m]

Thank you for your answer. I am indeed planing to try Salmon with this data, and see if I get to similar conclusion. I found some articles / posts claiming pseudoalignment tools are better, but both were based on human data. So I don't know if this claim holds for organisms with much higher intraspecies sequence variability (where tweaking alignment parameters is sometimes necessary).

[h.mon]

Pay attention because Salmon- like eXpress - doesn't want sorted input.