[BioC] Limma, SAM and Fold Change Calculation

Gordon Smyth smyth at wehi.EDU.AU
Sat Jul 14 05:22:38 CEST 2007


BTW, here's an amusing example of the potential contradictions 
introduced by the SAM fold change definition:

 > y <- matrix(rnorm(100*20),100,20)
 > y[1,] <- rep(c(-2,2),each=10); y[1,1] <- y[1,1]-0.1
 > rownames(y) <- paste("Gene",1:100)
 > y[1,1] <- 6
 > out <- sam(y,cl)
 > summary(out,0.4)
SAM Analysis for the Two-Class Unpaired Case Assuming Unequal Variances

s0 = 0.8  (The 100 % quantile of the s values.)

  Number of permutations: 100

  MEAN number of falsely called variables is computed.

Delta: 0.4
cutlow: -Inf
cutup: 2
p0: 0.866
Identified Genes: 1
Falsely Called Genes: 0
FDR: 0


Identified Genes (using Delta = 0.4):

   Row d.value stdev  rawp q.value R.fold   Name
1   1       2   0.8 1e-04 0.00866  0.604 Gene 1

Note that Gene1 is significantly up-regulated (as indicated by a 
positive d.value=regularized t), but the fold change is less one, 
i.e., in the opposite direction!

Best wishes
Gordon

At 12:55 AM 14/07/2007, Holger Schwender wrote:
>Hi Peter, hi Gordon,
>
>the reason for the fold change computation in SAM is simply that 
>Tusher et al. also use these fold changes. In 2002, SAM was part of 
>my diploma thesis and I tried to reproduce the results of Tusher et 
>al. and the results that are provided by the Excel SAM version. 
>Since then many things in sam have changed, mainly to make sam 
>faster, more user-friendly and more flexible. But I didn't change 
>the computation of the fold change to actually avoid the question 
>why the fold changes differ between the siggenes SAM and the Excel SAM.
>
>Admittedly, "because Tusher et al. do so" is not a very good reason, 
>but I have stopped counting how often I was ask why this and that 
>differs between Excel SAM and siggenes SAM.
>
>Would it help you if I add another slot to the SAM object that 
>contains the numerator of the test statistics? Note that you can 
>compute these numerators by
>
> > sam.out at d * (sam.out at s0 + sam.out at s)
>
>if sam.out is the output of sam.
>
>Best,
>Holger



More information about the Bioconductor mailing list