gamlss (version 5.1-4)

prof.term: 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. model=gamlss(y~cs(x,df=this),sigma.fo=~cs(x,df=3),data=abdom), where this indicates the (hyper)parameter to be profiled

criterion

whether global deviance ("GD") or information criterion ("GAIC") is profiled. The default is global deviance criterion="GD"

penalty

The penalty value if information criterion is used in criterion, default penalty=2.5

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. min=1

max

the maximum value for the parameter e.g. max=20

step

how often to evaluate the global deviance (defines the step length of the grid for the parameter) e.g. step=1

length

if the step is left NULL then length is considered for evaluating the grid for the parameter. It has a default value of 11

xlabel

if a label for the axis is required

plot

whether to plot, plot=TRUE the resulting profile deviance (or GAIC)

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

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.

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.

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.

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, http://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 http://www.gamlss.org/).

See Also

gamlss, prof.dev

Examples

# NOT RUN {
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)
# }