Question: mixed model using limma-voom
gravatar for Yanzhu Lin
2.0 years ago by
Yanzhu Lin120
Yanzhu Lin120 wrote:


I have questions with limma-voom RNA-Seq data analysis. My project has three factor: A (16 levels), B(three levels) and C(two levels), where A is random effect while B and C are fixed effects.

I would like to test the main effects of B and C, as well as A. To test the effect B, I used the following coding.


Group <- factor(paste(Design$A,Design$B,Design$C,Design$plateID,sep="."))
design <- model.matrix(~B+C)
y <- DGEList(counts=T,group=Group) ## T is the count data
y <- calcNormFactors(y,method="RLE")                  
v <- voom(y,design,plot=TRUE)

corfit <- duplicateCorrelation(v,design,block=Design$A)

fitRan <- lmFit(v,design,block=Design$A,correlation=corfit$consensus)
fitRan <- eBayes(fitRan)
topTable(fitRan,coef=c(2,3)) ## for testing the effect of factor B, which has three levels.

My questions:

1. Is the coding above correct for testing the main effect of B?

2. How to test the random effect A?






ADD COMMENTlink modified 2.0 years ago • written 2.0 years ago by Yanzhu Lin120

Hi Gordon,

Thank you very much for your comments and suggestion. I have filtered the very lowly expressed genes.

ADD REPLYlink written 2.0 years ago by Yanzhu Lin120
gravatar for Gordon Smyth
2.0 years ago by
Gordon Smyth35k
Walter and Eliza Hall Institute of Medical Research, Melbourne, Australia
Gordon Smyth35k wrote:

Yes, the code is correct for testing for B. (Although usually a filtering step is needed. Have you filtered out very lowly expressed genes?)

You can't test for factor A. It doesn't really make sense to test for a random effect.

ADD COMMENTlink modified 2.0 years ago • written 2.0 years ago by Gordon Smyth35k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 341 users visited in the last hour