Hi Richard, For t-tests, this is what I do: do.ttest <- function(y, x1=1, x2=2) { # split y according to categories in pData # For example you may have cell.line as a cateogry in pData # because the function is called with esApply, it knows the values of pData # therefore knows that cell.line is one of the categories # defaults for x1 and x2 are set as the first two categories # y[[1]]= cell line a, y[[2]] = cell line b, y[[3]] = cell line c, y[[4]] = cell line d exprs.values <- split(y, cell.line) # Do a t-test between the expression values of cell line a and cell line b # Return the p-values t.test(exprs.values[[x1]], exprs.values[[x2]])$p.value } # Function is called as follows # Changing x1 and x2 means that you can change which two cell lines are analysed # 1 means analyse eset row-wise > t.res <- esApply(eset,1, do.ttest, x1=1, x2=2) # t.res holds a list of probe set ids and their p-values For anova: calc.anova.pval <- function(x) { # [[1]] takes us to the summary information # [4:5] are columns in the summary data that contain the F-value and the P-value # [2] selects the P-value and [1] says the 1st row of P-value # Not too sure on the details of the numbers but it works for me! # Again I call it with esApply, so it knows what the pData categories are summary(aov(x ~cell.line))[[1]][4:5][[2]][1] } # Call the function and store the anova p-values as a list with the probeset ids > anova.pvalues <- esApply(eset, 1, calc.anova.pval) # Another way to do ANOVA, but you don't get the p-values back # Returns a logical vector which indicates whether the p-value from the ANOVA # is less than a specified value # Declare the filter # Test the expression levels in eset according to cell.line # Return true for all those with a p-value <= 0.01 > anova.filter <- Anova(eset$cell.line, p=0.01) # Apply the filter, but you don't get the raw p-values, therefore can't filter according to p-value > signif.anova.01 <- genefilter(exprs(eset), filterfun(anova.filter)) None of the above have included multiple testing, which is necessary if you are analysing large numbers of probesets - See the Bioconductor multtest package. Hope this helps/makes sense/is right! claire -- Claire Wilson Bioinformatics group Paterson Institute for Cancer Research Christies Hospital NHS Trust Wilmslow Road, Withington Manchester M20 4BX tel: +44 (0)161 446 8218 url: http://bioinf.picr.man.ac.uk/ -------------------------------------------------------- This email is confidential and intended solely for the use of th... {{dropped}}

An innocent remark about the sequence of indexing you are "not sure about but that works for you". There have been a "very similar" piece of code circulating around. It works for others. (if time to kill, try to google on 'summary aov "[1]][4:5][[2]][1]"'). Hopin' it helps, Laurent On Wed, May 21, 2003 at 03:18:43PM +0100, Claire Wilson wrote: > Hi Richard, > > For t-tests, this is what I do: > > do.ttest <- function(y, x1=1, x2=2) { > # split y according to categories in pData > # For example you may have cell.line as a cateogry in pData > # because the function is called with esApply, it knows the values of pData > # therefore knows that cell.line is one of the categories > # defaults for x1 and x2 are set as the first two categories > # y[[1]]= cell line a, y[[2]] = cell line b, y[[3]] = cell line c, y[[4]] = cell line d > exprs.values <- split(y, cell.line) > # Do a t-test between the expression values of cell line a and cell line b > # Return the p-values > t.test(exprs.values[[x1]], exprs.values[[x2]])$p.value > } > > # Function is called as follows > # Changing x1 and x2 means that you can change which two cell lines are analysed > # 1 means analyse eset row-wise > > t.res <- esApply(eset,1, do.ttest, x1=1, x2=2) > # t.res holds a list of probe set ids and their p-values > > For anova: > calc.anova.pval <- function(x) { > # [[1]] takes us to the summary information > # [4:5] are columns in the summary data that contain the F-value and the P-value > # [2] selects the P-value and [1] says the 1st row of P-value > # Not too sure on the details of the numbers but it works for me! > # Again I call it with esApply, so it knows what the pData categories are > summary(aov(x ~cell.line))[[1]][4:5][[2]][1] > } > > # Call the function and store the anova p-values as a list with the probeset ids > > anova.pvalues <- esApply(eset, 1, calc.anova.pval) > > # Another way to do ANOVA, but you don't get the p-values back > # Returns a logical vector which indicates whether the p-value from the ANOVA > # is less than a specified value > # Declare the filter > # Test the expression levels in eset according to cell.line > # Return true for all those with a p-value <= 0.01 > > anova.filter <- Anova(eset$cell.line, p=0.01) > > # Apply the filter, but you don't get the raw p-values, therefore can't filter according to p-value > > signif.anova.01 <- genefilter(exprs(eset), filterfun(anova.filter)) > > None of the above have included multiple testing, which is necessary if you are analysing large numbers of probesets - See the Bioconductor multtest package. > > Hope this helps/makes sense/is right! > > claire > -- > Claire Wilson > Bioinformatics group > Paterson Institute for Cancer Research > Christies Hospital NHS Trust > Wilmslow Road, > Withington > Manchester > M20 4BX > tel: +44 (0)161 446 8218 > url: http://bioinf.picr.man.ac.uk/ > > -------------------------------------------------------- > > > This email is confidential and intended solely for the use of th... {{dropped}} > > _______________________________________________ > Bioconductor mailing list > Bioconductor@stat.math.ethz.ch > https://www.stat.math.ethz.ch/mailman/listinfo/bioconductor -- -------------------------------------------------------------- currently at the National Yang-Ming University in Taipei, Taiwan -------------------------------------------------------------- Laurent Gautier CBS, Building 208, DTU PhD. Student DK-2800 Lyngby,Denmark tel: +45 45 25 24 89 http://www.cbs.dtu.dk/laurent