Fit a Nonlinear Mixed Effects Model
DESCRIPTION:
A method for the generic function nlme() for objects inheriting from class formula.
USAGE: 
nlme(object, fixed, random = fixed, cluster, data = sys.parent(),
     start, est.method = c("ML", "RML"),  re.block, re.structure =
     c("unstructured", "diagonal", "identity", "compsymm", "ar1"),
     re.paramtr = c("matrixlog", "logcholesky", "cholesky", "spherical",
     "givens"), serial.structure = c("identity", "ar1", "ar1.continuous",
     "compsymm", "ar2", "ma1", "ma2", "arma11"), serial.covariate = NULL,
     serial.covariate.transformation = c("rank.within.cluster", "none",
     "round", "global.rank"), var.function = c("identity", "power",
     "expon", "cte.power"), var.covariate = NULL, var.estimate = T,
     control, na.action, na.pattern, verbose = F)
REQUIRED ARGUMENTS:
object:
a nonlinear model formula with the response on the left of a ~ operator, and an expression involving parameters and covariates on the right. If data is given, all names used in the formula should be defined as parameters or variables in the data frame.
fixed:
a list of formulas of the parameter fixed effects. Usually there will be a fixed effects formula for each parameter. A . on the right hand side of a formula indicates a single fixed effect for that parameter. The formulas are evaluated as linear model formulas in the frame (data).
cluster:
an expression or formula object, specifying the experimental units over which the random effects vary. If cluster is given as a formula, it will have no left side to the ~expression.
start:
a list of initial estimates for the parameters. The fixed component is required. The order of initial estimates in that component corresponds to the order in the fixed argument. Optional components are random, D (scaled variance-covariance matrix of the random effects), theta (the factorized form of the scaled variance-covariance matrix of the random effects), alpha (the serial structure parameters), and delta (the variance function parameters). If present, the random component should be a matrix with as many rows as there are clusters and as many columns as there are random effects. Theta has to be consistent with the chosen structure and parametrization (see re.structure and paramtr below). See serial.structure and var.function below for definitions of alpha and delta respectively.
OPTIONAL ARGUMENTS:
random:
a list of formulas of the parameter random effects. There need not be a random effect formula for each parameter. A . on the right hand side of a formula indicates a single random effect for that parameter. The formula are evaluated as linear model formulas in the frame (data). By default random is equal to fixed.

Note: random effects are always assumed to have mean zero. A nonzero mean can be specified by including an identical term in the fixed effects part of the model.

