F-tests for factorial effects - limma
1
0
Entering edit mode
Naomi Altman ★ 6.0k
@naomi-altman-380
Last seen 2.9 years ago
United States
I am analyzing a 2-factor factorial Affy experiment, with 3 d.f. for each factor. I would like to get the F-tests for the main effects and interactions using limma. I have computed all the contrasts, and got the t-tests (both unadjusted and eBayes). I do know how to combine these into F-tests "by hand" but I could not figure out if there was a simple way to do this using limma. I had a look at FStat (classifyTestsF). There seems to be a problem, in that the matrix tstat is not premultiplied by the contrast matrix when the F-statistics are computed. So, if the contrasts are not full-rank, an error is generated (instead of the F-statistics) because nrow(Q) != ncol(tstat).. (See the lines below). if (fstat.only) { fstat <- drop((tstat%*% Q)^2 %*% array(1, c(r, 1))) attr(fstat, "df1") <- r attr(fstat, "df2") <- df return(fstat) } I figured that before I fiddled with the code, I would check to make sure that I didn't miss an existing routine to do this. Thanks in advance. Naomi S. Altman 814-865-3791 (voice) Associate Professor Bioinformatics Consulting Center Dept. of Statistics 814-863-7114 (fax) Penn State University 814-865-1348 (Statistics) University Park, PA 16802-2111
affy limma affy limma • 826 views
ADD COMMENT
0
Entering edit mode
@gordon-smyth
Last seen 1 hour ago
WEHI, Melbourne, Australia
> Date: Tue, 21 Dec 2004 17:21:19 -0500 > From: Naomi Altman <naomi@stat.psu.edu> > Subject: [BioC] F-tests for factorial effects - limma > To: bioconductor@stat.math.ethz.ch > > I am analyzing a 2-factor factorial Affy experiment, with 3 d.f. for each > factor. > > I would like to get the F-tests for the main effects and interactions using > limma. > > I have computed all the contrasts, and got the t-tests (both unadjusted and > eBayes). I do know how to combine these into F-tests "by hand" but I could > not figure out if there was a simple way to do this using limma. limma doesn't have any easy way to deal with main effects and interactions, at least not with main effects, interactions are actually simpler. I haven't implemented this because I've never been able to figure out what one would do with these things in a microarray context. To compute F-tests for main effects and interaction, the easiest way would probably be to compute the SS for main effects and interactions by non-limma means, then use shrinkVar() to adjust the residual mean squares, i.e., the F-statistic denominators. If you only want F-tests for interactions, the following code would work: X <- model.matrix(~a*b) fit <- lmFit(eset, X) p <- ncol(X) cont.ia <- diag(p)[,attr(X,"assign")==3] fit.ia <- eBayes(contrasts.fit(fit, cont.ia)) Now fit.ia contains the F-statistic and p-values for the interaction in fit.ia$F and fit.ia$F.p.value. > I had a look at FStat (classifyTestsF). There seems to be a problem, in > that the matrix tstat is not premultiplied by the contrast matrix when the > F-statistics are computed. So, if the contrasts are not full-rank, an > error is generated (instead of the F-statistics) because nrow(Q) != > ncol(tstat).. (See the lines below). No, the code is correct. FStat is quite happy with non full rank contrasts but the contrast matrix must be applied using contrasts.fit() before entering FStat(). You should not expect to see a contrast matrix inside the classifyTestsF() code. Gordon > if (fstat.only) { > fstat <- drop((tstat%*% Q)^2 %*% array(1, c(r, 1))) > attr(fstat, "df1") <- r > attr(fstat, "df2") <- df > return(fstat) > } > > I figured that before I fiddled with the code, I would check to make sure > that I didn't miss an existing routine to do this. > > Thanks in advance. > > Naomi S. Altman 814-865-3791 (voice) > Associate Professor > Bioinformatics Consulting Center > Dept. of Statistics 814-863-7114 (fax) > Penn State University 814-865-1348 (Statistics) > University Park, PA 16802-2111
ADD COMMENT

Login before adding your answer.

Traffic: 777 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6