[BioC] Bioconductor

Claire Wilson ClaireWilson at picr.man.ac.uk
Wed May 21 16:18:43 MEST 2003


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}}



More information about the Bioconductor mailing list