prof.term {gamlss} | R Documentation |
Plotting the Profile: deviance or information criterion for one of the terms (or hyper-parameters) in a GAMLSS model
Description
This functions plots the profile deviance for a chosen parameter included in the linear predictor of any of the mu
,sigma
, nu
or tau
models so profile confidence intervals can be obtained. In can also be used to plot the profile of a specified information criterion for any hyper-parameter when smooth additive terms are used.
Usage
prof.term(model = NULL, criterion = c("GD", "GAIC"), penalty = 2.5,
other = NULL, min = NULL, max = NULL, step = NULL,
length = 7, xlabel = NULL, plot = TRUE, perc = 95,
start.prev = TRUE, col="darkgreen")
Arguments
model |
this is a GAMLSS model, e.g. |
criterion |
whether global deviance ("GD") or information criterion ("GAIC") is profiled. The default is global deviance |
penalty |
The penalty value if information criterion is used in |
other |
this can be used to evaluate an expression before the actual fitting of the model (Make sure that those expressions are well define in the global environment) |
min |
the minimum value for the parameter e.g. |
max |
the maximum value for the parameter e.g. |
step |
how often to evaluate the global deviance (defines the step length of the grid for the parameter) e.g. |
length |
if the step is left NULL then |
xlabel |
if a label for the axis is required |
plot |
whether to plot, |
perc |
what % confidence interval is required |
start.prev |
whether to start from the previous fitted model parameters values or not (default is TRUE) |
col |
the color of the profile line |
Details
This function can be use to provide likelihood based confidence intervals for a parameter involved in terms in the linear predictor(s).
These confidence intervals are more accurate than the ones obtained from the parameters' standard errors.
The function can also be used to plot a profile information criterion (with a given penalty) against a hyper-parameter. This can be used to check the uniqueness in hyper-parameter determination using for example find.df
.
Value
Return a profile plot (if the argument plot=TRUE
) and an ProfLikelihood.gamlss
object if saved. The object contains:
values |
the values at the grid where the parameter was evaluated |
fun |
the function which approximates the points using splines |
min |
the minimum values in the grid |
max |
the maximum values in the grid |
max.value |
the value of the parameter maximising the Profile deviance (or GAIC) |
CI |
the profile confidence interval (if global deviance is used) |
criterion |
which criterion was used |
Warning
A dense grid (i.e. small step) evaluation of the global deviance can take a long time, so start with a sparse grid (i.e. large step) and decrease gradually the step length for more accuracy.
Author(s)
Mikis Stasinopoulos and Bob Rigby
References
Rigby, R. A. and Stasinopoulos D. M. (2005). Generalized additive models for location, scale and shape,(with discussion), Appl. Statist., 54, part 3, pp 507-554.
Rigby, R. A., Stasinopoulos, D. M., Heller, G. Z., and De Bastiani, F. (2019) Distributions for modeling location, scale, and shape: Using GAMLSS in R, Chapman and Hall/CRC. An older version can be found in https://www.gamlss.com/.
Stasinopoulos D. M. Rigby R.A. (2007) Generalized additive models for location scale and shape (GAMLSS) in R. Journal of Statistical Software, Vol. 23, Issue 7, Dec 2007, https://www.jstatsoft.org/v23/i07/.
Stasinopoulos D. M., Rigby R.A., Heller G., Voudouris V., and De Bastiani F., (2017) Flexible Regression and Smoothing: Using GAMLSS in R, Chapman and Hall/CRC.
(see also https://www.gamlss.com/).
See Also
Examples
data(aids)
# fitting a linear model
gamlss(y~x+qrt,family=NBI,data=aids)
# testing the linear beta parameter
mod<-quote(gamlss(y ~ offset(this * x) + qrt, data = aids, family = NBI))
prof.term(mod, min=0.06, max=0.11)
# find the hyper parameter using cubic splines smoothing
mod1<-quote(gamlss(y ~ cs(x,df=this) + qrt, data = aids, family = NBI))
prof.term(mod1, min=1, max=15, step=1, criterion="GAIC", penalty=log(45))
# find a break point in x
mod2 <- quote(gamlss(y ~ x+I((x>this)*(x-this))+qrt,family=NBI,data=aids))
prof.term(mod2, min=1, max=45, step=1, criterion="GD")
rm(mod,mod1,mod2)