R/AllGenerics.R
In StoreyLab/edge: Extraction of Differential Gene Expression
#' Performs F-test (likelihood ratio test using Normal likelihood)
#'
#' \code{lrt} performs a generalized likelihood ratio test using the full and
#' null models.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}. Optional.
#' @param nullDistn \code{character}: either "normal" or "bootstrap", If
#' "normal" then the p-values are calculated using the F distribution. If
#' "bootstrap" then a bootstrap algorithm is implemented to simulate
#' statistics from the null distribution. In the "bootstrap" case, empirical
#' p-values are calculated using the observed and null statistics (see
#' \code{\link{empPvals}}). Default is "normal".
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#' @param bs.its \code{integer}: number of null statistics generated (only
#' applicable for "bootstrap" method). Default is 100.
#' @param seed \code{integer}: set the seed value. Default is NULL.
#' @param verbose \code{boolean}: print iterations for bootstrap method.
#' Default is TRUE.
#' @param mod.F \code{boolean}: Moderated F-test, recommended for experiments
#' with a small sample size. Default is FALSE.
#' @param ... Additional arguments for \code{\link{apply_qvalue}} and
#' \code{\link{empPvals}} function.
#'
#' @details \code{lrt} fits the full and null models to each gene using the
#' function \code{\link{fit_models}} and then performs a likelihood ratio test.
#' The user has the option to calculate p-values a Normal distribution
#' assumption or through a bootstrap algorithm. If \code{nullDistn} is
#' "bootstrap" then empirical p-values will be determined from the
#' \code{\link{qvalue}} package (see \code{\link{empPvals}}).
#'
#' @author John Storey, Andrew Bass
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # lrt method
#' de_lrt <- lrt(de_obj, nullDistn = "normal")
#'
#' # to generate p-values from bootstrap
#' de_lrt <- lrt(de_obj, nullDistn = "bootstrap", bs.its = 30)
#'
#' # input a deFit object directly
#' de_fit <- fit_models(de_obj, stat.type = "lrt")
#' de_lrt <- lrt(de_obj, de.fit = de_fit)
#'
#' # summarize object
#' summary(de_lrt)
#'
#' @references
#' Storey JD, Xiao W, Leek JT, Tompkins RG, and Davis RW. (2005) Significance
#' analysis of time course microarray experiments. Proceedings of the National
#' Academy of Sciences, 102: 12837-12842.
#'
#' \url{http://en.wikipedia.org/wiki/Likelihood-ratio_test}
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{build_models}},
#' \code{\link{odp}}
#'
#' @export
setGeneric("lrt", function(object, de.fit,
nullDistn = c("normal","bootstrap"), weights = NULL,
bs.its = 100, seed = NULL, verbose = TRUE,
mod.F = FALSE, ...)
standardGeneric("lrt"))
#' The optimal discovery procedure
#'
#' \code{odp} performs the optimal discovery procedure, which is a framework for
#' optimally performing many hypothesis tests in a high-dimensional study. When
#' testing whether a feature is significant, the optimal discovery procedure
#' uses information across all features when testing for significance.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}. Optional.
#' @param odp.parms \code{list}: parameters for each cluster. See
#' \code{\link{kl_clust}}.
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#' @param bs.its \code{numeric}: number of null bootstrap iterations. Default
#' is 100.
#' @param n.mods \code{integer}: number of clusters used in
#' \code{\link{kl_clust}}. Default is 50.
#' @param seed \code{integer}: set the seed value. Default is NULL.
#' @param verbose \code{boolean}: print iterations for bootstrap method.
#' Default is TRUE.
#' @param ... Additional arguments for \code{\link{qvalue}} and
#' \code{\link{empPvals}}.
#'
#'
#' @details
#' The full ODP estimator computationally grows quadratically with respect to
#' the number of genes. This becomes computationally taxing at a certain point.
#' Therefore, an alternative method called mODP is used which has been shown to
#' provide results that are very similar. mODP utilizes a clustering algorithm
#' where genes are assigned to a cluster based on the Kullback-Leiber distance.
#' Each gene is assigned an module-average parameter to calculate the ODP score
#' and it reduces the computations time to approximately linear (see Woo, Leek
#' and Storey 2010). If the number of clusters is equal to the number of genes
#' then the original ODP is implemented. Depending on the number of hypothesis
#' tests, this can take some time.
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov,
#' null.model = null_model, full.model = full_model)
#'
#' # odp method
#' de_odp <- odp(de_obj, bs.its = 30)
#'
#' # input a deFit object or ODP parameters ... not necessary
#' de_fit <- fit_models(de_obj, stat.type = "odp")
#' de_clust <- kl_clust(de_obj, n.mods = 10)
#' de_odp <- odp(de_obj, de.fit = de_fit, odp.parms = de_clust,
#' bs.its = 30)
#'
#' # summarize object
#' summary(de_odp)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Woo S, Leek JT, Storey JD (2010) A computationally efficient modular
#' optimal discovery procedure. Bioinformatics, 27(4): 509-515.
#'
#' @author John Storey, Jeffrey Leek, Andrew Bass
#'
#' @seealso \code{\link{kl_clust}}, \code{\link{build_models}} and
#' \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("odp", function(object, de.fit, odp.parms = NULL, weights = NULL, bs.its = 100,
n.mods = 50, seed = NULL, verbose = TRUE, ...)
standardGeneric("odp"))
#' Modular optimal discovery procedure (mODP)
#'
#' \code{kl_clust} is an implementation of mODP that assigns genes to modules
#' based on of the Kullback-Leibler distance.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}.
#' @param n.mods \code{integer}: number of modules (i.e., clusters).
#'
#' @details mODP utilizes a k-means clustering algorithm where genes are
#' assigned to a cluster based on the Kullback-Leiber distance. Each gene is
#' assigned an module-average parameter to calculate the ODP score (See Woo,
#' Leek and Storey 2010 for more details). The mODP and full ODP produce nearly
#' exact results but mODP has the advantage of being computationally
#' faster.
#'
#' @note The results are generally insensitive to the number of modules after a
#' certain threshold of about n.mods>=50 in our experience. It is recommended
#' that users experiment with the number of modules. If the number of modules
#' is equal to the number of genes then the original ODP is implemented.
#'
#' @return
#' A list with the following slots:
#' \itemize{
#' \item {mu.full: cluster averaged fitted values from full model.}
#' \item {mu.null: cluster averaged fitted values from null model.}
#' \item {sig.full: cluster standard deviations from full model.}
#' \item {sig.null: cluster standard deviations from null model.}
#' \item {n.per.mod: total members in each cluster.}
#' \item {clustMembers: cluster membership for each gene.}
#' }
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # mODP method
#' de_clust <- kl_clust(de_obj)
#'
#' # change the number of clusters
#' de_clust <- kl_clust(de_obj, n.mods = 10)
#'
#' # input a deFit object
#' de_fit <- fit_models(de_obj, stat.type = "odp")
#' de_clust <- kl_clust(de_obj, de.fit = de_fit)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Woo S, Leek JT, Storey JD (2010) A computationally efficient modular optimal
#' discovery procedure. Bioinformatics, 27(4): 509-515.
#'
#' @author John Storey, Jeffrey Leek
#'
#' @seealso \code{\link{odp}}, \code{\link{fit_models}}
#'
#' @exportMethod kl_clust
setGeneric("kl_clust", function(object, de.fit = NULL, n.mods = 50)
standardGeneric("kl_clust"))
#' Linear regression of the null and full models
#'
#' \code{fit_models} fits a model matrix to each gene by using the least
#' squares method. Model fits can be either statistic type "odp" (optimal
#' discovery procedure) or "lrt" (likelihood ratio test).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param stat.type \code{character}: type of statistic to be used. Either
#' "lrt" or "odp". Default is "lrt".
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#'
#' @details
#' If "odp" method is implemented then the null model is removed from the full
#' model (see Storey 2007). Otherwise, the statistic type has no affect on the
#' model fit.
#'
#' @note \code{fit_models} does not have to be called by the user to use
#' \code{\link{odp}}, \code{\link{lrt}} or \code{\link{kl_clust}} as it is an
#' optional input and is implemented in the methods. The
#' \code{\linkS4class{deFit}} object can be created by the user if a different
#' statistical implementation is required.
#'
#' @return \code{\linkS4class{deFit}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # retrieve statistics from linear regression for each gene
#' fit_lrt <- fit_models(de_obj, stat.type = "lrt") # lrt method
#' fit_odp <- fit_models(de_obj, stat.type = "odp") # odp method
#'
#' # summarize object
#' summary(fit_odp)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Storey JD, Xiao W, Leek JT, Tompkins RG, and Davis RW. (2005) Significance
#' analysis of time course microarray experiments. Proceedings of the National
#' Academy of Sciences, 102: 12837-12842.
#'
#' @seealso \code{\linkS4class{deFit}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey
#' @exportMethod fit_models
setGeneric("fit_models",
function(object, stat.type = c("lrt", "odp"), weights = NULL) {
standardGeneric("fit_models")
})
#' Create a deSet object from an ExpressionSet
#'
#' Creates a \code{\linkS4class{deSet}} object that extends the
#' \code{\link{ExpressionSet}} object.
#'
#' @param object \code{S4 object}: \code{\link{ExpressionSet}}
#' @param full.model \code{formula}: full model containing the both the
#' adjustment and the biological variables for the experiment.
#' @param null.model \code{formula}: null model containing the adjustment
#' variables for the experiment.
#' @param individual \code{factor}: information on repeated samples in
#' experiment.
#'
#' @note It is essential that the null and full models have the same variables
#' as the ExpressionSet phenoType column names.
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#' pDat <- as(cov, "AnnotatedDataFrame")
#' exp_set <- ExpressionSet(assayData = kidexpr, phenoData = pDat)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- deSet(exp_set, null.model = null_model,
#' full.model = full_model)
#'
#' # optionally add individuals to experiment, in this case there are 36
#' # individuals that were sampled twice
#' indSamples <- as.factor(rep(1:36, each = 2))
#' de_obj <- deSet(exp_set, null.model = null_model,
#' full.model = full_model, ind = indSamples)
#' summary(de_obj)
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#'
#' @export
setGeneric("deSet", function(object, full.model, null.model,
individual=NULL) standardGeneric("deSet"))
#' Estimate the q-values for a given set of p-values
#'
#' Runs \code{\link{qvalue}} on a \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param ... Additional arguments for \code{\link{qvalue}}
#'
#' @return \code{\linkS4class{deSet}} object with slots updated by \code{\link{qvalue}}
#' calculations.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # Run lrt (or odp) and apply_qvalue
#' de_lrt <- lrt(de_obj)
#' de_lrt <- apply_qvalue(de_lrt, fdr.level = 0.05,
#' pi0.method = "bootstrap", adj=1.2)
#' summary(de_lrt)
#'
#' @references
#' Storey JD and Tibshirani R. (2003) Statistical significance for
#' genome-wide studies. Proceedings of the National Academy of Sciences,
#' 100: 9440-9445
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#'
#' @export
setGeneric("apply_qvalue", function(object, ...)
standardGeneric("apply_qvalue"))
#' Estimate surrogate variables
#'
#' Runs \code{\link{sva}} on the null and full models in
#' \code{\linkS4class{deSet}}. See \code{\link{sva}} for additional details.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param ... Additional arguments for \code{\link{sva}}
#'
#' @return \code{\linkS4class{deSet}} object where the surrogate variables
#' estimated by \code{\link{sva}} are added to the full model and null model
#' matrices.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run surrogate variable analysis
#' de_sva <- apply_sva(de_obj)
#'
#' # run odp/lrt with surrogate variables added
#' de_odp <- odp(de_sva, bs.its = 30)
#' summary(de_odp)
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @references
#' Leek JT, Storey JD (2007) Capturing Heterogeneity in Gene Expression
#' Studies by Surrogate Variable Analysis. PLoS Genet 3(9): e161.
#' doi:10.1371/journal.pgen.0030161
#'
#' Leek JT and Storey JD. (2008) A general framework for multiple testing
#' dependence. Proceedings of the National Academy of Sciences, 105: 18718-
#' 18723.
#'
#' @author John Storey, Jeffrey Leek, Andrew Bass
#' @export
setGeneric("apply_sva", function(object, ...)
standardGeneric("apply_sva"))
#' Supervised normalization of data in edge
#'
#' Runs \code{snm} on a deSet object based on the null and full models in
#' \code{\linkS4class{deSet}}. See \code{\link{snm}} for additional details
#' on the algorithm.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param int.var \code{data frame}: intensity-dependent effects (see
#' \code{\link{snm}} for details)
#' @param ... Additional arguments for \code{\link{snm}}
#'
#' @return \code{apply_snm} returns a \code{\linkS4class{deSet}} object where
#' assayData (the expression data) that has been passed to apply_snm is replaced
#' with the normalized data that \code{\link{snm}} returns. Specifically,
#' \code{exprs(object)} is replaced by \code{$norm.dat} from \code{\link{snm}},
#' where \code{object} is the \code{\link{deSet}} object.
#'
#' @references
#' Mechan BH, Nelson PS, Storey JD. Supervised normalization of microarrays.
#' Bioinformatics 2010;26:1308-1315.
#'
#' @examples
#' # simulate data
#' library(snm)
#' singleChannel <- sim.singleChannel(12345)
#' data <- singleChannel$raw.data
#'
#' # create deSet object using build_models (can use ExpressionSet see manual)
#' cov <- data.frame(grp = singleChannel$bio.var[,2])
#' full_model <- ~grp
#' null_model <- ~1
#'
#' # create deSet object using build_models
#' de_obj <- build_models(data = data, cov = cov, full.model = full_model,
#' null.model = null_model)
#'
#' # run snm using intensity-dependent adjustment variable
#' de_snm <- apply_snm(de_obj, int.var = singleChannel$int.var,
#' verbose = FALSE, num.iter = 1)
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#' @export
setGeneric("apply_snm", function(object, int.var=NULL, ...)
standardGeneric("apply_snm"))
#' Non-Parametric Jackstraw for Principal Component Analysis (PCA)
#'
#' Estimates statistical significance of association between variables and
#' their principal components (PCs).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param r1 a numeric vector of principal components of interest. Choose a subset of r significant PCs to be used.
#' @param r a number (a positive integer) of significant principal components.
#' @param s a number (a positive integer) of synthetic null variables. Out of m variables, s variables are independently permuted.
#' @param B a number (a positive integer) of resampling iterations. There will be a total of s*B null statistics.
#' @param covariate a data matrix of covariates with corresponding n observations.
#' @param verbose a logical indicator as to whether to print the progress.
#' @param seed a seed for the random number generator.
#'
#' @details
#' This function computes m p-values of linear association between m variables
#' and their PCs. Its resampling strategy accounts for the over-fitting
#' characteristics due to direct computation of PCs from the observed data
#' and protects against an anti-conservative bias.
#'
#' Provide the \code{\linkS4class{deSet}},
#' with m variables as rows and n observations as columns. Given that there are
#' r significant PCs, this function tests for linear association between m
#' varibles and their r PCs.
#'
#' You could specify a subset of significant PCs
#' that you are interested in r1. If PC is given, then this function computes
#' statistical significance of association between m variables and PC, while
#' adjusting for other PCs (i.e., significant PCs that are not your interest).
#' For example, if you want to identify variables associated with 1st and 2nd
#' PCs, when your data contains three significant PCs, set r=3 and r1=c(1,2).
#'
#' Please take a careful look at your data and use appropriate graphical and
#' statistical criteria to determine a number of significant PCs, r. The number
#' of significant PCs depends on the data structure and the context. In a case
#' when you fail to specify r, it will be estimated from a permutation test
#' (Buja and Eyuboglu, 1992) using a function \code{\link{permutationPA}}.
#'
#' If s is not supplied, s is set to about 10% of m variables. If B is not
#' supplied, B is set to m*10/s.
#'
#' @return \code{apply_jackstraw} returns a \code{list} containing the following
#' slots:
#' \itemize{
#' \item{\code{p.value} the m p-values of association tests between variables
#' and their principal components}
#' \item{\code{obs.stat} the observed F-test statistics}
#' \item{\code{null.stat} the s*B null F-test statistics}
#' }
#'
#'
#' @references
#' Chung and Storey (2013) Statistical Significance of
#' Variables Driving Systematic Variation in
#' High-Dimensional Data. arXiv:1308.6013 [stat.ME]
#' \url{http://arxiv.org/abs/1308.6013}
#'
#'More information available at \url{http://ncc.name/}
#'
#' @examples
# import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#' ## apply the jackstraw
#' out = apply_jackstraw(de_obj, r1=1, r=1)
#' ## Use optional arguments
#' ## For example, set s and B for a balance between speed of the algorithm and accuracy of p-values
#' ## out = apply_jackstraw(dat, r1=1, r=1, s=10, B=1000, seed=5678)
#'
#' @seealso \code{\link{permutationPA}}
#'
#' @author Neo Christopher Chung \email{nc@@princeton.edu}
#' @import jackstraw
#' @export
setGeneric("apply_jackstraw", function(object, r1 = NULL, r = NULL, s = NULL, B = NULL,
covariate = NULL, verbose = TRUE, seed = NULL)
standardGeneric("apply_jackstraw"))
#' Full model equation
#'
#' These generic functions access and set the full model for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{formula}: The experiment design for the full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract out the full model equation
#' mod_full <- fullModel(de_obj)
#'
#' # change the full model in the experiment
#' fullModel(de_obj) <- ~sex + ns(age, df = 2)
#'
#'
#' @return the formula for the full model.
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("fullModel", function(object) standardGeneric("fullModel"))
#' @rdname fullModel
#' @export
setGeneric("fullModel<-", function(object, value) {
standardGeneric("fullModel<-")
})
#' Null model equation from deSet object
#'
#' These generic functions access and set the null model for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{formula}: The experiment design for the null model.
#'
#' @return \code{nullModel} returns the formula for the null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the null model equation
#' mod_null <- nullModel(de_obj)
#'
#' # change null model in experiment but must update full model
#' nullModel(de_obj) <- ~1
#' fullModel(de_obj) <- ~1 + ns(age, df=4)
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @keywords nullModel, nullModel<-
#'
#' @exportMethod nullModel
setGeneric("nullModel", function(object) standardGeneric("nullModel"))
#' @rdname nullModel
#' @export
setGeneric("nullModel<-", function(object, value) {
standardGeneric("nullModel<-")
})
#' Matrix representation of null model
#'
#' These generic functions access and set the null matrix for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{matrix}: null model matrix where columns are covariates
#' and rows are observations
#'
#' @return \code{nullMatrix} returns the value of the null model matrix.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the null model as a matrix
#' mat_null <- nullMatrix(de_obj)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{fullModel}} and
#' \code{\link{fullModel}}
#'
#' @export
setGeneric("nullMatrix", function(object) standardGeneric("nullMatrix"))
#' @rdname nullMatrix
#' @export
setGeneric("nullMatrix<-", function(object, value) {
standardGeneric("nullMatrix<-")
})
#' Matrix representation of full model
#'
#' These generic functions access and set the full matrix for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{matrix}: full model matrix where the columns are the
#' covariates and rows are observations
#'
#' @return \code{fullMatrix} returns the value of the full model matrix.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the full model equation as a matrix
#' mat_full <- fullMatrix(de_obj)
#' @author Andrew Bass, John Storey
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{fullModel}}
#'
#' @export
setGeneric("fullMatrix", function(object) standardGeneric("fullMatrix"))
#' @rdname fullMatrix
#' @export
setGeneric("fullMatrix<-", function(object, value) {
standardGeneric("fullMatrix<-")
})
#' Access/set qvalue slot
#'
#' These generic functions access and set the \code{qvalue} object in the
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value S3 \code{object}: \code{\link{qvalue}}
#'
#' @return \code{qvalueObj} returns a \code{\link{qvalue}} object.
#'
#' @examples
#' # import data
#' library(splines)
#' library(qvalue)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run the odp method
#' de_odp <- odp(de_obj, bs.its = 20)
#'
#' # extract out significance results
#' qval_obj <- qvalueObj(de_odp)
#'
#' # run qvalue and assign it to deSet slot
#' pvals <- qval_obj$pvalues
#' qval_new <- qvalue(pvals, pfdr = TRUE, fdr.level = 0.1)
#' qvalueObj(de_odp) <- qval_new
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{lrt}}, \code{\link{odp}} and
#' \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("qvalueObj", function(object) standardGeneric("qvalueObj"))
#' @rdname qvalueObj
#' @export
setGeneric("qvalueObj<-", function(object, value) {
standardGeneric("qvalueObj<-")
})
#' Individuals sampled in experiment
#'
#' These generic functions access and set the individual slot in
#' \code{\linkS4class{deSet}}.
#'
#' @param object \code{\linkS4class{deSet}}
#' @param value \code{factor}: Identifies which samples correspond to which
#' individuals. Important if the same individuals are sampled multiple times
#' in a longitudinal fashion.
#'
#' @return \code{individual} returns information regarding dinstinct individuals
#' sampled in the experiment.
#'
#' @examples
#' library(splines)
#' # import data
#' data(endotoxin)
#' ind <- endotoxin$ind
#' time <- endotoxin$time
#' class <- endotoxin$class
#' endoexpr <- endotoxin$endoexpr
#' cov <- data.frame(individual = ind, time = time, class = class)
#'
#' # create ExpressionSet object
#' pDat <- as(cov, "AnnotatedDataFrame")
#' exp_set <- ExpressionSet(assayData = endoexpr, phenoData = pDat)
#'
#' # formulate null and full models in experiement
#' # note: interaction term is a way of taking into account group effects
#' mNull <- ~ns(time, df=4, intercept = FALSE)
#' mFull <- ~ns(time, df=4, intercept = FALSE) +
#' ns(time, df=4, intercept = FALSE):class + class
#'
#' # create deSet object
#' de_obj <- deSet(exp_set, full.model = mFull, null.model = mNull,
#' individual = ind)
#'
#' # extract out the individuals factor
#' ind_exp <- individual(de_obj)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("individual", function(object) standardGeneric("individual"))
#' @rdname individual
#' @export
setGeneric("individual<-", function(object, value) {
standardGeneric("individual<-")
})
#' Regression coefficients from full model fit
#'
#' Access the full model fitted coefficients of a
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @return \code{betaCoef} returns the regression coefficients for the full
#' model fit.
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract beta coefficients
#' beta <- betaCoef(de_fit)
#'
#' @export
setGeneric("betaCoef", function(object) standardGeneric("betaCoef"))
#' Statistic type used in analysis
#'
#' Access the statistic type in a \code{\linkS4class{deFit}} object. Can
#' either be the optimal discovery procedure (odp) or the likelihood ratio
#' test (lrt).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract the statistic type of model fits
#' stat_type <- sType(de_fit)
#'
#' @return \code{sType} returns the statistic type- either "odp" or "lrt".
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}, \code{\linkS4class{deFit}} and
#' \code{\linkS4class{deSet}}
#'
#' @keywords sType
#'
#' @exportMethod sType
setGeneric("sType", function(object) standardGeneric("sType"))
#' Fitted data from the full model
#'
#' Access the fitted data from the full model in a
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage fitFull(object)
#'
#' @return \code{fitFull} returns a matrix of fitted values from full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract fitted values for full model
#' fitted_full <- fitFull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("fitFull", function(object) standardGeneric("fitFull"))
#' Fitted data from the null model
#'
#' Access the fitted data from the null model in an
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage fitNull(object)
#'
#' @return \code{fitNull} returns a matrix of fitted values from null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract fitted values from null model
#' fitted_null <- fitNull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("fitNull", function(object) standardGeneric("fitNull"))
#' Residuals of full model fit
#'
#' Access the fitted full model residuals in an \code{\linkS4class{deFit}}
#' object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage resFull(object)
#'
#' @return \code{resFull} returns a matrix of residuals from full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract out the full residuals from the model fit
#' res_full <- resFull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("resFull", function(object) standardGeneric("resFull"))
#' Residuals of null model fit
#'
#' Access the fitted null model residuals in an \code{\linkS4class{deFit}}
#' object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage resNull(object)
#'
#' @return \code{resNull} returns a matrix of residuals from null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract out the null residuals from the model fits
#' res_null <- resNull(de_fit)
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("resNull", function(object) standardGeneric("resNull"))
#' Summary of deFit and deSet
#'
#' Summary of \code{\linkS4class{deFit}} and \code{\linkS4class{deSet}} objects.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param \dots additional parameters
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # get summary
#' summary(de_obj)
#'
#' # run odp and summarize
#' de_odp <- odp(de_obj, bs.its= 20)
#' summary(de_odp)
#' @author John Storey, Andrew Bass
#'
#' @return
#' Summary of \code{\linkS4class{deSet}} object
#'
#' @keywords summary
#'
#' @export summary
setGeneric("summary")
#' Show function for deFit and deSet
#'
#' Show function for \code{\linkS4class{deFit}} and \code{\linkS4class{deSet}}
#' objects.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param \dots additional parameters
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # get summary
#' summary(de_obj)
#'
#' # run odp and summarize
#' de_odp <- odp(de_obj, bs.its= 20)
#' de_odp
#' @author John Storey, Andrew Bass
#'
#' @return
#' Nothing of interest
#'
#' @export
setGeneric("show")
StoreyLab/edge documentation built on May 9, 2019, 3:09 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#' Performs F-test (likelihood ratio test using Normal likelihood)
#'
#' \code{lrt} performs a generalized likelihood ratio test using the full and
#' null models.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}. Optional.
#' @param nullDistn \code{character}: either "normal" or "bootstrap", If
#' "normal" then the p-values are calculated using the F distribution. If
#' "bootstrap" then a bootstrap algorithm is implemented to simulate
#' statistics from the null distribution. In the "bootstrap" case, empirical
#' p-values are calculated using the observed and null statistics (see
#' \code{\link{empPvals}}). Default is "normal".
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#' @param bs.its \code{integer}: number of null statistics generated (only
#' applicable for "bootstrap" method). Default is 100.
#' @param seed \code{integer}: set the seed value. Default is NULL.
#' @param verbose \code{boolean}: print iterations for bootstrap method.
#' Default is TRUE.
#' @param mod.F \code{boolean}: Moderated F-test, recommended for experiments
#' with a small sample size. Default is FALSE.
#' @param ... Additional arguments for \code{\link{apply_qvalue}} and
#' \code{\link{empPvals}} function.
#'
#' @details \code{lrt} fits the full and null models to each gene using the
#' function \code{\link{fit_models}} and then performs a likelihood ratio test.
#' The user has the option to calculate p-values a Normal distribution
#' assumption or through a bootstrap algorithm. If \code{nullDistn} is
#' "bootstrap" then empirical p-values will be determined from the
#' \code{\link{qvalue}} package (see \code{\link{empPvals}}).
#'
#' @author John Storey, Andrew Bass
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # lrt method
#' de_lrt <- lrt(de_obj, nullDistn = "normal")
#'
#' # to generate p-values from bootstrap
#' de_lrt <- lrt(de_obj, nullDistn = "bootstrap", bs.its = 30)
#'
#' # input a deFit object directly
#' de_fit <- fit_models(de_obj, stat.type = "lrt")
#' de_lrt <- lrt(de_obj, de.fit = de_fit)
#'
#' # summarize object
#' summary(de_lrt)
#'
#' @references
#' Storey JD, Xiao W, Leek JT, Tompkins RG, and Davis RW. (2005) Significance
#' analysis of time course microarray experiments. Proceedings of the National
#' Academy of Sciences, 102: 12837-12842.
#'
#' \url{http://en.wikipedia.org/wiki/Likelihood-ratio_test}
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{build_models}},
#' \code{\link{odp}}
#'
#' @export
setGeneric("lrt", function(object, de.fit,
nullDistn = c("normal","bootstrap"), weights = NULL,
bs.its = 100, seed = NULL, verbose = TRUE,
mod.F = FALSE, ...)
standardGeneric("lrt"))
#' The optimal discovery procedure
#'
#' \code{odp} performs the optimal discovery procedure, which is a framework for
#' optimally performing many hypothesis tests in a high-dimensional study. When
#' testing whether a feature is significant, the optimal discovery procedure
#' uses information across all features when testing for significance.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}. Optional.
#' @param odp.parms \code{list}: parameters for each cluster. See
#' \code{\link{kl_clust}}.
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#' @param bs.its \code{numeric}: number of null bootstrap iterations. Default
#' is 100.
#' @param n.mods \code{integer}: number of clusters used in
#' \code{\link{kl_clust}}. Default is 50.
#' @param seed \code{integer}: set the seed value. Default is NULL.
#' @param verbose \code{boolean}: print iterations for bootstrap method.
#' Default is TRUE.
#' @param ... Additional arguments for \code{\link{qvalue}} and
#' \code{\link{empPvals}}.
#'
#'
#' @details
#' The full ODP estimator computationally grows quadratically with respect to
#' the number of genes. This becomes computationally taxing at a certain point.
#' Therefore, an alternative method called mODP is used which has been shown to
#' provide results that are very similar. mODP utilizes a clustering algorithm
#' where genes are assigned to a cluster based on the Kullback-Leiber distance.
#' Each gene is assigned an module-average parameter to calculate the ODP score
#' and it reduces the computations time to approximately linear (see Woo, Leek
#' and Storey 2010). If the number of clusters is equal to the number of genes
#' then the original ODP is implemented. Depending on the number of hypothesis
#' tests, this can take some time.
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov,
#' null.model = null_model, full.model = full_model)
#'
#' # odp method
#' de_odp <- odp(de_obj, bs.its = 30)
#'
#' # input a deFit object or ODP parameters ... not necessary
#' de_fit <- fit_models(de_obj, stat.type = "odp")
#' de_clust <- kl_clust(de_obj, n.mods = 10)
#' de_odp <- odp(de_obj, de.fit = de_fit, odp.parms = de_clust,
#' bs.its = 30)
#'
#' # summarize object
#' summary(de_odp)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Woo S, Leek JT, Storey JD (2010) A computationally efficient modular
#' optimal discovery procedure. Bioinformatics, 27(4): 509-515.
#'
#' @author John Storey, Jeffrey Leek, Andrew Bass
#'
#' @seealso \code{\link{kl_clust}}, \code{\link{build_models}} and
#' \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("odp", function(object, de.fit, odp.parms = NULL, weights = NULL, bs.its = 100,
n.mods = 50, seed = NULL, verbose = TRUE, ...)
standardGeneric("odp"))
#' Modular optimal discovery procedure (mODP)
#'
#' \code{kl_clust} is an implementation of mODP that assigns genes to modules
#' based on of the Kullback-Leibler distance.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param de.fit \code{S4 object}: \code{\linkS4class{deFit}}.
#' @param n.mods \code{integer}: number of modules (i.e., clusters).
#'
#' @details mODP utilizes a k-means clustering algorithm where genes are
#' assigned to a cluster based on the Kullback-Leiber distance. Each gene is
#' assigned an module-average parameter to calculate the ODP score (See Woo,
#' Leek and Storey 2010 for more details). The mODP and full ODP produce nearly
#' exact results but mODP has the advantage of being computationally
#' faster.
#'
#' @note The results are generally insensitive to the number of modules after a
#' certain threshold of about n.mods>=50 in our experience. It is recommended
#' that users experiment with the number of modules. If the number of modules
#' is equal to the number of genes then the original ODP is implemented.
#'
#' @return
#' A list with the following slots:
#' \itemize{
#' \item {mu.full: cluster averaged fitted values from full model.}
#' \item {mu.null: cluster averaged fitted values from null model.}
#' \item {sig.full: cluster standard deviations from full model.}
#' \item {sig.null: cluster standard deviations from null model.}
#' \item {n.per.mod: total members in each cluster.}
#' \item {clustMembers: cluster membership for each gene.}
#' }
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # mODP method
#' de_clust <- kl_clust(de_obj)
#'
#' # change the number of clusters
#' de_clust <- kl_clust(de_obj, n.mods = 10)
#'
#' # input a deFit object
#' de_fit <- fit_models(de_obj, stat.type = "odp")
#' de_clust <- kl_clust(de_obj, de.fit = de_fit)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Woo S, Leek JT, Storey JD (2010) A computationally efficient modular optimal
#' discovery procedure. Bioinformatics, 27(4): 509-515.
#'
#' @author John Storey, Jeffrey Leek
#'
#' @seealso \code{\link{odp}}, \code{\link{fit_models}}
#'
#' @exportMethod kl_clust
setGeneric("kl_clust", function(object, de.fit = NULL, n.mods = 50)
standardGeneric("kl_clust"))
#' Linear regression of the null and full models
#'
#' \code{fit_models} fits a model matrix to each gene by using the least
#' squares method. Model fits can be either statistic type "odp" (optimal
#' discovery procedure) or "lrt" (likelihood ratio test).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}.
#' @param stat.type \code{character}: type of statistic to be used. Either
#' "lrt" or "odp". Default is "lrt".
#' @param weights \code{matrix}: weights for each observation. Default is NULL.
#'
#' @details
#' If "odp" method is implemented then the null model is removed from the full
#' model (see Storey 2007). Otherwise, the statistic type has no affect on the
#' model fit.
#'
#' @note \code{fit_models} does not have to be called by the user to use
#' \code{\link{odp}}, \code{\link{lrt}} or \code{\link{kl_clust}} as it is an
#' optional input and is implemented in the methods. The
#' \code{\linkS4class{deFit}} object can be created by the user if a different
#' statistical implementation is required.
#'
#' @return \code{\linkS4class{deFit}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # retrieve statistics from linear regression for each gene
#' fit_lrt <- fit_models(de_obj, stat.type = "lrt") # lrt method
#' fit_odp <- fit_models(de_obj, stat.type = "odp") # odp method
#'
#' # summarize object
#' summary(fit_odp)
#'
#' @references
#' Storey JD. (2007) The optimal discovery procedure: A new approach to
#' simultaneous significance testing. Journal of the Royal Statistical
#' Society, Series B, 69: 347-368.
#'
#' Storey JD, Dai JY, and Leek JT. (2007) The optimal discovery procedure for
#' large-scale significance testing, with applications to comparative
#' microarray experiments. Biostatistics, 8: 414-432.
#'
#' Storey JD, Xiao W, Leek JT, Tompkins RG, and Davis RW. (2005) Significance
#' analysis of time course microarray experiments. Proceedings of the National
#' Academy of Sciences, 102: 12837-12842.
#'
#' @seealso \code{\linkS4class{deFit}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey
#' @exportMethod fit_models
setGeneric("fit_models",
function(object, stat.type = c("lrt", "odp"), weights = NULL) {
standardGeneric("fit_models")
})
#' Create a deSet object from an ExpressionSet
#'
#' Creates a \code{\linkS4class{deSet}} object that extends the
#' \code{\link{ExpressionSet}} object.
#'
#' @param object \code{S4 object}: \code{\link{ExpressionSet}}
#' @param full.model \code{formula}: full model containing the both the
#' adjustment and the biological variables for the experiment.
#' @param null.model \code{formula}: null model containing the adjustment
#' variables for the experiment.
#' @param individual \code{factor}: information on repeated samples in
#' experiment.
#'
#' @note It is essential that the null and full models have the same variables
#' as the ExpressionSet phenoType column names.
#'
#' @return \code{\linkS4class{deSet}} object
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#' pDat <- as(cov, "AnnotatedDataFrame")
#' exp_set <- ExpressionSet(assayData = kidexpr, phenoData = pDat)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- deSet(exp_set, null.model = null_model,
#' full.model = full_model)
#'
#' # optionally add individuals to experiment, in this case there are 36
#' # individuals that were sampled twice
#' indSamples <- as.factor(rep(1:36, each = 2))
#' de_obj <- deSet(exp_set, null.model = null_model,
#' full.model = full_model, ind = indSamples)
#' summary(de_obj)
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#'
#' @export
setGeneric("deSet", function(object, full.model, null.model,
individual=NULL) standardGeneric("deSet"))
#' Estimate the q-values for a given set of p-values
#'
#' Runs \code{\link{qvalue}} on a \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param ... Additional arguments for \code{\link{qvalue}}
#'
#' @return \code{\linkS4class{deSet}} object with slots updated by \code{\link{qvalue}}
#' calculations.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # Run lrt (or odp) and apply_qvalue
#' de_lrt <- lrt(de_obj)
#' de_lrt <- apply_qvalue(de_lrt, fdr.level = 0.05,
#' pi0.method = "bootstrap", adj=1.2)
#' summary(de_lrt)
#'
#' @references
#' Storey JD and Tibshirani R. (2003) Statistical significance for
#' genome-wide studies. Proceedings of the National Academy of Sciences,
#' 100: 9440-9445
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#'
#' @export
setGeneric("apply_qvalue", function(object, ...)
standardGeneric("apply_qvalue"))
#' Estimate surrogate variables
#'
#' Runs \code{\link{sva}} on the null and full models in
#' \code{\linkS4class{deSet}}. See \code{\link{sva}} for additional details.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param ... Additional arguments for \code{\link{sva}}
#'
#' @return \code{\linkS4class{deSet}} object where the surrogate variables
#' estimated by \code{\link{sva}} are added to the full model and null model
#' matrices.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run surrogate variable analysis
#' de_sva <- apply_sva(de_obj)
#'
#' # run odp/lrt with surrogate variables added
#' de_odp <- odp(de_sva, bs.its = 30)
#' summary(de_odp)
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @references
#' Leek JT, Storey JD (2007) Capturing Heterogeneity in Gene Expression
#' Studies by Surrogate Variable Analysis. PLoS Genet 3(9): e161.
#' doi:10.1371/journal.pgen.0030161
#'
#' Leek JT and Storey JD. (2008) A general framework for multiple testing
#' dependence. Proceedings of the National Academy of Sciences, 105: 18718-
#' 18723.
#'
#' @author John Storey, Jeffrey Leek, Andrew Bass
#' @export
setGeneric("apply_sva", function(object, ...)
standardGeneric("apply_sva"))
#' Supervised normalization of data in edge
#'
#' Runs \code{snm} on a deSet object based on the null and full models in
#' \code{\linkS4class{deSet}}. See \code{\link{snm}} for additional details
#' on the algorithm.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param int.var \code{data frame}: intensity-dependent effects (see
#' \code{\link{snm}} for details)
#' @param ... Additional arguments for \code{\link{snm}}
#'
#' @return \code{apply_snm} returns a \code{\linkS4class{deSet}} object where
#' assayData (the expression data) that has been passed to apply_snm is replaced
#' with the normalized data that \code{\link{snm}} returns. Specifically,
#' \code{exprs(object)} is replaced by \code{$norm.dat} from \code{\link{snm}},
#' where \code{object} is the \code{\link{deSet}} object.
#'
#' @references
#' Mechan BH, Nelson PS, Storey JD. Supervised normalization of microarrays.
#' Bioinformatics 2010;26:1308-1315.
#'
#' @examples
#' # simulate data
#' library(snm)
#' singleChannel <- sim.singleChannel(12345)
#' data <- singleChannel$raw.data
#'
#' # create deSet object using build_models (can use ExpressionSet see manual)
#' cov <- data.frame(grp = singleChannel$bio.var[,2])
#' full_model <- ~grp
#' null_model <- ~1
#'
#' # create deSet object using build_models
#' de_obj <- build_models(data = data, cov = cov, full.model = full_model,
#' null.model = null_model)
#'
#' # run snm using intensity-dependent adjustment variable
#' de_snm <- apply_snm(de_obj, int.var = singleChannel$int.var,
#' verbose = FALSE, num.iter = 1)
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{odp}} and
#' \code{\link{lrt}}
#'
#' @author John Storey, Andrew Bass
#' @export
setGeneric("apply_snm", function(object, int.var=NULL, ...)
standardGeneric("apply_snm"))
#' Non-Parametric Jackstraw for Principal Component Analysis (PCA)
#'
#' Estimates statistical significance of association between variables and
#' their principal components (PCs).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param r1 a numeric vector of principal components of interest. Choose a subset of r significant PCs to be used.
#' @param r a number (a positive integer) of significant principal components.
#' @param s a number (a positive integer) of synthetic null variables. Out of m variables, s variables are independently permuted.
#' @param B a number (a positive integer) of resampling iterations. There will be a total of s*B null statistics.
#' @param covariate a data matrix of covariates with corresponding n observations.
#' @param verbose a logical indicator as to whether to print the progress.
#' @param seed a seed for the random number generator.
#'
#' @details
#' This function computes m p-values of linear association between m variables
#' and their PCs. Its resampling strategy accounts for the over-fitting
#' characteristics due to direct computation of PCs from the observed data
#' and protects against an anti-conservative bias.
#'
#' Provide the \code{\linkS4class{deSet}},
#' with m variables as rows and n observations as columns. Given that there are
#' r significant PCs, this function tests for linear association between m
#' varibles and their r PCs.
#'
#' You could specify a subset of significant PCs
#' that you are interested in r1. If PC is given, then this function computes
#' statistical significance of association between m variables and PC, while
#' adjusting for other PCs (i.e., significant PCs that are not your interest).
#' For example, if you want to identify variables associated with 1st and 2nd
#' PCs, when your data contains three significant PCs, set r=3 and r1=c(1,2).
#'
#' Please take a careful look at your data and use appropriate graphical and
#' statistical criteria to determine a number of significant PCs, r. The number
#' of significant PCs depends on the data structure and the context. In a case
#' when you fail to specify r, it will be estimated from a permutation test
#' (Buja and Eyuboglu, 1992) using a function \code{\link{permutationPA}}.
#'
#' If s is not supplied, s is set to about 10% of m variables. If B is not
#' supplied, B is set to m*10/s.
#'
#' @return \code{apply_jackstraw} returns a \code{list} containing the following
#' slots:
#' \itemize{
#' \item{\code{p.value} the m p-values of association tests between variables
#' and their principal components}
#' \item{\code{obs.stat} the observed F-test statistics}
#' \item{\code{null.stat} the s*B null F-test statistics}
#' }
#'
#'
#' @references
#' Chung and Storey (2013) Statistical Significance of
#' Variables Driving Systematic Variation in
#' High-Dimensional Data. arXiv:1308.6013 [stat.ME]
#' \url{http://arxiv.org/abs/1308.6013}
#'
#'More information available at \url{http://ncc.name/}
#'
#' @examples
# import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#' ## apply the jackstraw
#' out = apply_jackstraw(de_obj, r1=1, r=1)
#' ## Use optional arguments
#' ## For example, set s and B for a balance between speed of the algorithm and accuracy of p-values
#' ## out = apply_jackstraw(dat, r1=1, r=1, s=10, B=1000, seed=5678)
#'
#' @seealso \code{\link{permutationPA}}
#'
#' @author Neo Christopher Chung \email{nc@@princeton.edu}
#' @import jackstraw
#' @export
setGeneric("apply_jackstraw", function(object, r1 = NULL, r = NULL, s = NULL, B = NULL,
covariate = NULL, verbose = TRUE, seed = NULL)
standardGeneric("apply_jackstraw"))
#' Full model equation
#'
#' These generic functions access and set the full model for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{formula}: The experiment design for the full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract out the full model equation
#' mod_full <- fullModel(de_obj)
#'
#' # change the full model in the experiment
#' fullModel(de_obj) <- ~sex + ns(age, df = 2)
#'
#'
#' @return the formula for the full model.
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("fullModel", function(object) standardGeneric("fullModel"))
#' @rdname fullModel
#' @export
setGeneric("fullModel<-", function(object, value) {
standardGeneric("fullModel<-")
})
#' Null model equation from deSet object
#'
#' These generic functions access and set the null model for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{formula}: The experiment design for the null model.
#'
#' @return \code{nullModel} returns the formula for the null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the null model equation
#' mod_null <- nullModel(de_obj)
#'
#' # change null model in experiment but must update full model
#' nullModel(de_obj) <- ~1
#' fullModel(de_obj) <- ~1 + ns(age, df=4)
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @keywords nullModel, nullModel<-
#'
#' @exportMethod nullModel
setGeneric("nullModel", function(object) standardGeneric("nullModel"))
#' @rdname nullModel
#' @export
setGeneric("nullModel<-", function(object, value) {
standardGeneric("nullModel<-")
})
#' Matrix representation of null model
#'
#' These generic functions access and set the null matrix for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{matrix}: null model matrix where columns are covariates
#' and rows are observations
#'
#' @return \code{nullMatrix} returns the value of the null model matrix.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the null model as a matrix
#' mat_null <- nullMatrix(de_obj)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{fullModel}} and
#' \code{\link{fullModel}}
#'
#' @export
setGeneric("nullMatrix", function(object) standardGeneric("nullMatrix"))
#' @rdname nullMatrix
#' @export
setGeneric("nullMatrix<-", function(object, value) {
standardGeneric("nullMatrix<-")
})
#' Matrix representation of full model
#'
#' These generic functions access and set the full matrix for
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value \code{matrix}: full model matrix where the columns are the
#' covariates and rows are observations
#'
#' @return \code{fullMatrix} returns the value of the full model matrix.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # extract the full model equation as a matrix
#' mat_full <- fullMatrix(de_obj)
#' @author Andrew Bass, John Storey
#'
#' @seealso \code{\linkS4class{deSet}}, \code{\link{fullModel}}
#'
#' @export
setGeneric("fullMatrix", function(object) standardGeneric("fullMatrix"))
#' @rdname fullMatrix
#' @export
setGeneric("fullMatrix<-", function(object, value) {
standardGeneric("fullMatrix<-")
})
#' Access/set qvalue slot
#'
#' These generic functions access and set the \code{qvalue} object in the
#' \code{\linkS4class{deSet}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param value S3 \code{object}: \code{\link{qvalue}}
#'
#' @return \code{qvalueObj} returns a \code{\link{qvalue}} object.
#'
#' @examples
#' # import data
#' library(splines)
#' library(qvalue)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run the odp method
#' de_odp <- odp(de_obj, bs.its = 20)
#'
#' # extract out significance results
#' qval_obj <- qvalueObj(de_odp)
#'
#' # run qvalue and assign it to deSet slot
#' pvals <- qval_obj$pvalues
#' qval_new <- qvalue(pvals, pfdr = TRUE, fdr.level = 0.1)
#' qvalueObj(de_odp) <- qval_new
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{lrt}}, \code{\link{odp}} and
#' \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("qvalueObj", function(object) standardGeneric("qvalueObj"))
#' @rdname qvalueObj
#' @export
setGeneric("qvalueObj<-", function(object, value) {
standardGeneric("qvalueObj<-")
})
#' Individuals sampled in experiment
#'
#' These generic functions access and set the individual slot in
#' \code{\linkS4class{deSet}}.
#'
#' @param object \code{\linkS4class{deSet}}
#' @param value \code{factor}: Identifies which samples correspond to which
#' individuals. Important if the same individuals are sampled multiple times
#' in a longitudinal fashion.
#'
#' @return \code{individual} returns information regarding dinstinct individuals
#' sampled in the experiment.
#'
#' @examples
#' library(splines)
#' # import data
#' data(endotoxin)
#' ind <- endotoxin$ind
#' time <- endotoxin$time
#' class <- endotoxin$class
#' endoexpr <- endotoxin$endoexpr
#' cov <- data.frame(individual = ind, time = time, class = class)
#'
#' # create ExpressionSet object
#' pDat <- as(cov, "AnnotatedDataFrame")
#' exp_set <- ExpressionSet(assayData = endoexpr, phenoData = pDat)
#'
#' # formulate null and full models in experiement
#' # note: interaction term is a way of taking into account group effects
#' mNull <- ~ns(time, df=4, intercept = FALSE)
#' mFull <- ~ns(time, df=4, intercept = FALSE) +
#' ns(time, df=4, intercept = FALSE):class + class
#'
#' # create deSet object
#' de_obj <- deSet(exp_set, full.model = mFull, null.model = mNull,
#' individual = ind)
#'
#' # extract out the individuals factor
#' ind_exp <- individual(de_obj)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\linkS4class{deSet}}
#'
#' @export
setGeneric("individual", function(object) standardGeneric("individual"))
#' @rdname individual
#' @export
setGeneric("individual<-", function(object, value) {
standardGeneric("individual<-")
})
#' Regression coefficients from full model fit
#'
#' Access the full model fitted coefficients of a
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @return \code{betaCoef} returns the regression coefficients for the full
#' model fit.
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract beta coefficients
#' beta <- betaCoef(de_fit)
#'
#' @export
setGeneric("betaCoef", function(object) standardGeneric("betaCoef"))
#' Statistic type used in analysis
#'
#' Access the statistic type in a \code{\linkS4class{deFit}} object. Can
#' either be the optimal discovery procedure (odp) or the likelihood ratio
#' test (lrt).
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract the statistic type of model fits
#' stat_type <- sType(de_fit)
#'
#' @return \code{sType} returns the statistic type- either "odp" or "lrt".
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}, \code{\linkS4class{deFit}} and
#' \code{\linkS4class{deSet}}
#'
#' @keywords sType
#'
#' @exportMethod sType
setGeneric("sType", function(object) standardGeneric("sType"))
#' Fitted data from the full model
#'
#' Access the fitted data from the full model in a
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage fitFull(object)
#'
#' @return \code{fitFull} returns a matrix of fitted values from full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract fitted values for full model
#' fitted_full <- fitFull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("fitFull", function(object) standardGeneric("fitFull"))
#' Fitted data from the null model
#'
#' Access the fitted data from the null model in an
#' \code{\linkS4class{deFit}} object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage fitNull(object)
#'
#' @return \code{fitNull} returns a matrix of fitted values from null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract fitted values from null model
#' fitted_null <- fitNull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("fitNull", function(object) standardGeneric("fitNull"))
#' Residuals of full model fit
#'
#' Access the fitted full model residuals in an \code{\linkS4class{deFit}}
#' object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage resFull(object)
#'
#' @return \code{resFull} returns a matrix of residuals from full model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract out the full residuals from the model fit
#' res_full <- resFull(de_fit)
#'
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("resFull", function(object) standardGeneric("resFull"))
#' Residuals of null model fit
#'
#' Access the fitted null model residuals in an \code{\linkS4class{deFit}}
#' object.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deFit}}
#'
#' @usage resNull(object)
#'
#' @return \code{resNull} returns a matrix of residuals from null model.
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # run fit_models to get model fits
#' de_fit <- fit_models(de_obj)
#'
#' # extract out the null residuals from the model fits
#' res_null <- resNull(de_fit)
#' @author John Storey, Andrew Bass
#'
#' @seealso \code{\link{fit_models}}
#'
#' @export
setGeneric("resNull", function(object) standardGeneric("resNull"))
#' Summary of deFit and deSet
#'
#' Summary of \code{\linkS4class{deFit}} and \code{\linkS4class{deSet}} objects.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param \dots additional parameters
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # get summary
#' summary(de_obj)
#'
#' # run odp and summarize
#' de_odp <- odp(de_obj, bs.its= 20)
#' summary(de_odp)
#' @author John Storey, Andrew Bass
#'
#' @return
#' Summary of \code{\linkS4class{deSet}} object
#'
#' @keywords summary
#'
#' @export summary
setGeneric("summary")
#' Show function for deFit and deSet
#'
#' Show function for \code{\linkS4class{deFit}} and \code{\linkS4class{deSet}}
#' objects.
#'
#' @param object \code{S4 object}: \code{\linkS4class{deSet}}
#' @param \dots additional parameters
#'
#' @examples
#' # import data
#' library(splines)
#' data(kidney)
#' age <- kidney$age
#' sex <- kidney$sex
#' kidexpr <- kidney$kidexpr
#' cov <- data.frame(sex = sex, age = age)
#'
#' # create models
#' null_model <- ~sex
#' full_model <- ~sex + ns(age, df = 4)
#'
#' # create deSet object from data
#' de_obj <- build_models(data = kidexpr, cov = cov, null.model = null_model,
#' full.model = full_model)
#'
#' # get summary
#' summary(de_obj)
#'
#' # run odp and summarize
#' de_odp <- odp(de_obj, bs.its= 20)
#' de_odp
#' @author John Storey, Andrew Bass
#'
#' @return
#' Nothing of interest
#'
#' @export
setGeneric("show")
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.