A Multilevel Structural Equation Model for Dyadic Data

Jason T. Newsom

Portland State Unversity

I begin by giving a brief overview of latent growth models and multilevel regression (i.e., hierarchical linear models). I've assumed some familiarity with both of these techniques, but I've summarized them to illustrate their parallels in the case of repeated measures. I then proceed by proposing a (new, I think) structural equation model approach to hierarchical regression in the case of dyadic data. All comments and questions are welcome and encouraged ().

Latent Growth Curve Models

Structural equation modeling (SEM) can be used to estimate individual growth curves by using repeated measures as indicators of two latent variables, an intercept variable (0) and a slope variable (1), called “latent growth curve models”. The interpretation of the intercept variable depends on how the loadings on the slope factor are fixed. For instance, one approach to defining the slope variable is to fix loadings to values 0, 1, 2, 3,4, . . . t-1. In this case, the intercept latent variable represents the initial value, because the first loading on 1 is set to 0. One can also “center “ these loadings by setting the middle time point to 0 (e.g., -3.,-2,-1,0,1,2,3), giving the intercept factor the value of the average score across all time points [see here for a Figure showing the general model specification (with estimated loadings on the slope factor, discussed later)].

Mathematically, the latent growth curve model is represented by the following set of formulas:

level-1 equation (measurement model):

/ (1.1)

level-2 equations (structural model):

/ (1.2)
/ (1.3)

yit is the dependent variable. The subscripts i and t indicate a measurement within an individual, i, for each time point, t. 0 is a latent variable that represents the level-1 intercept, 1 is a latent variable that represents the relationship between the time code and the dependent variable (i.e., the growth trajectory), it are the loadings for each time point on the intercept latent variable (0) and the slope latent variable (1). (The intercept, , associated with each loading is assumed to be zero and is not shown above.) For simplicity, no level-2 predictors are presented in (1.2) or (1.3), but could be included as predictors of the intercept or slope variables. In the level-2 equations, 0 and 1 are the intercepts or average value of 0 and 1 respectively, and 0 and 1 are error terms.

More traditionally, the structural model would be represented by grouping each variable into matrices:

/ (1.4)
/ (1.5)

In these equations,  is a 2 X t matrix representing the relationship between 2 latent variables, 0 and 1, and t indicators, one for each time point. The column in  which corresponds to 0 is comprised of all 1s, because each loading for this variable is set equal to 1 to define it as the intercept.  is the matrix of measurement errors for each indicator at each time point.  is a 2 X 1 vector containing latent means for the intercept and slope, representing the average intercept across individuals and the average slope (i.e., trajectory) across individuals.  is the error term. The variance of the intercepts and slopes, 0 and 1, are obtained by estimation of the  matrix.

Multilevel Regression Models (Hierarchical Linear Models: HLM)

Multilevel regression models estimate predictive relationships when the data are nested or hierarchically structured, as in the case of students nested within schools. The statistical model used for hierarchically structured data is the same statistical model used for longitudinal analysis of individual growth curves. With growth curve models, longitudinal data measurements are considered to be nested within individuals. In general, a multilevel regression with a single level-1 predictor and no level-2 predictors can be written with two sets of equations:

level-1 equation:

/ (1.6)

level-2 equations:

/ (1.7)
/ (1.8)

The first equation is the familiar regression equation, with r representing error or unexplained variance. The subscripts i and g indicate whether the value is for each individual or each group. In the second equation, the intercept values for each group serve as the dependent variable. For simplicity sake, there are no predictors in equation (1.7) or (1.8). In (1.7), 00 is the intercept (mean of all group intercepts), and u0 is the error or remaining variance. u0 can also be interpreted as the variance of the intercept values across groups. Since the intercepts represent adjusted means for each group (i.e., adjusting or controlling for the effects of xi), u0 is the variance of the adjusted means for each group. In the third equation (1.8), the slopes from the level-1 equation, 1g, serve as dependent variables for each group. 01 is the intercept in this equation, and represents the average of all slopes, 1g, interpreted as the average effect of xi on the dependent variable across all groups. u1 is the error term or the variance of the slopes across groups (i.e., the variability in the relationship between x and y across the groups). By substituting equations (1.7) and (1.8) into equation (1.6), the HLM model can be expressed as,