data:
a data frame in which to interpret the variables named in fixed, random, cluster, serial.covariate, and var.covariate. Default is the calling frame.
est.method:
estimation method; if equal to "ML", Maximum Likelihood is used, otherwise if "RML" is specified, Residual Maximum Likelihood is used. Partial matching of arguments is used, so only the first character needs to be provided. Default is "ML".
re.block:
a list indicating how the random effects should be blocked. Random effects pertaining to different blocks are assumed to be independent. The random effects can be referenced either by their names, or the order in which they appear in the random formula. Elements in the re.block list are vectors containing the names or numbers of the random effects. Within a vector all elements have to be of the same type (i.e. all names or all numbers). By default all random effects are included in the same block.
re.structure:
a character vector describing the variance-covariance structure in each block of the random effects. Predefined structure names are "unstructured" for a general variance-covariance matrix with q(q+1)/2 parameters (q = the number of random effects for the block); "diagonal" for independent random effects with possibly different variances (q parameters); "identity" for independent random effects with the same variance (1 parameter); "compsymm" for a compound symmetry structure, random effects with common variance and common covariance (2 parameters); "ar1" for a common variance and AR(1) autocorrelation structure (2 parameters). Users can define their own variance-covariance structure named as a class with methods for the generic functions lme.re.factor() and lme.re.param(). The generic function lme.re.factor() is used to obtain the factorization of D from the parameter vector theta, and lme.re.param() is used to obtain theta from D (see lme.re.factor.ar1() and lme.re.param.ar1() for examples of the use of the generic functions). If re.structure is of length 1 the name will be used for all blocks. Partial matching against the predefined names is used, so only the first character needs to be provided. Default is "unstructured".
re.paramtr:
the parameterization to be used for the unstructured scaled variance-covariance matrices specified in re.structure; possible values are "cholesky" for the Cholesky decomposition, "logcholesky" for Cholesky using logs of the diagonal elements, "spherical" for spherical coordinates of columns of the Cholesky decomposition, "matrixlog" for the matrix logarithm, and "givens" for a Givens rotation approach to defining eigenvalues and eigenvectors. These variance-covariance parametrizations are described in Pinheiro and Bates (1996) - see the references section. Partial matching of arguments is used, so only the first character needs to be provided. Default is "matrixlog".
serial.structure:
the variance-covariance structure to be used for the within-cluster errors. Predefined structure names are "identity" for independent errors with constant variance; "ar1" for an AR(1) time series structure with respect to an integer valued time covariate (see serial.covariate below); "ar1.continuous" for a continuous time version of an AR(1) structure with the correlation between two observations d units of time apart being given by (a ^ d) where a is a parameter assumed to be one-unit lag correlation; "compsymm" for a compound symmetry structure, within-cluster errors with common variance and common covariance; "ar2" for an AR(2) structure with respect to an integer valued time covariate; "ma1" for an MA(1) structure with respect to an integer valued time covariate; "ma2" for an MA(2) structure with respect to an integer valued time covariate; and "arma11" for an ARMA(1,1) structure with respect to an integer valued time covariate. Users can define their own serial structure named as a class with methods for the generic functions lme.serial.factor(), lme.serial.transform(), and lme.serial.untransform(). The generic function lme.serial.factor() is used to obtain the correlation matrix of the within-cluster errors for a given cluster, lme.serial.transform() is used to transform the original parameter vector into an unrestricted form (used in the optimization algorithms), and lme.serial.untransform() does the opposite operation (see lme.serial.factor.ar1(), lme.serial.transform.ar1(), and lme.serial.untransform.ar1() for examples of such functions). Default is "identity".
serial.covariate:
a formula object specifying the time covariate to be used with serial.structure. This argument has the same characteristics as random. If serial.structure is not equal to "identity" and no serial.covariate is provided, the time covariate will be set to the order of the observations within each cluster.
serial.covariate.transformation:
a transformation to be applied to the time covariate defined in serial.covariate. Possible values are "rank.within.cluster" for the ranks of the time covariate within each cluster; "none" for no transformation; "round" for rounding the time covariate to the nearest integer; and "global.rank" for the overall ranks of the time covariate. Default is "rank.within.cluster".
var.function:
the variance function structure to be used for the within-cluster errors. Predefined structure names are "identity" for constant within-cluster error variance, "power" for a variance function that is proportional to some power of the absolute value of the predictions (var ~ |yhat|^(2a)), "expon" for a variance function that varies exponentially with the absolute value of the predictions (var ~ exp(a*|yhat|)); and "cte.power" for a variance function that varies with some power of the absolute value of the predictions, but approaches a constant value when the predictions go to zero (var ~ a + |yhat|^(2b)). Users can define their own variance functions named as a class with methods for the generic functions lme.varfunc.transform(), lme.varfunc.untransform(), and lme.varfunc.weights(). The generic functions lme.varfunc.transform() and lme.varfunc.untransform() are used to transform the original parameters into an unrestricted form (used in the optimization algorithms) and vice-versa. lme.varfunc.weights() is used to obtain the variance function weights, defined as the inverse of the scaled error standard deviations. This function has the general form proposed by Davidian and Giltinan (1995) - see references below: 1/sqrt(g(yhat, parameters, covariates)), where yhat denotes the predictions, parameters is the parameter vector, and covariates are optional covariates defined in var.covariate (see below). See lme.varfunc.transform.power(), lme.varfunc.untransform.power(), and lme.varfunc.weights.power() for examples of such functions). Default is "identity".
var.covariate:
a formula object specifying the covariates to be used in the calculation of the variance function weights. These covariates will be passed to the generic function lme.varfunc.weights() as the model.frame corresponding to var.covariate, using the data frame data. This argument has the same characteristics as random.
var.estimate:
a logical variable indicating whether the parameters in the variance function structure should be estimated. If equal to TRUE the variance function parameters will be estimated together with the other parameters in the model; if FALSE the variance function parameters will be kept fixed at their initial values (which may be given in the start list). Default is TRUE.
control:
a list of control values for the nonlinear estimation algorithm.
subset:
expression saying which subset of the rows of the data should be used in the fit. This can be a logical vector, or a numeric vector indicating which observation numbers are to be included, or a character vector of the row names to be included. All observations are included by default.
na.action:
a function to filter missing data. This is applied to the model.frame after any subset argument has been used. The default (with na.fail) is to create an error if any missing values are found. A possible alternative is na.omit, which deletes observations that contain one or more missing values.
na.pattern:
an expression or formula object, specifying which returned values are to be regarded as missing.
verbose:
if TRUE, information on the evolution of the iterative algorithm (number of iterations, convergence criterion, etc.) will be printed. Default is FALSE.
VALUE:
an object of class nlme representing the fit. Generic functions such as print(), plot(), predict(), and summary() have methods to show the results of the fit. See nlme.object for the components of the fit. The functions residuals(), coefficients() and fitted.values() can be used to extract some of its components.
DETAILS:
The computational methods are described in Lindstrom, M.J. and Bates, D.M. (1990). The variance-covariance parametrizations are described in Pinheiro, J.C. and Bates., D.M. (1996). The different time series structures available for the serial.structure argument are described in Box, G.E.P., Jenkins, G.M., and Reinsel G.C. (1994). The use of variance functions for linear and nonlinear mixed effects models is presented in detail in Davidian, M. and Giltinan, D.M. (1995).
REFERENCES:
Lindstrom, M.J. and Bates, D.M. (1990). Nonlinear Mixed Effects Models for Repeated Measures Data. Biometrics 46, 673-687.

