RPKM and multiple reads, tophat and cufflinks

[mgogol]

My understanding of the RPKM calculation in Tophat is that it includes multi-reads that match < 40 times in the genome (by default).

It seems like cufflinks does something more complex involving some kind of allocation of multiple reads.

http://cufflinks.cbcb.umd.edu/manual.html

Cufflinks models the sequencing process by asking what the probability is of observing each read, given a set of transcripts and a set of abundances. The program then multiplies these probabilities to compute the overall likelihood that one would observe the reads in the experiment, given the proposed abundances on the transcripts. Because Cufflinks' statistical model is linear, the likelihood function has a unique maximum value, and Cufflinks finds it with a numerical optimization algorithm.

I know we're supposed to start using cufflinks' RPKM now, but I'd like to understand tophat's as well. Does anyone know if my description is correct?

[Cole Trapnell]

Correct - both TopHat and Cufflinks include multireads in their RPKM calculations

[sandmann]

I would like to use Cufflinks to assemble SOLiD reads, so unfortunately I can't use bowtie. Are there any recommendations about the maximum number of mappings any read should not exceed to be included in Cufflink's analysis ? In other words, should I discard reads that map e.g. to >40 locations in the genome before running Cufflinks ?

Any advice is greatly appreciated.

[Cole Trapnell]

You should try leaving them in, but make sure their mapping quality is set to 0. Cufflinks will use them to assemble transcripts, but they will contribute little to nothing to RPKM's of transcripts that contain them. Leaving them out could create voids in read coverage where there are transcribed repeats (as happens not uncommonly in UTRs, etc), which could break up the gene models.

[HTS]

Please allow me to ask two related questions:

1. If I do use TopHat and Bowtie to do the alignment, will the mapping qualities of multireads be automatically adjusted? From the Bowtie manual it seems that mapping quality is always set to 255 at the moment.

2. Is the default value of 40 a good one for mammalian genomes? I have been using 10 in my own calculation so far and I am a bit curious about this default value.

Thanks a lot!

-- Leo

[Cole Trapnell]

Please allow me to ask two related questions:

1. If I do use TopHat and Bowtie to do the alignment, will the mapping qualities of multireads be automatically adjusted? From the Bowtie manual it seems that mapping quality is always set to 255 at the moment.

Yes, TopHat calculates its own mapping quality. After ranking the alignments for each read, TopHat outputs all the alignments in the top rank. If there are more than one, then the probability that the read is incorrect (i.e. the mapping quality) is (1 - 1/p), assuming one of the alignments in the top rank *is* the correct one.

2. Is the default value of 40 a good one for mammalian genomes? I have been using 10 in my own calculation so far and I am a bit curious about this default value.

It should be fine - Reducing it might cut down on your running time a bit, but not drastically in my opinion.

[carmeyeii]

Hi,

VEEERY basic question: When calculating FPKM, what is the "correct" M value to use: Total unique mapped reads (i.e., the awk-ed file of accepted hits), or ALL the hits, i.e, all the alignments that TopHat output? What does cufflinks take as the M in FPKM?

Thanks for your help.

Carmen