HLM is a statistical technique for analyzing longitudinal data that have a hierarchical structure. For example, the longitudinal data can be a record of students' grades collected over time. Data collected this way are often in a hierarchical structure, e.g., students are grouped in classes within schools. Longitudinal data are sometimes called repeated-measure data. To analyze them, psychologists often use repeated-measure Analysis Of Variance (ANOVA). HLM is more powerful than repeated-measure ANOVA. One advantage is that HLM is more flexible, and it allows the analyst to use complex models to describe changes of outcomes over time. In this section we will learn how to use SPSS to analyze the mathematics achievement data in Bryk and Raudenbush (1992). A simplified background of the study is given first, followed by an example in SPSS. Downloadable SPSS data and syntax files are at the bottom of this page.

Mathematical Achievement in Bryk and Raudenbush (1992)

Bryk and Raudenbush wrote an influential book on HLM. In that book they used HLM to describe how socioeconomic status (SES) affects students' math achievement. The data were based on the mathematics achievement of over 7,000 students in 160 schools: 90 public and 70 Catholic. The basic ideas are simple: Students from low-SES families generally do less well than students from high-SES families. But there are subtle differences. It would not be surprising if we observe that students from low-SES families could excel equally well if they are in a nurturing environment with equal opportunities. Bryk and Raudenbush showed the subtle interactions between SES and the environment (that is, the school the students go to). HLM is particularly useful for this type of analysis.

Their main findings were that:

• Students in Catholic schools tend to have higher average SES, but
• High SES affects students in Catholic and public schools differently. High SES has a stronger impact on math achievemnt for students in public schools than for students in Catholic schools.

The second finding is of particular interest. But we will not go into its implications right now. Instead we will go over how to use SPSS to carry out Bryk and Raudenbush's analysis and verify their results.

A Hierarchy of Statistical Models

The statistical procedure is called hierarchical linear modeling because it involves two levels of regression-like models. The level-1 model describes how SES affects math achievement in a simple regression equation.

Y_ij = β_0j + β_1j (SES)_ij + ε_ij.

The equation states that math achievement is a function of SES plus unexplained random variations. The subscripts in the equation may appear intimidating to beginners. But they simply say that math achievement (Y_ij) for student i at school j depends on many predictors, or terms. The first predictor is an average math achievement for each school β_0j. Here β_0j represents the average math achievement score when SES is zero. This may appear nonsense at first, because it is hard to think of a family with zero SES. We will take care of it later. The β_0j parameter is referred to as the INTERCEPT. The other level-1 parameter, β_1j, represents the increase in math achieve per unit increase in SES. It is referred to as the SLOPE. Suppose that β_1j equals 3.2, then it means that on average, a unit increase in SES is estimated to boost math achievement by 3.2 points.

Note that the first subscript of β_0j, 0, represents a fixed number, an average; and the second subscript, j, represents several parameters, in fact, it represents one parameter for each of the j schools. Thus β_0j represents j averages, one for each school. The second predictor is the SES (β_1j (SES)_ij). Here β_1j represents that there is one slope coefficient for each of the j schools. The last predictor is random error residuals for person i and school j (ε_ij). Algerbraically, this is not very different from saying that a worker gets paid $8 per every hour's work, plus a basic pay of $20 for showing up every day. Here $8 per hour pay is analogous to the β_1j parameter, and $20 is analogous to β_0j. After you get used to the notations, you will soon feel less intimidated by them.

The INTERCEPT problem can be taken care of by rescaling. If for each school we calculate the average SES scores across students and call it MEANSES. Then MEANSES represents the overall SES for each school. Moreover, if we subtract MEANSES from each student's unique SES measure and call it CSES, then CSES represents by how much the student's SES differs from the school average. Zero CSES means that the student's SES is at the school average, not null SES any more. This subtraction procedure is also referred to as "centering". It does not affect the statistical tests. After the centering, the intercept parameter β_0j represents the estimated math achievement at school average SES. So that the level-1 equation now becomes:

Y_ij = β_0j + β_1j (CSES)_ij + γ_ij,
where CSES stands for centered SES.

The method is hierarchical because the parameters in the level-1 model (the equation above) are further dissected in a second level. In the second level β_0j and β_1j are further modeled by the type of schools the students go to. Bryk and Raudenbush use the second level of models to ask whether or not the impact of SES on math achievement is different between public and Catholic schools. The model equations are:

β_0j = γ₀₀ + γ₀₁ (MEANSES) + γ₀₂(Catholic) + u_0i, and
β_1j = γ₁₀ + γ₁₁ (MEANSES) + γ₁₂ (Catholic) + u_1i.

Bryk and Raudenbush call this the "intercept- and slope-as-outcomes" model. The intercept and slope in the level-1 model become outcome variables in the lower-level model. The level-1 intercept (the average SES after centering for each school) parameter β_0j now is a function of an average intercept for all schools (γ₀₀), average SES in each school (γ₀₁), an estimated difference between Catholic and public schools (γ₀₂), and random variations among students within the schools (u_0i). It is important to note that the variable MEANSES is calcuated for each school (averaged across students), a detail that will become imporant later in this page. A statistically significant γ₀₁ then suggests that a Catholic school with an average SES is expected to have a higher math achievement score than that of a public school with an average SES.

The slope parameter β_1j represents the "growth rate" in math achievement due to SES. At level-2, it is also a function of an average growth rate for all schools, average SES in the student's school, an estimated difference between Catholic and public schools, and random variations. A statistically significant γ₁₂ suggests that SES exerts differential impact on math achievement for students in Catholic schools than for students in public schools. Suppose that γ₁₂ turns out to be -2.3 (also with a p < 0.05), then for each unit increase in SES, its boost on math achievement is 2.3 points less in Catholic schools than in public schools. In other words, SES is less important in Catholic schools than in public schools.

SPSS Syntax

When we combine the two levels together and substitute β_0j and β_1j with the level-2 parameters, we end up with the following long equation:

Y_ij = [γ₀₀ + γ₀₁ (MEANSES) + γ₀₂(Catholic) + u_0i] + (CESE) [γ₁₀ + γ₁₁ (MEANSES) + γ₁₂ (Catholic) + u_1i] + γ_ij, and by distributing (CSES) into the terms within the square brackets, we get

Y_ij = γ₀₀ + γ₀₁ (MEANSES) + γ₀₂(Catholic) + u_0i + γ₁₀ (CSES) + γ₁₁ (MEANSES)*(CSES) + γ₁₂ (Catholic)*(CSES) + u_1i (CSES) + γ_ij.

The terms in bold font, u_0i + u_1i (CSES) + γ_ij, deserves special attention. They are called the "random effects" (as opposed to the "fixed effects"---all the other terms in the long equation, more on them later). The first term, u_0i represents the math achievement score for a student whose SES status is at the school average. The term u_1i (CSES) represents that there is a unique growth rate of CSES on math achievement for each student. And finally, γ_ij represents the residual variations for all students across all schools, variations not explained by any of the parameters in this complex model. In HLM jargon these terms and coefficients are called "random effects". This complicated concenpt is best explained with the graph below.

Random vs. Fixed Effects

The graph below explains the differences between a random and a fixed effect in a hierarchical linear model. There are 20 small squares in the graph, each represents a randomly selected school from the total of 160 schools in Bryk and Raudenbush's example. The top 10 squares are schools selected from the 70 Catholic schools and the bottom 10 are selected from the 90 public schools. On top of each panel there is a unique school id. In each square panel, we plot the relationship between math achievement scores and SES. Also, we fit a regression line like the level-1 model above.

