The HPMIXED Procedure

Model Assumptions

The following sections provide an overview of the approach used by the HPMIXED procedure for likelihood-based analysis of linear mixed models with sparse matrix technique. Additional theory and examples are provided in Littell et al. (1996); Verbeke and Molenberghs (1997, 2000); Brown and Prescott (1999).

The HPMIXED procedure fits models generally of the form

$\mb {y} = \bX \bbeta + \bZ \bgamma +\bepsilon$

Models of this form contain both fixed-effects parameters, $\bbeta$ , and random-effects parameters, $\bgamma$ ; hence, they are called mixed models. See Henderson (1990) and Searle, Casella, and McCulloch (1992) for historical developments of the mixed model. Note that the matrix $\bZ$ can contain either continuous or dummy variables, just like $\bX$ .

So far this is the same general form of model fit by the MIXED procedure. The difference between the models handled by the two procedures lies in the assumptions about the distributions of $\bgamma$ and $\bepsilon$ . For both procedures a key assumption is that $\bgamma$ and $\bepsilon$ are normally distributed with

$\begin{eqnarray*} \mbox{E}\left[ \begin{array}{c} \bgamma \\ \bepsilon \end{array} \right] & = & \left[\begin{array}{c} \Strong{0} \\ \Strong{0} \end{array} \right] \\ \mbox{Var}\left[ \begin{array}{c} \bgamma \\ \bepsilon \end{array} \right] & = & \left[\begin{array}{cc} \bG & \Strong{0} \\ \Strong{0} & \bR \end{array} \right] \end{eqnarray*}$

The two procedures differ in their assumptions about the variance matrices $\bG$ and $\bR$ for $\bgamma$ and $\bepsilon$ , respectively. The MIXED procedure allows a variety of different structures for both $\bG$ and $\bR$ ; while in HPMIXED procedure, $\bR$ is always assumed to be of the form $\bR = \bI \sigma ^2$ , and the structures available for modeling $\bG$ are only a small subset of the structures offered by the MIXED procedure.

Estimates of fixed effects and predictions for random effects are obtained by solving the so-called mixed model equations:

$\left[ \begin{array}{cc} \bX ’\bX /\sigma ^2 & \bX ’\bZ /\sigma ^2 \\ \bZ ’\bX /\sigma ^2 & \bZ ’\bZ /\sigma ^2 + \bG ^{-1} \end{array}\right] \left[\begin{array}{c} \widehat{\bbeta } \\ \widehat{\bgamma } \end{array}\right] = \left[\begin{array}{c} \bX ’\mb {y}/\sigma ^2 \\ \bZ ’\mb {y}/\sigma ^2 \end{array}\right]$

Let $\bC$ denote the coefficient matrix of the mixed model equations:

$\bC = \left[ \begin{array}{cc} \bX ’\bX /\sigma ^2 & \bX ’\bZ /\sigma ^2 \\ \bZ ’\bX /\sigma ^2 & \bZ ’\bZ /\sigma ^2 + \bG ^{-1} \end{array}\right]$

Under the assumptions given previously for the moments of $\bgamma$ and $\bepsilon$ , the variance of $\mb {y}$ is $\bV = \bZ \bG \bZ ’ + \bI \sigma ^2$ . You can model $\bV$ by setting up the random-effects design matrix $\bZ$ and by specifying covariance structures for $\bG$ . Let $\btheta$ be a vector of all unknown parameters in $\bG$ . Then the general form of the restricted likelihood function for the mixed models that the HPMIXED procedure can fit is

$\begin{equation*} L(\btheta ,\sigma ^2) = -2\log l = (n-p) \log (2 \pi ) + \log |\bC | + \log |\bG | + n \log (\sigma ^2) + \mb {y}’\bP \mb {y} \end{equation*}$

where

$\bP = \bV ^{-1} - \bV ^{-1}\bX (\bX ’\bV ^{-1}\bX )^{-}\bX ’\bV ^{-1}$

and p is the rank of $\bX$ . The HPMIXED procedure minimizes $L(\btheta ,\sigma ^2)$ over all unknown parameters in $\btheta$ and $\sigma ^2$ by using nonlinear optimization algorithms.