For each gene in sca, fits the hurdle model in formula
(linear for et>0), logistic for et==0 vs et>0.
Return an object of class ZlmFit
containing slots giving the coefficients, variance-covariance matrices, etc.
After each gene, optionally run the function on the fit named by 'hook'
zlm.SingleCellAssay(formula, sca, method = "bayesglm", silent = TRUE, ebayes = TRUE, ebayesControl = NULL, force = FALSE, hook = NULL, parallel = TRUE, LMlike, onlyCoef = FALSE, ...)
fit
after each gene.option(mc.cores)>1
then multiple cores will be used in fitting.fitArgsC
and fitArgsD
, or coefPrior
.a object of class ZlmFit
with methods to extract coefficients, etc.
The empirical bayes regularization of the gene variance assumes that the precision (1/variance) is drawn from a
gamma distribution with unknown parameters.
These parameters are estimated by considering the distribution of sample variances over all genes.
The procedure used for this is determined from
ebayesControl
, a named list with components 'method' (one of 'MOM' or 'MLE') and 'model' (one of 'H0' or 'H1')
method MOM uses a method-of-moments estimator, while MLE using the marginal likelihood.
H0 model estimates the precisions using the intercept alone in each gene, while H1 fits the full model specified by formula
ebayes, glmlike-class, ZlmFit-class, BayesGLMlike-class
data(vbetaFA) zlmVbeta <- zlm.SingleCellAssay(~ Stim.Condition, subset(vbetaFA, ncells==1)[1:10,])#> #>slotNames(zlmVbeta)#> [1] "coefC" "coefD" "vcovC" #> [4] "vcovD" "LMlike" "sca" #> [7] "deviance" "loglik" "df.null" #> [10] "df.resid" "dispersion" "dispersionNoshrink" #> [13] "priorDOF" "priorVar" "converged" #> [16] "hookOut"#A matrix of coefficients coef(zlmVbeta, 'D')['CCL2',]#> (Intercept) Stim.ConditionUnstim #> -3.8329217 -0.5108005#An array of covariance matrices vcov(zlmVbeta, 'D')[,,'CCL2']#> (Intercept) Stim.ConditionUnstim #> (Intercept) 0.1439196 -0.1254838 #> Stim.ConditionUnstim -0.1254838 0.9182853#> , , metric = lambda #> #> test.type #> primerid cont disc hurdle #> B3GAT1 0.9617702324 0.038250068 1.000020 #> BAX 7.2211565188 3.645901878 10.867058 #> BCL2 0.3766814067 2.202291748 2.578973 #> CCL2 0.8414775522 0.284135226 1.125613 #> CCL3 NA 3.548463195 NA #> CCL4 NA 2.012308210 NA #> CCL5 0.1746468538 0.862093478 1.036740 #> CCR2 5.3383734489 2.308408187 7.646782 #> CCR4 2.0437612666 0.003737042 2.047498 #> CCR5 0.0005534473 2.811952306 2.812506 #> #> , , metric = df #> #> test.type #> primerid cont disc hurdle #> B3GAT1 1 1 2 #> BAX 1 1 2 #> BCL2 1 1 2 #> CCL2 1 1 2 #> CCL3 1 1 2 #> CCL4 1 1 2 #> CCL5 1 1 2 #> CCR2 1 1 2 #> CCR4 1 1 2 #> CCR5 1 1 2 #> #> , , metric = Pr(>Chisq) #> #> test.type #> primerid cont disc hurdle #> B3GAT1 0.326741272 0.84494185 0.606524503 #> BAX 0.007204927 0.05620735 0.004367654 #> BCL2 0.539384664 0.13780573 0.275412150 #> CCL2 0.358974532 0.59400356 0.569608276 #> CCL3 NA 0.05960063 NA #> CCL4 NA 0.15602777 NA #> CCL5 0.676014596 0.35315351 0.595490308 #> CCR2 0.020860929 0.12867576 0.021853574 #> CCR4 0.152831343 0.95125460 0.359245545 #> CCR5 0.981231129 0.09356445 0.245059834 #>