Note that the regression lines are not the same. Some schools have a steep regression line (e.g., schools 1461 and 5619, two public schools in which SES strongly affects the math achievement of the students in the schools), while others have a relatively flat regression line (e.g., schools 9359 and 5404, two Catholic schools in which SES exerts a weak impact on math achievement). If we allow each school to have a unique slope and intercept in its regression line, then the INTERCEPT for each of the j schools (β_0j) can very randomly, so can the SLOPE for each school (β_1j). One of the greatest innovations in HLM is that these random variations between schools are further explained in the level-2 model equations, by MEANSES (the average SES for each school, averaged across individual student's SES within that school) and by the variable Catholic (a dummy variable that separates Catholic and public schools). Therefore, the "random" factors are the parameters that are considered to change across the unit of analysis, which in this case are the schools. All other parameters are considered static across schools, or they are "fixed". As seen in the figure, the "random" versus "fixed" distinction depends on the data analyst's hypotheses on how the variables behave across the units of analysis.

SPSS Syntax for HLM

The SPSS syntax below mimics the equations above. FIXED effects and RANDOM effects are separately specified. The grouping of the subjects (and of MEANSES) are at the SCHOOL level. Other parts of the SPSS syntax are mostly for controlling how SPSS will fit the HLM model. Details on their use can be found on the PDF file called the SPSS Syntax Guide (accessible from the Help menu).

* ****************
* RECODE variable sector (1 = public; 2 = Catholic) to a dummy
* variable Catholic.
* ****************
RECODE
  sector
  (2=1)  (1=0)  INTO  Catholic .
VARIABLE LABELS Catholic 'Catholic school'.
EXECUTE .

MIXED
  mathach  WITH Catholic meanses cses
  /FIXED = meanses cses cses*meanses Catholic Catholic*cses  | SSTYPE(3)
  /RANDOM INTERCEPT cses | SUBJECT(school) COVTYPE(UN)
  /METHOD = REML
  /CRITERIA = CIN(95) MXITER(100) MXSTEP(5) SCORING(1)
  SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE)
  LCONVERGE(0, ABSOLUTE)
  PCONVERGE(0.000001, ABSOLUTE)
  /PRINT = CORB COVB G SOLUTION .

The menu-driven, point-and-click method is not recommended. As of SPSS version 12 for Windows, the point-and-click method seems to omit the /RANDOM component of the syntax. The /RANDOM command is important because it says the unit of analysis is SUBJECT(school), and it allows one INTERCEPT and CSES slope for each school. Hence the following syntax:

    /RANDOM INTERCEPT cses | SUBJECT(school) COVTYPE(UN)

Downloads

The SPSS data can be downloaded below. You can verify the SPSS output with the results on pp. 70 - 75 of Bryk and Raudenbush (1992). NOTE: The data file is in ASCII text format, not in SPSS binary format. It is a good idea to share data in plain text because plain text is not restricted by the particular SPSS version you have. If your SPSS is an older version than mine, you can't read the binary .sav file I create. But you can always use the syntax file below to read plain data into any version of SPSS (including SPSS on a Macintosh and on a UNIX-based operating system).

Download SPSS data. NOTE: Use your "Save Target As" option (right mouse button on a Windows machine running IE) to save the file on your hard disk.

Download SPSS syntax. NOTE: This syntax file contains a "DATA LIST FILE = 'Bryk.sav'" command that processes the plain text data into SPSS. Download both the syntax and the Bryk.sav files into the same working directory, laungh SPSS, read the syntax into SPSS. You should have the full dataset when you Run the syntax.

R and nlme() for HLM

Finally, the same analysis can be carried out in R (www.r-project.org). The graph above is produced by the xyplot() function in R. R is an open-source computer language designed for statistical programming. R can be downloaded for free. R comes with a basic system and some required packages. Additional packages can be downloaded and incorporated into the base R system. The nlme library (Non-Linear Mixed-Effect Model) is one of the add-on packages. The nlme library is designed to analyze and visualize hierarchical data. The Bryk and Raudenbush example is in the MathAchieve dataset in the nlme() package by Jose Pinheiro and Douglas Bates. The example can be modeled in R with the following syntax, which is analogous to the SPSS syntax above.

library(nlme)


     lme(mathach ~ meanses*cses + sector*cses,
            random = ~ cses | school, data=Bryk)

Details on how to use R to analyze Bryk and Raudenbush's example can be found in John Fox's chapter (a link to it is in the bibliography below) and are not duplicated here.

References

Bryk, A. S. & Raudenbush, S. W. (1992). Hierarchical Linear Models : Applications and Data Analysis Methods. Newbury Park, CA: Sage Publications, Inc.

There is a second edition by Raudenbush and Bryk (2001).

The SPSS data came from the MathAchieve dataset in the nlme library of the statistical package R (www.r-project.org). The nlme library is written by Pinheiro and Bates and made freely available to the R community. I also relied on the foreign() function in R (written by Thomas Lumley) to port the MathAchieve dataset from R into the downloadable SPSS files above. Kudos to all.

A more sophisticated treatment of hierarchical linear model can be found in John Fox's chapter online (http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-mixed-models.pdf, last accessed December 6, 2005.).