This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How do I run DE for TPM values (not CPM)?

Can someone explain to me how to run DE using TPM instead of CPM, please? All the DE guides I'm seeing only use CPM values, but I need to work with TPM values. Providing either a guide of code to use or a manual/tutorial providing one would be most appreciated.

Also, note that I'm working with a .gz file that I downloaded onto RStudio with data.table::fread, which seems like it would affect how I would need to code the DE process in terms of what kinds of objects I'm using.

rstudio differential-expression tpm
• 3,540 views • link
updated by Gordon Smyth Gordon Smyth • written by John • John

2 answers

TPM and CPM are highly similar scales, so if you're following (say) a tutorial for DE with CPMs using LIMMA (which I would recommend if you're stuck with this pre-transformed data), then you can use TPM values wherever you see the CPMs and you should be good to go.

• 0 views • link
written by LChart LChart

Okay, so even if the tutorial gives a code like

logCPM <- cpm(dge, log=TRUE, prior.count=3)

It should still be fine to just put in my tpm file where CPM should go? (Note here that they use a cpm function, but there isn't a tpm function, to my knowledge.)

And as for the second part of what I asked, here is what I have so far:

dge <- DGEList(matchedgeneTPM)
dge <- calcNormFactors(dge)
dge$samples
fnames <- colnames(females)
mnames <- colnames(males)
group <- interaction(fnames, mnames)
fnames <- colnames(females)
mnames <- colnames(males)
group <- interaction(fnames, mnames)
plotMDS(dge, col = as.numeric(group))
mm <- model.matrix(~0 + group)
fit <- lmFit(dge, mm)

"matchedgeneTPM" is the table of TPM's I have (I already converted to log2(TPM+1) beforehand) and "females" and "males" are the list of female and male sample ID's I also constructed beforehand. Everything works fine until the last line, which gives the error "Error in getEAWP(object) : data object isn't of a recognized data class". Would you be able to explain how to fix this error, please?

• 0 views • link
written by John • John

The cpm() function takes in a matrix of counts. So does calcNormFactors.

Just take your TPM, log it, and pass it through limma.

• 0 views • link
written by LChart LChart

That's what I'm trying to do but again, I'm getting an error on the last line and I'm asking how to fix it.

• 0 views • link
written by John • John

you're mixing two packages, edgeR ("dgelist") and limma ("lmFit"). Just put your log(1+TPM) matrix as the first argument to lmFit.

• 0 views • link
written by LChart LChart

Thanks, but it's giving me the error that expression object should be numeric and that there are 2 non-numeric columns (because the first two columns are the labels of the genes and gene ID's I'm working with). Here is my code:

fit <- lmFit(matchedgeneTPM)[,-c(1,2)]

How should I reformat the code? And apologies, I know this should be simple but I'm new to this kind of stuff.

Also, do I need the stuff with group and mm?

• 0 views • link
written by John • John

See "Differential expression analysis starting from TPM data" https://support.bioconductor.org/p/98820/

The short answer is that you can analyse log(TPM+1) data passably well using limma with arrayWeights and trend, but you will nevertheless pay a heavy price in terms of statistical power compared to an analysis with actual read counts and library sizes. The use of arrayWeights() in limma is a way to try to estimate the library sizes that the TPMs were originally computed from.

• 0 views • link
written by Gordon Smyth Gordon Smyth

Log in to answer this question.

More posts like this