hi Dvir,

Simon is pointing out that the RSEM output you pasted is calculated on the
transcript level, as the rows have a transcript_id column. DESeq operates
on gene-level counts, which involve collapsing all the possible transcripts
of a gene into one feature.  For more information on the details of
appropriate gene-level counting for DESeq, you can see this page of the
HTSeq documentation (
http://www-huber.embl.de/users/anders/HTSeq/doc/count.html ):

"the features are typically genes, where each gene is considered here as
the union of all its exons"

Suppose you want to do gene-level DE analysis, and you have 1 gene with 10
transcripts. The count of reads that hit this gene is somehow distributed
by the RSEM software over the 10 transcripts. Since we do not know exactly
how the reads are distributed by the RSEM software, we cannot answer how to
arrive at a sensible gene-level count.  I would recommend contacting the
RSEM developers to see if it is appropriate to add, for example the counts
in these files from the 10 transcripts to arrive at a gene-level count
which conforms to the definition from the HTSeq link above.

Mike


On Mon, Mar 25, 2013 at 12:41 PM, Dvir Netanely <dvir.tau@gmail.com> wrote:

> Hi Simon,
>
> Thanks again for your help.
>
> The data files I'm using are the GENE LEVEL and UNNORMALIZED and named
> "...expression_rsem_gene.txt". Exon level and normalized versions of the
> RSEM output also exist but I do not use them.
> These files contain a 'raw_count' column which I thought should be okay as
> input for DESeq and EdgeR.
>
> Do you still think it's calculated on the exon level ?
>
> Many thanks,
> Dvir
>
> -----Original Message-----
> From: Simon Anders [mailto:anders@embl.de]
> Sent: Thursday, March 21, 2013 3:20 PM
> To: bioconductor@r-project.org
> Subject: Re: [BioC] Applying DESeq on RSEM output
>
> Hi Dvir
>
> On 20/03/13 14:15, dvir.tau@gmail.com wrote:
> > I'm running DESeq and EdgeR on RNA-Seq data that was already processed
> > with RSEM (downloaded from TCGA web site).
> >
> > Since these methods require the raw read counts I'm using the
> > raw_count column of the RSEM output but I'm not sure this is the right
> > thing to do (is it the actual raw count required ?)
>
> The real issue is not that your counts are not integer, but that RSEM gives
> you counts per isoform rather than per gene. Now, if you have two very
> similar isoforms, RSEM will be unable to decide which isoform to assign a
> read to and just spread them proportionally over both. Hence, even if only
> one of the two isoforms is differentially expressed, you will incorrectly
> see differential expression for both isoforms.
>
> This is why the output of isoform quantification methods such as RSEM of
> MMSeq are not suitable as input for differential expression tests.
>
> At the very minimum, you need also the information about the uncertainty of
> the assignments of reads to isoforms. In fact, RSEM provides this
> information if you run it in its Bayesian mode, but this seems to be hardly
> ever done in practice.
>
> If you really need to perform differential expression analysis on a level
> finer than whole gene expression, you should either use a tool for
> differential exon usage testing, such as our DEXSeq package, or one that
> combines isoform abundance estimation and testing for differences in a
> unified framework, such as BitSeq. In both cases, you will need the SAM
> files.
>
> If you are fine with staying on the gene level for your analysis, you need
> to get counts per gene, not per isoform. I am not familiar enough with
> RSEM,
> though, to tell you whether adding up the counts from all the isoforms per
> gene would be a good idea.
>
>    Simon
>
>
>
> -----
> No virus found in this message.
> Checked by AVG - www.avg.com
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor@r-project.org
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives:
> http://news.gmane.org/gmane.science.biology.informatics.conductor
>

	[[alternative HTML version deleted]]

