Fits the Bell regression model to overdispersed count data.
Arguments
- formula
an object of class "formula" (or one that can be coerced to that class): a symbolic description of the model to be fitted.
- data
an optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the variables in the model. If not found in data, the variables are taken from environment(formula), typically the environment from which ypbp is called.
- approach
approach to be used to fit the model (mle: maximum likelihood; bayes: Bayesian approach).
- hessian
hessian logical; If TRUE (default), the hessian matrix is returned when approach="mle".
- link1
assumed link function for degenerate distribution (logit, probit, cloglog, cauchy); default is logit.
- link2
assumed link function for count distribution (log, sqrt or identiy); default is log.
- hyperpars
a list containing the hyperparameters associated with the prior distribution of the regression coefficients; if not specified then default choice is hyperpars = c(mu_psi = 0, sigma_psi = 10, mu_beta = 0, sigma_beta = 10).
- ...
further arguments passed to either
rstan::optimizing
orrstan::sampling
.
Examples
# \donttest{
# ML approach:
data(cells)
mle <- zibellreg(cells ~ smoker+gender|smoker+gender, data = cells, approach = "mle")
summary(mle)
#> Call:
#> zibellreg(formula = cells ~ smoker + gender | smoker + gender,
#> data = cells, approach = "mle")
#>
#> Zero-inflated regression coefficients:
#> Estimate StdErr z.value p.value
#> (Intercept) -1.95159 0.84476 -2.3102 0.020876 *
#> smoker 2.17588 0.82300 2.6439 0.008197 **
#> gender -0.49634 0.42068 -1.1798 0.238062
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#>
#> Count regression coefficients:
#> Estimate StdErr z.value p.value
#> (Intercept) 0.716811 0.179865 3.9853 6.74e-05 ***
#> smoker -0.611907 0.183417 -3.3361 0.0008495 ***
#> gender 0.036084 0.177493 0.2033 0.8389021
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> ---
#> logLik = -610.3234 AIC = 1232.647
# Bayesian approach:
bayes <- zibellreg(cells ~ 1|smoker+gender, data = cells, approach = "bayes", refresh = FALSE)
summary(bayes)
#> Call:
#> zibellreg(formula = cells ~ 1 | smoker + gender, data = cells,
#> approach = "bayes", refresh = FALSE)
#>
#> Zero-inflated regression coefficients:
#> mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff
#> (Intercept) -1.208 0.051 0.694 -1.997 -1.349 -1.122 -0.929 -0.641 187.654
#> Rhat
#> (Intercept) 1.019
#>
#> Count regression coefficients:
#> mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff
#> (Intercept) 0.717 0.003 0.150 0.424 0.616 0.716 0.815 1.012 2692.681
#> smoker -1.077 0.003 0.144 -1.356 -1.174 -1.075 -0.979 -0.802 2615.948
#> gender 0.173 0.003 0.144 -0.105 0.073 0.172 0.270 0.458 2754.854
#> Rhat
#> (Intercept) 1.000
#> smoker 1.001
#> gender 1.001
#> ---
#> Inference for Stan model: zibellreg.
#> 4 chains, each with iter=2000; warmup=1000; thin=1;
#> post-warmup draws per chain=1000, total post-warmup draws=4000.
#>
# }