calculateDMN {mia} | R Documentation |
These functions are accessors for functions implemented in the
DirichletMultinomial
package
calculateDMN(x, ...) ## S4 method for signature 'ANY' calculateDMN( x, k = 1, BPPARAM = SerialParam(), seed = runif(1, 0, .Machine$integer.max), ... ) ## S4 method for signature 'SummarizedExperiment' calculateDMN( x, abund_values = exprs_values, exprs_values = "counts", transposed = FALSE, ... ) runDMN(x, name = "DMN", ...) getDMN(x, name = "DMN", ...) ## S4 method for signature 'SummarizedExperiment' getDMN(x, name = "DMN") bestDMNFit(x, name = "DMN", type = c("laplace", "AIC", "BIC"), ...) ## S4 method for signature 'SummarizedExperiment' bestDMNFit(x, name = "DMN", type = c("laplace", "AIC", "BIC")) getBestDMNFit(x, name = "DMN", type = c("laplace", "AIC", "BIC"), ...) ## S4 method for signature 'SummarizedExperiment' getBestDMNFit(x, name = "DMN", type = c("laplace", "AIC", "BIC")) calculateDMNgroup(x, ...) ## S4 method for signature 'ANY' calculateDMNgroup( x, variable, k = 1, seed = runif(1, 0, .Machine$integer.max), ... ) ## S4 method for signature 'SummarizedExperiment' calculateDMNgroup( x, variable, abund_values = exprs_values, exprs_values = "counts", transposed = FALSE, ... ) performDMNgroupCV(x, ...) ## S4 method for signature 'ANY' performDMNgroupCV( x, variable, k = 1, seed = runif(1, 0, .Machine$integer.max), ... ) ## S4 method for signature 'SummarizedExperiment' performDMNgroupCV( x, variable, abund_values = exprs_values, exprs_values = "counts", transposed = FALSE, ... )
x |
a numeric matrix with samples as rows or a
|
... |
optional arguments not used. |
k |
the number of Dirichlet components to fit. See
|
BPPARAM |
A
|
seed |
random number seed. See
|
abund_values |
a single |
exprs_values |
a single |
transposed |
Logical scalar, is x transposed with samples in rows? |
name |
the name to store the result in
|
type |
the type of measure used for the goodness of fit. One of ‘laplace’, ‘AIC’ or ‘BIC’. |
variable |
a variable from |
calculateDMN
and getDMN
return a list of DMN
objects,
one element for each value of k provided.
bestDMNFit
returns the index for the best fit and getBestDMNFit
returns a single DMN
object.
calculateDMNgroup
returns a
DMNGroup
object
performDMNgroupCV
returns a data.frame
DMN-class
,
DMNGroup-class
,
dmn
,
dmngroup
,
cvdmngroup
,
accessors for DMN objects
fl <- system.file(package="DirichletMultinomial", "extdata", "Twins.csv") counts <- as.matrix(read.csv(fl, row.names=1)) fl <- system.file(package="DirichletMultinomial", "extdata", "TwinStudy.t") pheno0 <- scan(fl) lvls <- c("Lean", "Obese", "Overwt") pheno <- factor(lvls[pheno0 + 1], levels=lvls) colData <- DataFrame(pheno = pheno) se <- SummarizedExperiment(assays = list(counts = counts), colData = colData) # dmn <- calculateDMN(se) dmn[[1L]] # since this take a bit of resources to calculate for k > 1, the data is # loaded ## Not run: se <- runDMN(se, name = "DMN", k = 1:7) ## End(Not run) data(dmn_se) names(metadata(dmn_se)) # return a list of DMN objects getDMN(dmn_se) # return, which objects fits best bestDMNFit(dmn_se, type = "laplace") # return the model, which fits best getBestDMNFit(dmn_se, type = "laplace")