LMS Quantile Regression with a Box-Cox transformation to a Gamma Distribution
LMS quantile regression with the Box-Cox transformation to the gamma distribution.
lms.bcg(percentiles = c(25, 50, 75), zero = c("lambda", "sigma"),
llambda = "identitylink", lmu = "identitylink", lsigma = "loglink",
idf.mu = 4, idf.sigma = 2, ilambda = 1, isigma = NULL)Given a value of the covariate, this function applies a Box-Cox
transformation to the response to best obtain a gamma distribution.
The parameters
chosen to do this are estimated by maximum likelihood or penalized
maximum likelihood.
Similar details can be found at lms.bcn.
An object of class "vglmff" (see vglmff-class).
The object is used by modelling functions such as vglm,
rrvglm
and vgam.
This VGAM family function comes with the same
warnings as lms.bcn.
Also, the expected value of the second derivative with
respect to lambda may be incorrect (my calculations do
not agree with the Lopatatzidis and Green manuscript.)
Similar notes can be found at lms.bcn.
Thomas W. Yee
Lopatatzidis A. and Green, P. J. (unpublished manuscript). Semiparametric quantile regression using the gamma distribution.
Yee, T. W. (2004). Quantile regression via vector generalized additive models. Statistics in Medicine, 23, 2295–2315.
# This converges, but deplot(fit) and qtplot(fit) do not work
fit0 <- vglm(BMI ~ sm.bs(age, df = 4), lms.bcg, data = bmi.nz, trace = TRUE)
coef(fit0, matrix = TRUE)
## Not run:
par(mfrow = c(1, 1))
plotvgam(fit0, se = TRUE) # Plot mu function (only)
## End(Not run)
# Use a trick: fit0 is used for initial values for fit1.
fit1 <- vgam(BMI ~ s(age, df = c(4, 2)), etastart = predict(fit0),
lms.bcg(zero = 1), bmi.nz, trace = TRUE)
# Difficult to get a model that converges.
# Here, we prematurely stop iterations because it fails near the solution.
fit2 <- vgam(BMI ~ s(age, df = c(4, 2)), maxit = 4,
lms.bcg(zero = 1, ilam = 3), bmi.nz, trace = TRUE)
summary(fit1)
head(predict(fit1))
head(fitted(fit1))
head(bmi.nz)
# Person 1 is near the lower quartile of BMI amongst people his age
head(cdf(fit1))
## Not run:
# Quantile plot
par(bty = "l", mar=c(5, 4, 4, 3) + 0.1, xpd = TRUE)
qtplot(fit1, percentiles=c(5, 50, 90, 99), main = "Quantiles",
xlim = c(15, 90), las = 1, ylab = "BMI", lwd = 2, lcol = 4)
# Density plot
ygrid <- seq(15, 43, len = 100) # BMI ranges
par(mfrow = c(1, 1), lwd = 2)
(aa <- deplot(fit1, x0 = 20, y = ygrid, xlab = "BMI", col = "black",
main = "Density functions at Age = 20 (black), 42 (red) and 55 (blue)"))
aa <- deplot(fit1, x0=42, y=ygrid, add=TRUE, llty=2, col="red")
aa <- deplot(fit1, x0=55, y=ygrid, add=TRUE, llty=4, col="blue", Attach=TRUE)
aa@post$deplot # Contains density function values
## End(Not run)Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.