Variance Components
DESCRIPTION:
Returns an object of class "varcomp" which provides estimates of variance components, coefficients and the random variables. The default estimates are MINQUE0, but REML, maximum likelihood and Winsorized estimates are alternatives.
USAGE: 
varcomp(formula, data=<<see below>>, method="minque0",
    start.var=<<see below>>, tol=<<see below>>, nderiv = 1...)
REQUIRED ARGUMENTS:
formula:
a formula describing the model.
OPTIONAL ARGUMENTS:
data:
if given, this should be a data frame. If this is missing, then the objects should exist somewhere on the search list. The data used in formula must contain at least one factor that is random (see is.random).
method:
a character string (or a vector of two strings) giving the estimation technique to be used. The choices are: "minque0", "reml", "ml", and "winsor". If "winsor" is the first string, then one of the other three may be given as a second string. You only need to give enough of the beginning of the string to uniquely identify it.
start.var:
vector of variances at which to start iterative estimation (REML or maximum likelihood). The default is the MINQUE0 estimate.
tol:
tolerance for iterative estimation. The default is .Machine$double.eps^0.25 for REML and .Machine$double.eps^0.5 for ML.
nderiv:
maximum number of derivatives to be used in the optimization (0 <= nderiv <= 2). In general, more derivative information means fewer iterations in the optimization, but more storage and more computation per iteration. All derivatives are calculated within varcomp. The default is 1.
...:
additional arguments for models such as subset and na.action.
VALUE:
an object of class "varcomp", which inherits from "lm". See varcomp.object for details.
DETAILS:
Models of the form y = Xa + Zb + e are fit by this function, where y is a univariate response, X is a known matrix, a is an unknown vector of (fixed) coefficients, Z is a known matrix representing one or more random effects, b is a vector of unknown random variables, and e is a vector of the unknown errors. There needs to be a partition of the columns of Z such that the submatrices, Zi, have only one non-zero entry in each row.

MINQUE0 estimates are computed as in Hartley, Rao, and LaMotte (1978).

Variances components are obtained via unconstrained optimization using the NETLIB functions dmnf, dmng, and dmnh (Gay (1983), Dongarra and Grosse (1987)). The objective function and derivatives are computed using a version of the W-transformation (Hemmerle and Hartley (1973), Corbeil and Searle (1976)) that is numerically stable when the variance component estimates are small in magnitude (Fraley and Burns (1992)).

The Winsorization method produces robust estimates of the variance components. The data are "cleaned" using an initial robust estimate, and then a standard estimation procedure is performed on the cleaned data. See Burns (1992) for more information.

The raov function produces ANOVA estimates for balanced random models.

REFERENCES:
Burns, P. J. (1992). Winsorized REML estimates of variance components. Technical Report, Statistical Sciences, Inc.

Corbeil, R. R. and Searle, S. R. (1976). Restricted maximum likelihood (REML) estimation of variance components in the mixed model. Technometrics 18 31-38.

Dongarra, J. J. and Grosse, E. (1987). Distribution of mathematical software via electronic mail. Communications of the ACM 30 403-407.

Fraley, C. and Burns, P. J. (1995). Large-scale estimation of variance and covariance components. SIAM Journal on Scientific Computing Vol. 16 No. 1.

Gay, D. M. (1983). Algorithm 611. Subroutines for Unconstrained Minimization using a Model/ Trust-Region Approach. ACM Transactions on Mathematical Software, 9 369-383.

Hartley, H. O., Rao, J. N. K., and LaMotte, L. R. (1978). A simple "synthesis"-based mathod of variance component estimation. Biometrics 34 233-242.

Hemmerle, W. J. and Hartley, H. O. (1973). Computing maximum likelihood estimates for the mixed AOV model using the W-transformation. Technometrics 15 819-831.

BUGS:
It should be noted that varcomp will not work for models in which the error variance vanishes. For example, you cannot put a term in the model that is equivalent to the residuals. For example with the pigment data, the formula Moisture ~ Batch/Sample/Test will not work because the variable Test is equivalent to the residuals.
SEE ALSO:
raov , varcomp.object , is.random .
EXAMPLES: 
is.random(pigment) <- T # make all factors random
varcomp(Moisture ~ Batch/Sample, pigment)


varcomp(Moisture ~ Batch/Sample, pigment, method=c("winsor", "r"))




# if subject is a random factor and time.point is numeric, then
# the following model gives a random slope and intercept for each
# subject, plus an overall fixed slope and intercept


varcomp(response ~ time.point * subject)