perks.RdMaximum likelihood estimation of the 2-parameter Perks distribution.
Logical. Suppress a warning? Ignored for VGAM 0.9-7 and higher.
Parameter link functions applied to the
shape parameter shape,
scale parameter scale.
All parameters are treated as positive here
See Links for more choices.
The Perks distribution
has cumulative distribution function
$$F(y; \alpha, \beta) =
1 -
\left\{
\frac{1 + \alpha}{1 + \alpha e^{\beta y}}
\right\}^{1 / \beta}
$$
which leads to a probability density function
$$f(y; \alpha, \beta) =
\left[ 1 + \alpha \right]^{1 / \beta}
\alpha e^{\beta y} / (1 + \alpha e^{\beta y})^{1 + 1 / \beta}
$$
for \(\alpha > 0\),
\(\beta > 0\),
\(y > 0\).
Here, \(\beta\) is called the scale parameter scale,
and \(\alpha\) is called a shape parameter.
The moments for this distribution do
not appear to be available in closed form.
Simulated Fisher scoring is used and multiple responses are handled.
An object of class "vglmff"
(see vglmff-class).
The object is used by modelling functions such as
vglm, and vgam.
Perks, W. (1932). On some experiments in the graduation of mortality statistics. Journal of the Institute of Actuaries, 63, 12–40.
Richards, S. J. (2012). A handbook of parametric survival models for actuarial use. Scandinavian Actuarial Journal. 1–25.
A lot of care is needed because
this is a rather difficult distribution for parameter estimation.
If the self-starting initial values fail then try experimenting
with the initial value arguments, especially iscale.
Successful convergence depends on having very good initial values.
Also, monitor convergence by setting trace = TRUE.
if (FALSE) set.seed(123)
pdata <- data.frame(x2 = runif(nn <- 1000)) # x2 unused
pdata <- transform(pdata, eta1 = -1,
ceta1 = 1)
pdata <- transform(pdata, shape1 = exp(eta1),
scale1 = exp(ceta1))
pdata <- transform(pdata, y1 = rperks(nn, sh = shape1, sc = scale1))
fit1 <- vglm(y1 ~ 1, perks, data = pdata, trace = TRUE)
#> Iteration 1: loglikelihood = -1240.5237
#> Iteration 2: loglikelihood = -1240.4644
#> Iteration 3: loglikelihood = -1240.4614
#> Iteration 4: loglikelihood = -1240.4612
#> Iteration 5: loglikelihood = -1240.4612
#> Iteration 6: loglikelihood = -1240.4612
coef(fit1, matrix = TRUE)
#> loglink(scale) loglink(shape)
#> (Intercept) 1.242033 -1.195504
summary(fit1)
#>
#> Call:
#> vglm(formula = y1 ~ 1, family = perks, data = pdata, trace = TRUE)
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept):1 1.2420 0.1494 8.312 < 2e-16 ***
#> (Intercept):2 -1.1955 0.1972 -6.063 1.34e-09 ***
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Names of linear predictors: loglink(scale), loglink(shape)
#>
#> Log-likelihood: -1240.461 on 1998 degrees of freedom
#>
#> Number of Fisher scoring iterations: 6
#>
#> No Hauck-Donner effect found in any of the estimates
#>
# \dontrun{}