The function HYPERPO2()
defines the
hyper Poisson distribution (with mu as mean)
a two parameter distribution,
for a gamlss.family
object to be used in GAMLSS
fitting using the function gamlss()
.
HYPERPO2(mu.link = "log", sigma.link = "log")
Returns a gamlss.family
object which can be used
to fit a hyper-Poisson distribution version 2
in the gamlss()
function.
The hyper-Poisson distribution with parameters \(\mu\) and \(\sigma\) has a support 0, 1, 2, ...
Note: in this implementation the parameter \(\mu\) is the mean of the distribution and \(\sigma\) corresponds to the dispersion parameter. If you fit a model with this parameterization, the time will increase because an internal procedure to convert \(\mu\) to \(\lambda\) parameter.
Sáez-Castillo, A. J., & Conde-Sánchez, A. (2013). A hyper-Poisson regression model for overdispersed and underdispersed count data. Computational Statistics & Data Analysis, 61, 148-157.
# Example 1
# Generating some random values with
# known mu and sigma
set.seed(1234)
y <- rHYPERPO2(n=200, mu=4, sigma=1.5)
# Fitting the model
library(gamlss)
mod1 <- gamlss(y~1, sigma.fo=~1, family=HYPERPO2,
control=gamlss.control(n.cyc=500, trace=TRUE))
#> GAMLSS-RS iteration 1: Global Deviance = 858.2628
#> GAMLSS-RS iteration 2: Global Deviance = 858.2627
# Extracting the fitted values for mu and sigma
# using the inverse link function
exp(coef(mod1, what="mu"))
#> (Intercept)
#> 3.954995
exp(coef(mod1, what="sigma"))
#> (Intercept)
#> 1.52176
# Example 2
# Generating random values under some model
# \donttest{
# A function to simulate a data set with Y ~ HYPERPO2
gendat <- function(n) {
x1 <- runif(n)
x2 <- runif(n)
mu <- exp(1.21 - 3 * x1) # 0.75 approximately
sigma <- exp(1.26 - 2 * x2) # 1.30 approximately
y <- rHYPERPO2(n=n, mu=mu, sigma=sigma)
data.frame(y=y, x1=x1, x2=x2)
}
set.seed(12345)
dat <- gendat(n=200)
mod2 <- gamlss(y~x1, sigma.fo=~x2, family=HYPERPO2, data=dat,
control=gamlss.control(n.cyc=500, trace=TRUE))
#> GAMLSS-RS iteration 1: Global Deviance = 434.3557
#> GAMLSS-RS iteration 2: Global Deviance = 434.3483
#> GAMLSS-RS iteration 3: Global Deviance = 434.3483
summary(mod2)
#> Warning: summary: vcov has failed, option qr is used instead
#> ******************************************************************
#> Family: c("HYPERPO2", "Hyper-Poisson-2")
#>
#> Call: gamlss(formula = y ~ x1, sigma.formula = ~x2, family = HYPERPO2,
#> data = dat, control = gamlss.control(n.cyc = 500, trace = TRUE))
#>
#> Fitting method: RS()
#>
#> ------------------------------------------------------------------
#> Mu link function: log
#> Mu Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 1.1581 0.1121 10.33 <2e-16 ***
#> x1 -2.9576 0.2941 -10.06 <2e-16 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> Sigma link function: log
#> Sigma Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 1.5962 0.8132 1.963 0.0511 .
#> x2 -2.9907 1.5643 -1.912 0.0573 .
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> No. of observations in the fit: 200
#> Degrees of Freedom for the fit: 4
#> Residual Deg. of Freedom: 196
#> at cycle: 3
#>
#> Global Deviance: 434.3483
#> AIC: 442.3483
#> SBC: 455.5416
#> ******************************************************************
# }