Pinheiro, J.C. and Bates., D.M. (1996). Unconstrained Parametrizations for Variance-Covariance Matrices. Statistics and Computing, to appear.

Box, G.E.P., Jenkins, G.M., and Reinsel G.C. (1994). Time Series Analysis: Forecasting and Control, 3rd Edition. Holden-Day.

Davidian, M. and Giltinan, D.M. (1995). Nonlinear Mixed Effects Models for Repeated Measurement Data. Chapman and Hall.

SEE ALSO:
nlme.nlsList , nlme.control , nlme.object .
EXAMPLES: 
# Example from Lindstrom and Bates (1990) Biometrics
Orange.fit <- nlme(circumference ~ A/(1 + exp((B - age)/C)),
                   fixed = list(A ~ ., B ~ ., C ~ .),
                   random = list(A ~ ., B ~ .), cluster = ~ Tree,
                   data = Orange, start = list(fixed = c(160, 700, 350)),
                   serial.structure = "ar2", var.function = "expon")
Orange.fit


# Returns the following:
Call:
  Model: circumference ~ A/(1 + exp((B - age)/C))
  Fixed: list(A ~ ., B ~ ., C ~ .)
 Random: list(A ~ ., B ~ .)
Cluster:  ~ Tree
   Data: Orange


Variance/Covariance Components Estimate(s):


  Structure: matrixlog
  Standard Deviation(s) of Random Effect(s)
        A        B
 35.94847 68.57305
 Correlation of Random Effects
          A
B 0.6613264


 Cluster Residual Variance: 27.79603


 Serial Correlation Structure: ar2
 Serial Correlation Parameter(s): -0.42132853  0.05348223


 Variance Function: expon
 Variance Function Parameter(s): 0.0034951


Fixed Effects Estimate(s):
         A        B        C
 189.5978 705.7406 342.2209


Number of Observations: 35
Number of Clusters: 5