fitFeatureModel {metagenomeSeq} | R Documentation |
Wrapper to actually run zero-inflated log-normal model given a MRexperiment object and model matrix. User can decide to shrink parameter estimates.
fitFeatureModel(obj, mod, coef = 2, B = 1, szero = FALSE, spos = TRUE)
obj |
A MRexperiment object with count data. |
mod |
The model for the count distribution. |
coef |
Coefficient of interest to grab log fold-changes. |
B |
Number of bootstraps to perform if >1. If >1 performs permutation test. |
szero |
TRUE/FALSE, shrink zero component parameters. |
spos |
TRUE/FALSE, shrink positive component parameters. |
A list of objects including:
call - the call made to fitFeatureModel
fitZeroLogNormal - list of parameter estimates for the zero-inflated log normal model
design - model matrix
taxa - taxa names
counts - count matrix
pvalues - calculated p-values
permuttedfits - permutted z-score estimates under the null
data(lungData) lungData = lungData[,-which(is.na(pData(lungData)$SmokingStatus))] lungData=filterData(lungData,present=30,depth=1) lungData <- cumNorm(lungData, p=.5) s <- normFactors(lungData) pd <- pData(lungData) mod <- model.matrix(~1+SmokingStatus, data=pd) lungres1 = fitFeatureModel(lungData,mod)