lqmm (version 1.0)

lqmm: Fitting Linear Quantile Mixed Models

Description

lqmm is used to fit linear quantile mixed models based on the asymmetric Laplace distribution.

Usage

lqmm(fixed, random, group, covariance = "pdDiag", iota = 0.5,
	nK = 11, type = "normal", data = sys.frame(sys.parent()),
	subset, weights, na.action = na.fail, control = list(),
	contrasts = NULL, forceK = FALSE, fit = TRUE)

Arguments

fixed

an object of class formula for fixed effects: a symbolic description of the model to be fitted.

random

a one-sided formula of the form ~x1 + x2 + ... + xn for random effects: a symbolic description of the model to be fitted.

group

grouping factor.

covariance

variance--covariance matrix of the random effects. Default is pdDiag (see details).

iota

the quantile(s) to be estimated.

nK

number of quadrature knots.

type

type of quadrature "c("normal","robust")" (see details).

data

an optional data frame containing the variables named in fixed, random and group. By default the variables are taken from the environment from which lqmm is called.

subset

an optional vector specifying a subset of observations to be used in the fitting process.

weights

an optional vector of weights to be used in the fitting process of the same length as the number of rows of data.

na.action

a function that indicates what should happen when the data contain NAs. The default action (na.fail) causes lqmm to print an error message and terminate if there are any incomplete observations.

control

list of control parameters of the fitting process. See lqmmControl.

contrasts

not yet implemented.

forceK

logical flag. If TRUE the number of quadrature knots if forced to be equal to nK. This option is used to bypass a control check that otherwise automatically sets nK to 7 in case nK > 7 and the number of

fit

logical flag. If FALSE the function returns a list of arguments to be passed to lqmm.fit.

Value

• lqmm returns an object of class lqmm. The function summary is used to obtain and print a summary of the results. The generic accessor functions coefficients, predict and residuals extract various useful features of the value returned by lqmm. An object of class lqmm is a list containing the following components:
• theta
  a vector containing fixed regression coefficients and parameters of the variance--covariance matrix of the random effects. See cov.lqmm to extract the variance--covariance of the random effects from an "lqmm" object.
• theta_x,theta_z
  partition of theta: fixed regression coefficients (theta_x) and unique variance--covariance parameters (theta_z).
• scale
  the scale parameter.
• logLik
  the log--likelihood.
• opt
  details on optimization (see lqmm.fit.gs and lqmm.fit.df).
• call
  the matched call.
• nn
  column names of mmf.
• mm
  column names of mmr.
• nobs
  the number of observations.
• dim_theta
  the number of columns in mmf and mmr.
• dim_theta_z
  the length of theta_z.
• edf
  length of theta.
• rdf
  the number of residual degrees of freedom.
• df
  edf + 1 (scale parameter).
• iota
  the estimated quantile(s).
• mmf
  the model matrix -- fixed effects.
• mmr
  the model matrix -- random effects.
• y
  the model response.
• revOrder
  original order of observations (now ordered according to group).
• weights
  the likelihood weights used in the fitting process (a vector of 1's if weights is missing or NULL).
• group
  the grouping factor.
• ngroups
  the number of groups.
• QUAD
  the knots and weights of the quadrature.
• type
  the type of quadrature.
• InitialPar
  starting values for theta.
• control
  list of control parameters used for optimization (see lqmmControl).
• cov_name
  class of variance-covariance matrix for the random effects.
• mfArgs
  arguments for model.frame to return the full data frame.

Details

The function computes an estimate on the iota-th quantile function of the response, conditional on the covariates, as specified by the formula argument, and on random effects, as specified by the random argument. The quantile predictor is assumed to be linear. The function maximizes the (log)likelihood of the Laplace regression proposed by Geraci and Bottai (2007). The likelihood is numerically integrated via Gaussian quadrature techniques. The optimization algorithm is based on a Nelder-Mead algorithm (control = list(method = "df")) via optim. An alternative optimization algorithm is based on the gradient of the Laplace log--likelihood (control = list(method = "gs")) (Bottai, Orsini and Geraci, 2011). The scale parameter is optimized in a refinement step via optimize. Quadrature approaches include Gauss-Hermite ("normal") and Gauss-Laguerre ("robust") quadrature. Different standard classes of positive--definite matrices for the random effects can be specified: pdSymm general positive--definite matrix, with no additional structure; pdDiag diagonal; pdIdent multiple of an identity; pdCompSymm compound symmetry structure (constant diagonal and constant off--diagonal elements).

References

Geraci M and Bottai M (2007). Quantile regression for longitudinal data using the asymmetric Laplace distribution. Biostatistics 8(1), 140--154. Bottai M, Orsini N and Geraci M (2011). A Gradient Search Maximization Algorithm for Laplace Likelihood. Unpublished manuscript. Geraci M and Bottai M (1 June 2011). Linear Quantile Mixed Models. Unpublished manuscript.

See Also

lqm, summary.lqmm

Examples

# Test example
set.seed(123)

M <- 50
n <- 10
test <- data.frame(x = runif(n*M,0,1), group = rep(1:M,each=n))
test$y <- 10*test$x + rep(rnorm(M, 0, 2), each = n) + rchisq(n*M, 3)
fit.lqmm <- lqmm(fixed = y ~ x, random = ~ 1, group = group,	
	data = test, iota = 0.5, nK = 11, type = "normal")
fit.lqmm

#Call: lqmm(fixed = y ~ x, random = ~1, group = group, iota = 0.5, nK = 11, 
#    type = "normal", data = test)
#Quantile 0.5 
#Fixed effects:
#Intercept          x  
#    3.529      9.201  
#Random effects:
#Intercept 
#    3.312 
#Residual scale parameter: 0.8692 (standard deviation 2.459)
#Log-likelihood: -1178 
#Number of observations: 500 
#Number of groups: 50 


## Orthodont data
data(Orthodont)

# Random intercept model
fitOi.lqmm <- lqmm(distance ~ age, random = ~ 1, group = Subject,
	iota = c(0.1,0.5,0.9), data = Orthodont)
coef(fitOi.lqmm)

# Random slope model
fitOs.lqmm <- lqmm(distance ~ age, random = ~ age, group = Subject,
	iota = c(0.1,0.5,0.9), cov = "pdDiag", data = Orthodont)

# Extract estimates
cov.lqmm(fitOs.lqmm)
coef(fitOs.lqmm)
raneff.lqmm(fitOs.lqmm)

# AIC
AIC(fitOi.lqmm)
AIC(fitOs.lqmm)