/ (1.9)

or, by rearranging the terms,

/ (1.10)

If growth curve models are tested, the level-1 x variable is replaced by time codes, xt (e.g., 0,1,2,3...t-1). The dependent variable at each time point is regressed on the time code at level-1. Instead of individuals nested within groups, repeated measures are nested within individuals. Level 2 consists of individuals rather than groups.

Comparing the SEM and HLM growth models

Both approaches to latent growth models are essentially equivalent. One important difference between the two approaches is that the SEM approach accounts for measurement error at each time point (in the  matrix). The parallels between the SEM and HLM approaches can be seen by comparing their algebraic formulas.

SEM / HLM
Level-1 / (1.1) / (1.6)
Level-2 / (1.2)
(1.3) / (1.7)
(1.8)

These formulas are parallel, although this may not be fully apparent at first glance. In SEM, loadings are used in place of level-1 regression coefficients. In equation (1.1), the level-1 intercept is represented by the product term it0i, which refers to the loadings for latent intercept variable and the intercept variable itself. The  matrix is analogous to the X matrix in matrix regression in which the first column is a vector of 1's used to produce the intercept. By setting the loadings for the intercept (0i) to 1, the product of the loadings and the intercept (it0i ) of (1.1) is simply equivalent to the intercept term, 0, of equation (1.6). The next term in equation (1.1), it1i, representing the slope factor, can also be considered identical to the slope in equation (1.6) as long as the loadings in the  matrix are set to values that would be used as predictors in growth curve analysis, such as 0, 1, 2, . . . t-1. Here, it in (1.1) is equivalent to xt in (1.6). Because 's are equivalent to x's and the 's are equivalent to 's, it would make more sense to re-express equation (1.1) as,

/ (1.11)

By estimating the means and variances of 0i and 1i , we can obtain estimates of the average latent intercept and average latent slope and the extent to which they vary across individuals.

A Multilevel Structural Equation Model for Dyadic Data

Because growth models (repeated measures) and two-level hierarchical regression models are identical statistical models, as I've shown above, it is possible to specify a multilevel SEM for certain hierarchical data situations which uses the same model specifications as those used in the growth model case. I start by describing a dyadic data situation (e.g., couples, twins, mother-child dyads), for which this approach is clearly the simplest and possibly the best suited. I describe the model specification and give an example. I then discuss the possibilities of applying this technique to other data analytic situations. The generalization of the approach should be relatively straight forward for small groups sizes where the sample size is balanced (equal in all groups). Finally, I will discuss the possibilities of generalizing to situations in which group sizes are not equal but remain fairly small.

