Question: Standard error and effect size from Limma
1
3.5 years ago by
Vani20
United States
Vani20 wrote:

Is it possible to find the standard error and effect size of a dataset using limma's lmfit and toptable? Please advise.

modified 3.5 years ago by Gordon Smyth36k • written 3.5 years ago by Vani20
Answer: Standard error and effect size from Limma
5
3.5 years ago by
Gordon Smyth36k
Walter and Eliza Hall Institute of Medical Research, Melbourne, Australia
Gordon Smyth36k wrote:

The effect sizes are contained in fit$coefficients. The standard errors can be obtained from SE <- sqrt(fit$s2.post) * fit$stdev.unscaled This gives a matrix containing standard errors for every coefficient and every gene. ADD COMMENTlink modified 3.4 years ago • written 3.5 years ago by Gordon Smyth36k I am confused; shouldn't that be sqrt(fit$s2.post)? Similar to the line 255 of the toptable.R file in limma source package, the confidence intervals are calculated using

sqrt(eb$s2.post[top])*fit$stdev.unscaled[top,coef]*qt(alpha,df=eb$df.total[top]) ADD REPLYlink modified 3.4 years ago • written 3.4 years ago by January Weiner230 You are right. Now corrected. ADD REPLYlink modified 3.4 years ago • written 3.4 years ago by Gordon Smyth36k I'm sorry, just tried this code, but in my fit object I don't have 's2.post'. How do I get this? ADD REPLYlink written 17 months ago by s.w.vanderlaan20 If you've run eBayes() on the fit object, then you will haves s2.post. ADD REPLYlink written 17 months ago by Gordon Smyth36k However, I'm following this protocol here: https://molepi.github.io/DNAmArray_workflow/06_EWAS.html#correct_for_bias_and_inflation to run a EWAS. So, we start with limma(), but correct for inflation using bacon(). And thus, this is my code. designp <- model.matrix(metadata(Mvalues)$formula, data = colData(Mvalues))
datap <- assays(Mvalues)$data fitp <- limma::lmFit(datap, designp) Given that tstatp <-fitp$coef/fitp$stdev.unscaled/fitp$sigma had gotten me the T-statistic, I assumed fitp$stdev.unscaled/fitp$sigma was equal to standard error. And I assumed that fitp$coef would give me the effect sizes. If I don't run eBayes() where do I get the standard error from? Or is that effectively not possible? Thanks! ADD REPLYlink modified 17 months ago • written 17 months ago by s.w.vanderlaan20 If you want to ask questions about a non-standard workflow, you should start a new post with the appropriate tags. (In this case, it doesn't seem to be a Bioconductor package, so you might as well ask the authors directly.) The workflow in question is a bit bemusing as limma is run without EB shrinkage, which defeats the purpose - you might as well use lm.fit. Anyway, fitp$stdev.unscaled*fitp\$sigma is the standard error of the coefficient.

Ah, thank you for the answer, this was not immediately clear to me.

Answer: Standard error and effect size from Limma
4
3.5 years ago by
Denali
Steve Lianoglou12k wrote:

Yes.

Calling topTable(fit, ..., confint=TRUE) provides you with the effect size (logFC) and 95% confidence intervals (CI.L, CI.R), from which you can back calculate the standard error, if need be.

Cool thanks. Just a quick question: How would I calculate the standard error from the confidence intervals?

I'd say, (CI.R - CI.L)/3.92 (since the 95% confidence intervals should be plus/minus 1.96se)