makeOrderings {baySeq} | R Documentation |
Given a model structure as defined in the ‘@groups’ slot of a
countData
object containing more than one group, it is
often possible to define an ordering on the groups for a given genomic
event. To take a simple example, if the average expression of a gene
is higher in sample set A then in sample set B, then we might impose
an ordering A>B.
makeOrderings(cD, orderingFunction)
cD |
A |
orderingFunction |
A function defining the orderings. If not given, will be taken from the ‘@densityFunction’ slot of ‘cD’. See Details, and examples below. |
The orderingFunction takes 'dat' and 'observables' as arguments. 'dat' is equivalent to the ‘@data’ slot of the ‘cD’ object, and 'observables' the combined data in the ‘@sampleObservables’, ‘@rowObservables’ and ‘@cellObservables’ slots.
A countData
with populated ‘@orderings’ slot.
Thomas J. Hardcastle
Hardcastle T.J., and Kelly, K. baySeq: Empirical Bayesian Methods For Identifying Differential Expression In Sequence Count Data. BMC Bioinformatics (2010)
# load test data data(simData) # Create a countData object from test data replicates <- c("simA", "simA", "simA", "simA", "simA", "simB", "simB", "simB", "simB", "simB") groups <- list(NDE = c(1,1,1,1,1,1,1,1,1,1), DE = c(1,1,1,1,1,2,2,2,2,2)) CD <- new("countData", data = simData, replicates = replicates, groups = groups) libsizes(CD) <- getLibsizes(CD) # order on expression normalised for library size (scaling factor) and gene length CD <- makeOrderings(CD, orderingFunction = function(dat, observables) dat / observables$libsizes) # orderings calculated for DE group head(CD@orderings) # load test (paired) data data(pairData) # create a countData object from paired data pairCD <- new("countData", data = list(pairData[,1:4], pairData[,5:8]), replicates = c(1,1,2,2), groups = list(NDE = c(1,1,1,1), DE = c(1,1,2,2)), densityFunction = bbDensity) libsizes(pairCD) <- getLibsizes(pairCD) # order on (log-)ratio of pairs, with fudge-factor on zeros. pairCD <- makeOrderings(pairCD, orderingFunction = function(dat, observables) { data <- dat / observables$libsizes adjmin <- min(data[data > 0]) / 10 log(data[,,1] + adjmin) - log(data[,,2] + adjmin) }) # orderings calculated for DE group head(pairCD@orderings)