The usual multilevel approach has been to set up between and within covariance matrices, corresponding to level-1 and level-2 variables, and use a multigroup strategy to estimate the model (eg., Muthen, 1997; Muthen and Satorra, 1995). This approach can be implemented in any of the current SEM software packages by creating between and within covariance matrices, which are then analyzed in a multigroup structural model. This process has been recently automated in Mplus. This approach has a few limitations. First, except for Mplus, the procedure of constructing and reading separate covariance matrices can be cumbersome and analytically intensive, especially for larger models. Second, the multigroup approach is limited to separate within and between models. That is, one cannot analyze relationships between level-1 and level-2 variables. Multilevel regression, in contrast, allows for prediction of level-1 intercepts by level-2 predictors and for prediction of level-1 slopes by level-2 variables (called cross-level interactions"). Third, with small groups, nonconvergence of multilevel structural equation models can be a problem, especially with lower intraclass correlations (see Muthen & Satorra, 1995).

Data requirements. At minimum, one merely needs to have a single dependent measure obtained from each member of the dyad. Multiple indicators for each individual can also be used (see Specifications 2 and 3 below). For dyadic data, it would be optimal to have at least three indicators for each member of the couple. Members of each couple must also be nonexchangable. That is, there must be a basis for distinguishing members of each couple in an identical manner in all groups. Examples might include husbands and wives, mother and child, first born and second born, or caregiver and care recipient. The data set is set up in a so-called "repeated measures" format, in which each case in the data matrix contains information about the dyad. For instance, each record contains information about the husband and the wife, recorded under different variable names (e.g., y1h, y2h, y3h, y4w, y5w, y6w).

Example. To illustrate, I will use an example from a study I conducted recently examining interactions between spousal caregivers and care recipients. There are 118 couples (236 individuals), in which each member of the couple was interviewed separately. I examine five items from the Veit and Ware (1983) positive affect subscale of the Mental Health Inventory. Items such as "How much of the time have you felt the future look hopeful and promising?" on a 6-point scale of frequency of occurrence. Thus the analysis is based on 10 variables--5 items for caregivers and 5 items for care recipients.

Model specification 1. I first take the simplest case in which there is only one measure (i.e., indicator) of positive affect for caregivers and for care recipients (the measure was computed by averaging the five items for each). This model specification follows that of the growth curve model described above in the case in which there is only two time points tested. The first approach is depicted in Figure 1 below. Two latent variables are defined: a latent intercept, 0, and a latent slope, 1,. There is only one indicator for each of these latent variables. The loadings are fixed to specified values and the measurement errors are fixed to zero. The intercept variable, 0, is defined by fixing loadings on each of the two indicators to 1. The slope variable, 1, is defined by fixing loadings on the same two indicators to 0 and 1. The average intercept and average slope across couples (i.e., 00 and 01 in HLM notation) are obtained by estimating the mean structures of each (not estimated by default in SEM software packages). The variances of the slopes and intercepts across couples are indicated by the variances of 0 and 1 (i.e., the PSI matrix).


Results. Parallel analyses were conducted using Mplus (Muthen & Muthen, 1999) and HLM 5 (Raudenbush, Bryk, Cheong, & Congdon, 2000). In HLM, it is not possible to estimate variances (i.e., random effects) for the slopes and intercepts simultaneously with dyadic data. One can obtain estimates by running separate models fixing the random effect of either the intercept of the slope to zero, but these estimates will differ from a simultaneous estimation of both. Therefore, they are not presented here. Analyses are presented for dummy coding of the caregiving variable (0 and 1) and for group-mean centering. When dummy coding is used, the intercept represents the average score for caregivers (because they were coded as 0). When group-centering is used, the average intercept represents the grand mean for all couples (caregivers and care recipients combined). In multilevel regression, group-mean centering is achieved by subtracting each individual’s score from the mean of the dyad (this can be done automatically in the HLM software). To obtain the group-mean centered solution using the SEM approach, however, one simply sets the loadings of the slope variable to -.5 and +.5, rather than 0 and 1.

As can be seen in Table 1, the means (and their standard errors) for the intercept and slope variable are highly similar SEM method and the HLM method. The average intercept, approximately 4.2, represents the average positive affect for caregivers. Table 2 presents results when a group-centered approach is used. Notice that the average intercept obtained with centering differs little from that obtained with dummy coding. This is because there is very little difference between caregivers and care recipients on positive affect scores. Thus the average of caregivers and care recipients combined is similar to the average for caregivers only. The average slope (-.056 using the SEM method) represents the relationship between the difference variable (caregivers vs. care recipients) and the dependent variable, and is identical to a test between the caregiver and care recipient means on positive affect.

Table 1. Mean and variance estimates when the caregiving variable is uncentered/dummy coded (0,1).

Means (SE) / Variances
HLM / SEM / HLM / SEM
Intercept / 4.221 (.088) / 4.214 (.082) / NA / .791 (.103)
Slope / -.031 (.126) / -.056 (.099) / NA / 1.141 (.149)

Table 2. Mean and variance estimates when the caregiving variable is group-mean centered (-.5,+.5).

Means (SE) / Variances (SE)
HLM / SEM / HLM / SEM
Intercept / 4.205 (.072) / 4.186 (.069) / NA / .562 (.073)
Slope / -.031 (.103) / -.056 (.099) / NA / 1.141 (.149)

Model specification 2. The above approach has two shortcomings: it assumes no measurement error and only single indicators are used. Another specification is possible when there are several indicators. In this specification true latent variables are used for the intercept, 0, and a latent slope, 1,. Each of the ten indicators (5 caregiver items, 5 care recipient items) define the latent intercept and each loading is set equal to 1. The slope variable represents a difference variable, or dummy variable, distinguishing between the caregiver and care recipient. In Figure 2, I have coded caregivers as 0 and care recipients as 1. To accomplish this, the loadings for each of the caregiver items are set to 0 and each of the loadings for care recipients is set to 1 for the slope factor. The slope variable, 1, then represents the difference between caregivers and recipients on positive affect. In this specification, each of the five items are assumed to be equivalent indicators for caregivers and for care recipients (as is the case when computes an equally weighted composite score). Because of the 0-1 coding used for the latent slope variable, the intercept variable represents the latent mean for caregivers. One could also use an effects code (-1,+1 or -.5,+.5) for the slope variable, to center the variable within each group. In this case, the latent intercept would represent the average score for each couple. To account for intercorrelations between parallel items for caregivers and care recipients, I also estimate the correlations among measurement errors. The error correlations represent association between the errors over and above the variance accounted for by the the intercept and the difference factor. Ideally, one should test these correlations for significance. If they are not significant, the correlations can be set to zero. As in growth curve analysis, the individual intercepts in the measurement equation, , are set to zero

The model shown in Figure 2 below, corresponds to a simple multilevel regression model, with one level-1 predictor and no level-2 predictors as in equations (1.6) through (1.8).


Results. Table 3 and Table 4 present the results previously obtained from HLM and new results using Mplus and Specification 2. The most notable difference between results obtained with Specification 1 and Specification 2 is that the variances (or random effects) are smaller under Specification 2. This is the result of estimation of measurement error by the use of true latent variables for the intercept and slope.

Table 3. Mean and variance estimates when the caregiving variable is uncentered/dummy coded (0,1).

Means (SE) / Variances
HLM / SEM / HLM / SEM
Intercept / 4.221 (.088) / 4.253 (.075) / NA / .538 (.086)
Slope / -.031 (.126) / -.043 (.096) / NA / .824 (.142)

Similar analyses were conducted using group-mean centering of the caregiving variable. Using SEM, the equivalent model to the group-mean centered model in HLM constrains the loadings on 1to -.5 and +.5. With group-mean centering, the intercept value is the mean of both groups.

Table 4. Mean and variance estimates when the caregiving variable is group-mean centered (-.5,+.5).

Means (SE) / Variances (SE)
HLM / SEM / HLM / SEM
Intercept / 4.205 (.072) / 4.226 (.068) / NA / .483 (.071)
Slope / -.031 (.103) / -.107 (.095) / NA / .820 (.139)

Model specification 3. In Model Specification 2 above, loadings on 1 were assumed to be equal for caregivers and for care recipients. An alternative specification, using a second-order factor model for the slope variable would not require this assumption. Using this model specification, the loadings for the first order factor can be freely estimated, allowing for unequal loadings for each item. Two second-order factors, 0 and 1 can be used to represent the intercept and slope. As before, an uncentered, dummy coding (0,1) or group-mean centered coding (-.5,+.5) can be used to provide different interpretations for the intercept variable 0 . If dummy coding is used, the intercept represents the mean for the member of the dyad who is assigned 0 for the slope variable; and, if centered codes are used, the intercept represents the mean for the full sample. Several details of the specification are important. First, the intercepts for the first order latent variables (cg and cr) should be constrained to zero. Second, it is best to constrain loadings for the first order factors to be equal for parallel items. Third, the disturbances and the correlation for the first order factors should be set to zero. Fourth, the individual intercepts in the first order measurement equation, , are set to zero. Figure 3 below illustrates the model using dummy coding.