inv.binomial.RdEstimates the two parameters of an inverse binomial distribution by maximum likelihood estimation.
inv.binomial(lrho = "extlogitlink(min = 0.5, max = 1)",
llambda = "loglink", irho = NULL, ilambda = NULL, zero = NULL)Link function for the \(\rho\) and \(\lambda\)
parameters.
See Links for more choices.
Numeric. Optional initial values for \(\rho\) and \(\lambda\).
The inverse binomial distribution of Yanagimoto (1989) has density function $$f(y;\rho,\lambda) = \frac{ \lambda \,\Gamma(2y+\lambda) }{\Gamma(y+1) \, \Gamma(y+\lambda+1) } \{ \rho(1-\rho) \}^y \rho^{\lambda}$$ where \(y=0,1,2,\ldots\) and \(\frac12 < \rho < 1\), and \(\lambda > 0\). The first two moments exist for \(\rho>\frac12\); then the mean is \(\lambda (1-\rho) /(2 \rho-1)\) (returned as the fitted values) and the variance is \(\lambda \rho (1-\rho) /(2 \rho-1)^3\). The inverse binomial distribution is a special case of the generalized negative binomial distribution of Jain and Consul (1971). It holds that \(Var(Y) > E(Y)\) so that the inverse binomial distribution is overdispersed compared with the Poisson distribution.
An object of class "vglmff"
(see vglmff-class).
The object is used by modelling functions
such as vglm
and vgam.
Yanagimoto, T. (1989). The inverse binomial distribution as a statistical model. Communications in Statistics: Theory and Methods, 18, 3625–3633.
Jain, G. C. and Consul, P. C. (1971). A generalized negative binomial distribution. SIAM Journal on Applied Mathematics, 21, 501–513.
Jorgensen, B. (1997). The Theory of Dispersion Models. London: Chapman & Hall
This VGAM family function only works reasonably well with
intercept-only models.
Good initial values are needed; if convergence failure occurs
use irho and/or ilambda.
Some elements of the working weight matrices use the expected information matrix while other elements use the observed information matrix. Yet to do: using the mean and the reciprocal of \(\lambda\) results in an EIM that is diagonal.
idata <- data.frame(y = rnbinom(n <- 1000, mu = exp(3), size = exp(1)))
fit <- vglm(y ~ 1, inv.binomial, data = idata, trace = TRUE)
#> Iteration 1: loglikelihood = -246318.9
#> Iteration 2: loglikelihood = -245324.08
#> Iteration 3: loglikelihood = -244341.44
#> Iteration 4: loglikelihood = -243391.08
#> Iteration 5: loglikelihood = -242522.41
#> Iteration 6: loglikelihood = -241837.76
#> Iteration 7: loglikelihood = -241464.45
#> Iteration 8: loglikelihood = -241375.16
#> Iteration 9: loglikelihood = -241370.96
#> Iteration 10: loglikelihood = -241370.95
#> Iteration 11: loglikelihood = -241370.95
with(idata, c(mean(y), head(fitted(fit), 1)))
#> [1] 20.291 20.291
summary(fit)
#>
#> Call:
#> vglm(formula = y ~ 1, family = inv.binomial, data = idata, trace = TRUE)
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept):1 -5.64687 1.41462 -3.992 6.56e-05 ***
#> (Intercept):2 -1.94354 0.03156 -61.583 < 2e-16 ***
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Names of linear predictors: extlogitlink(rho, min = 0.5, max = 1),
#> loglink(lambda)
#>
#> Log-likelihood: -241371 on 1998 degrees of freedom
#>
#> Number of Fisher scoring iterations: 11
#>
#> Warning: Hauck-Donner effect detected in the following estimate(s):
#> '(Intercept):1'
#>
coef(fit, matrix = TRUE)
#> extlogitlink(rho, min = 0.5, max = 1) loglink(lambda)
#> (Intercept) -5.646867 -1.943542
Coef(fit)
#> rho lambda
#> 0.5017581 0.1431958
sum(weights(fit)) # Sum of the prior weights
#> [1] 1000
sum(weights(fit, type = "work")) # Sum of the working weights
#> [1] 1004.513