[BioC] Re: multiple comparisons in limma [was: Now I get it(almost)]

[Gordon Smyth]

At 05:48 AM 7/08/2004, Richard Friedman wrote:
>Gordon,
>
>         I now believe that I understand your answer. In order to do 
> adjust for both  multiple comparisons and
>multiple tests I use classifyTestsF() with method="fdr". If I understand 
>the documentation correctly, the fdr part
>refers to the multiple test correction across genes on top of the multiple 
>comparison adjustment across levels performed
>  if no method were to be specified.
>
>         Do I have it straight?

No, not quite. classifyTestsF() doesn't actually have a 'method' option, so 
you can't use it with method="fdr".

Please understand that making multiple adjustments simutaneously across 
genes and across contrasts (comparisons between levels for example) is 
still an unsolved problem, so it may be too much to expect neat solutions 
from limma or concise explanations from me.

One way that you could proceed is to compute an F-test p-value for each 
gene, testing for any differences betwen the levels of your factor, and use 
p.adjust() to apply "fdr" adjustment to those p-values across genes. Then 
for those genes which are chosen as showing some differences, you can 
examine the individual contrast p-values more closely to decide which 
levels are different.

Notice that the F-test p-value is computed automatically by eBayes() and is 
stored as fit$F.p.value in your MArrayLM linear model data object.

Here is some example code which applies "fdr" adjustment in a hierarchical 
fashion, first at the overall gene level, and then at the contrast level 
after deciding on the gene cutoff:

fit <- eBayes(fit)
selectedgenes <- p.adjust(fit$F.p.value, method="fdr") < 0.05
pmax <- pmin(fit$F.p.value[!selectedgenes], na.rm=TRUE)
results <- classifyTestsP(fit, p.value=pmax, method="fdr")

Gordon

>Thanks and best wishes,
>Rich
>On Aug 6, 2004, at 9:57 AM, Gordon Smyth wrote:
>
>>At 11:36 PM 6/08/2004, Richard Friedman wrote:
>>>Gordon,
>>>
>>>         Thank you for answering my questions. The last equation in your 
>>> paper makes intuitive sense to me.
>>>
>>>         I'm wondering if you can take the time to answer two  more 
>>> questions:
>>>
>>>1. Say I have the following case:
>>>
>>>Level A (3 replicates)
>>>Level B (2 replicates)
>>>Level C(1 replicate)
>>>Level D(1 replicate)
>>>
>>>Can I legitimately calculate a P value for the contrast Level A to level 
>>>C in the linear model even though I have only one replicate on Level C. 
>>>I am not talking about just Limma here. I am talking about the linear 
>>>model in general.
>>
>>Given assumption of common variance across levels, yes.
>>
>>>  Also,
>>>I realize that one replicate is poor experimental design. This is what I 
>>>was given to analyze.
>>>
>>>2. If I wished to apply a multiple test correction to the pvalues from 
>>>non-orthogonal contrasts, would the following procedure be legitimate::
>>>
>>>         1. Generate a pvalue for each contrast in the set of 
>>> nonothogonal contrasts for each gene  using classifyTestsF().
>>>         2. Correct the pvalues using a multiple test correction such as 
>>> FDR.
>>
>>Nothing special about this design. All usual things, e.g. in limma, apply.
>>
>>Gordon
>>
>>>I realize that no multiple-test correction is entirely satisfactory, I 
>>>just want to get an approximate estimate of the p-values for each 
>>>contrast as a guide to further experimentation and literature searching.
>>>
>>>Best wishes,
>>>Rich