Search
Question: mixed model using limma-voom
0
21 months ago by
Yanzhu Lin110
Yanzhu Lin110 wrote:

Hi ALL,

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) corfit$consensus

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?

Thanks,

Yanzhu

modified 21 months ago • written 21 months ago by Yanzhu Lin110

Hi Gordon,

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

0
21 months ago by
Gordon Smyth34k
Walter and Eliza Hall Institute of Medical Research, Melbourne, Australia
Gordon Smyth34k 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.