betaprime.RdEstimation of the two shape parameters of the beta-prime distribution by maximum likelihood estimation.
betaprime(lshape = "loglink", ishape1 = 2, ishape2 = NULL,
zero = NULL)Parameter link function applied to the two (positive) shape
parameters. See Links for more choices.
See CommonVGAMffArguments for more information.
The beta-prime distribution is given by $$f(y) = y^{shape1-1} (1+y)^{-shape1-shape2} / B(shape1,shape2)$$ for \(y > 0\). The shape parameters are positive, and here, \(B\) is the beta function. The mean of \(Y\) is \(shape1 / (shape2-1)\) provided \(shape2>1\); these are returned as the fitted values.
If \(Y\) has a \(Beta(shape1,shape2)\) distribution then \(Y/(1-Y)\) and \((1-Y)/Y\) have a \(Betaprime(shape1,shape2)\) and \(Betaprime(shape2,shape1)\) distribution respectively. Also, if \(Y_1\) has a \(gamma(shape1)\) distribution and \(Y_2\) has a \(gamma(shape2)\) distribution then \(Y_1/Y_2\) has a \(Betaprime(shape1,shape2)\) distribution.
An object of class "vglmff"
(see vglmff-class).
The object is used by modelling functions
such as vglm,
rrvglm
and vgam.
Johnson, N. L. and Kotz, S. and Balakrishnan, N. (1995). Chapter 25 of: Continuous Univariate Distributions, 2nd edition, Volume 2, New York: Wiley.
The response must have positive values only.
The beta-prime distribution is also known as the beta distribution of the second kind or the inverted beta distribution.
nn <- 1000
bdata <- data.frame(shape1 = exp(1), shape2 = exp(3))
bdata <- transform(bdata, yb = rbeta(nn, shape1, shape2))
bdata <- transform(bdata, y1 = (1-yb) / yb,
y2 = yb / (1-yb),
y3 = rgamma(nn, exp(3)) / rgamma(nn, exp(2)))
fit1 <- vglm(y1 ~ 1, betaprime, data = bdata, trace = TRUE)
#> Iteration 1: loglikelihood = -4741.6257
#> Iteration 2: loglikelihood = -4135.7461
#> Iteration 3: loglikelihood = -3673.4055
#> Iteration 4: loglikelihood = -3363.8237
#> Iteration 5: loglikelihood = -3225.2613
#> Iteration 6: loglikelihood = -3201.0174
#> Iteration 7: loglikelihood = -3200.3655
#> Iteration 8: loglikelihood = -3200.3651
#> Iteration 9: loglikelihood = -3200.3651
coef(fit1, matrix = TRUE)
#> loglink(shape1) loglink(shape2)
#> (Intercept) 3.025668 1.022632
fit2 <- vglm(y2 ~ 1, betaprime, data = bdata, trace = TRUE)
#> Iteration 1: loglikelihood = -8157.1713
#> Iteration 2: loglikelihood = -7157.5035
#> Iteration 3: loglikelihood = -6158.4056
#> Iteration 4: loglikelihood = -5160.8488
#> Iteration 5: loglikelihood = -4167.4262
#> Iteration 6: loglikelihood = -3184.8533
#> Iteration 7: loglikelihood = -2229.2428
#> Iteration 8: loglikelihood = -1333.0333
#> Iteration 9: loglikelihood = -542.51781
#> Iteration 10: loglikelihood = 105.37514
#> Iteration 11: loglikelihood = 600.01698
#> Iteration 12: loglikelihood = 936.28849
#> Iteration 13: loglikelihood = 1101.0387
#> Iteration 14: loglikelihood = 1136.3296
#> Iteration 15: loglikelihood = 1137.7288
#> Iteration 16: loglikelihood = 1137.7309
#> Iteration 17: loglikelihood = 1137.7309
coef(fit2, matrix = TRUE)
#> loglink(shape1) loglink(shape2)
#> (Intercept) 1.022632 3.025668
fit3 <- vglm(y3 ~ 1, betaprime, data = bdata, trace = TRUE)
#> Iteration 1: loglikelihood = -2204.3036
#> Iteration 2: loglikelihood = -1870.5549
#> Iteration 3: loglikelihood = -1687.7341
#> Iteration 4: loglikelihood = -1642.7626
#> Iteration 5: loglikelihood = -1640.5168
#> Iteration 6: loglikelihood = -1640.5116
#> Iteration 7: loglikelihood = -1640.5116
coef(fit3, matrix = TRUE)
#> loglink(shape1) loglink(shape2)
#> (Intercept) 2.98223 1.991868
# Compare the fitted values
with(bdata, mean(y3))
#> [1] 3.10824
head(fitted(fit3))
#> [,1]
#> [1,] 3.117571
#> [2,] 3.117571
#> [3,] 3.117571
#> [4,] 3.117571
#> [5,] 3.117571
#> [6,] 3.117571
Coef(fit3) # Useful for intercept-only models
#> shape1 shape2
#> 19.731765 7.329211