[BioC] limma roast syntax for overall anova

[Gordon K Smyth]

Oh, yes, that's right.  I had forgotten that limma now recognizes the 
special intercept column name "(Intercept)" and removes it automatically. 
So you already had the correct anova test.

Best wishes
Gordon


On Thu, 3 Jul 2014, Juliet Hannah wrote:

> Hi Gordon,
>
> Thanks for your help (and software!).
>
> I confusingly showed only the first few rows of the design matrix. I was
> trying to convey that I left the
> intercept in.
>
> limma gave me the message:
>
>> mytt <- topTable(fit,number=Inf)
> Removing intercept from test coefficients
>
> which gives me the same result as the more clear version you gave me
>
> mytt <- topTable(fit,coef=2:4,number=Inf).
>
> I will proceed with geneSetTest.
>
> Thanks,
>
> Juliet
>
>
>
>
>
> On Thu, Jul 3, 2014 at 8:00 PM, Gordon K Smyth <smyth at wehi.edu.au> wrote:
>
>> Dear Juliet,
>>
>>
>> On Thu Jul 3 21:09:46 CEST 2014, Juliet Hannah wrote:
>>
>>  All,
>>>
>>> I am testing an overall anova (4 groups; is there any difference across
>>> groups) using the following design matrix in limma
>>>
>>>  (Intercept) B C D
>>> 1           1 0 0 0
>>> 2           1 0 0 0
>>> 3           1 0 0 0
>>> 4           1 0 0 0
>>> 5           1 0 0 0
>>> 6           1 0 0 0
>>>
>>
>> Mmm, the columns seem to be all zero except for the first.  I assume that
>> you are actually just giving the first 6 rows of a larger matrix.
>>
>>
>>  I have used ROAST (in limma) in the past with a specific contrast, such as
>>>
>>> roast(expression_matrix,index,design_matrix,contrast=10,nrot=3000).
>>>
>>> In the current case, I didn't need any contrasts. In limma, I had used
>>>
>>> design <- model.matrix(~key$stage)
>>> colnames(design)[2:4] = c("B","C","D")
>>> fit <- lmFit(e, design)
>>> fit <- eBayes(fit)
>>> topTable(fit,adjust.method="BH")
>>>
>>> to obtain the overall F-test
>>>
>>
>> To get the overall anova F-test you need:
>>
>>   topTable(fit[,-1])
>>
>> or alternatively
>>
>>   topTable(fit, coef=2:4)
>>
>> which gives the same result.  Note that the intercept needs to be removed
>> from the test, as it always is in classical anova, and that's what the -1
>> does.
>>
>>
>>  (I understand this test is not that interesting, but it is requested of
>>> me).
>>>
>>> Using this setup, how can I specify this overall test within ROAST.
>>>
>>
>> That's a good question.  We haven't implemented F-tests for ROAST.  If we
>> get enough requests, we will consider doing so, but we haven't so far.
>>
>> If the meantime, you can input the F-statistics to the geneSetTest()
>> function:
>>
>>   fit2 <- fit[,-1]
>>   geneSetTest(index, fit2$F)
>>
>> etc.
>>
>> Best wishes
>> Gordon
>>
>>
>>
>>  Thanks for your help.
>>>
>>> Juliet

______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}