[BioC] limma and p-values

[Giovanni Coppola]

Hi Gordon,
thanks for your reply.
I am analyzing 12 Agilent slides (6 dye swaps)
My previous steps were the following:

MA <- normalizeWithinArrays(RG, method="loess")
design <- c(1,-1,1,-1,1,-1,1,-1,1,-1,1,-1)
fit <- lmFit(MA,design)
fit<- eBayes(fit)

#is fit$p.value a matrix of p.values?
 >fit$p.value[1:10]
  [1]        NA        NA 0.9931832 0.1992031 0.4111630 0.4880706        NA
  [8] 0.8261918 0.2254227 0.2039834

my outputs were as follows:
 > p.corr.none.ord
  [1] 5.478123e-05 7.317678e-05 8.723562e-05 9.087221e-05 1.133770e-04
  [6] 3.923843e-04 4.478298e-04 6.033205e-04 7.023724e-04 8.562249e-04
 > toptable.none$P.Value
  [1] 5.478123e-05 7.317678e-05 8.723562e-05 9.087221e-05 1.133770e-04
  [6] 3.923843e-04 4.478298e-04 6.033205e-04 7.023724e-04 8.562249e-04

 > p.corr.fdr.ord
  [1] 0.001140071 0.001521493 0.001812122 0.001885914 0.002350789 0.008128272
  [7] 0.009268255 0.012474751 0.014509433 0.017671378
 > toptable.fdr$P.Value
  [1] 0.4833943 0.4833943 0.4833943 0.4833943 0.4833943 0.9990538 0.9990538
  [8] 0.9990538 0.9990538 0.9990538

thank you for your time (and patience!)
Giovanni






At 04:04 AM 1/6/2005, you wrote:
> > Date: Wed, 05 Jan 2005 18:46:10 -0800
> > From: Giovanni Coppola <gcoppola at ucla.edu>
> > Subject: [BioC] limma and p-values
> > To: bioconductor at stat.math.ethz.ch
> > Message-ID: <6.1.2.0.2.20050105175459.04ccf508 at mail.ucla.edu>
> > Content-Type: text/plain; charset="us-ascii"; format=flowed
> >
> > Hi everybody,
> > my question is again related to p-values...
> > Is the topTable function using only 'p.adjust' to adjust the p.values?
>
>Yes it uses p.adjust(), and uses only p.adjust().
>
>The mistake you have made below is to apply p.adjust() to the entire 
>matrix of p-values
>fit$p-value whereas topTable() applies p.adjust() to each column separately.
>
> > #R2.0.1, limma 1.8.13
> >
> > fit<-eBayes(fit)
> > #After applying eBayes, fit$p.value contains the p-values associated with
> > the modified t-statistics
> >
> > p.corr.none <- p.adjust (fit$p.value, "none")
> > p.corr.none.ord<-sort(p.corr.none,decreasing=FALSE,na.last=NA)
> > p.corr.none.ord<-p.corr.none.ord[1:10]
> > toptable.none <-topTable(fit,number=10,adjust="none",sort.by="P")
> >
> > # at this point, 'p.corr.none' = 'toptable.none$P.Value' , whereas...
> >
> > p.corr.fdr <- p.adjust (fit$p.value, "fdr")
> > p.corr.fdr.ord<-sort(p.corr.fdr,decreasing=FALSE,na.last=NA)
> > p.corr.fdr.ord<-p.corr.fdr.ord[1:10]
> > toptable.fdr <-topTable(fit,number=10,adjust="fdr",sort.by="P")
> >
> > #... 'p.corr.fdr' and 'toptable.fdr$P.Value' are not the same, because
> > 'toptable.fdr$P.Value' is much worse...why?
> >
> > Thanks #and sorry for bringing up p-values again :-)
> > Giovanni
>
>It is generally helpful to give an example which is complete and 
>reproducible and which includes
>the output to prove your point.  In this case I can guess that you've 
>fitted a linear model with
>more than one coefficient.  Had there been only one coefficient, the 
>adjusted p-values would have
>turned out equal:
>
> > M <- matrix(rnorm(1000*10),1000,10)
> > X <- matrix(rnorm(10*1),10,1)
> > X <- X[,1]
> > fit <- lmFit(M,X)
> > fit <- eBayes(fit)
> > p.corr.fdr <- p.adjust (fit$p.value, "fdr")
> > p.corr.fdr.ord<-sort(p.corr.fdr,decreasing=FALSE,na.last=NA)
> > p.corr.fdr.ord<-p.corr.fdr.ord[1:10]
> > toptable.fdr <-topTable(fit,number=10,adjust="fdr",sort.by="P")
> > toptable.fdr
>              M         t   P.Value         B
>287 -1.1922028 -3.908555 0.1079313 -1.212784
>495 -0.9497453 -3.107945 0.7510069 -2.561172
>338 -0.9432092 -3.070181 0.7510069 -2.617739
>494  0.9171600  2.984433 0.7510069 -2.743777
>38   0.8098800  2.660003 0.9661470 -3.190147
>582  0.8134650  2.641185 0.9661470 -3.214544
>194 -0.7925727 -2.601763 0.9661470 -3.265120
>274 -0.7757093 -2.511355 0.9661470 -3.378368
>429 -0.7725337 -2.509163 0.9661470 -3.381067
>597 -0.7519274 -2.451572 0.9661470 -3.451150
> > p.corr.fdr.ord
>  [1] 0.1079313 0.7510069 0.7510069 0.7510069 0.9661470 0.9661470 
> 0.9661470 0.9661470 0.9661470
>[10] 0.9661470
>
>Gordon