Practice session8HLM多層線性模型講義

1、R-practice session 8cs&ss 560Marijtje van Duijn Winter 2006The commands used in this session are available as R syntax file (Session8.R) at the website.Data input and preparationWe continue with the data used in Snijders & Bosker. For a description see Example 4.1 (p. 46). We will estimate the multi

2、variate models for the language and math scores and reproduce the tables in chapter 13.Download the data file SBch13.csv from the class website. Also get the file session8.r and execute the -by now wellknown- commands under data preparation.Transforming the data for a multivariate multilevel analysi

3、sFor a multivariate multilevel analysis, the data need to be reorganized, comparable to a repeated measures analysis, where each measurement (of each dependent variable) has its own row.This is done (see session8.r) by first creating two matrices with just one dependent variable (and identical colum

4、ns for all other variables). Note that an indicator variable is added, numbering the matrices consecutively, to be able to distinguish the dependent variables later on. In general you need as many matrices as you have dependent variables. These matrices are then stacked and reordered, using the indi

5、cator variable, so that the correct structure is obtained. This takes a bit of moving around the data in new matrices and defining data frames (maybe this can be done more efficiently), but it works. The variable used to indicate the multiple testscores is called testNR. The data matrix with which w

6、e will work further is now twice as long as the one used before.Estimating a multivariate multilevel modelWe will first estimate the empty multivariate model as specified in table 13.1. This is not so hard, the only thing we have to remember is that it is now - technically - a three-level model, whe

7、re testNR indicates the lowest level, and pupilNR and schoolNR the second and third, and that we want to estimate separate variances for each value of testNR, which we therefore treat as an (ordered) factor.model13.1model13.1ML-lme(scorefactor(testNR)-1, +random=factor(testNR)-1lschoolNR/pupilNR,dat

8、a=langaritmv,method=ML) summary(model13.1ML)VarCorr(model13.1ML)There is one important difference with the table, and that is that in addition to the level 2 (between-pupils) within-schools covariance matrix elements a level 1 residual variance (a sort of measurement variance) is estimated. This num

9、ber needs to be added to the level 2 variances to get comparable results. Try to figure out how to get to the covariances reported in the table. (We will discuss in the lab, of course.) And note that the (population) correlation coefficients at the school level and at the pupil level are reported in

10、 the model summary.The correlations between the observed variables at pupil and at school level can serve as a reference to get some idea of explained variance at both levels.For the estimation of the model in table 13.2, the random structure is unchanged, but we have to be a little more careful in

11、the specification of the fixed part where, again, we want all estimates separately for the two outcome measures.How could one test whether the effect of for instance IQ on the test scores is different for the language and arithmetic tests?A model with a random slope for IQ did not converge, unfortun

12、ately, but in principle it is possible to estimate random slopes.We continue the inspection of residuals. We explored level 1 last week, and will now continue with level 2, following section 9.6.2 in the book.The first idea is to extend the concept of Cooks distance to a multilevel model. Cooks dist

13、ance measures the influence of an observation by investigating how much the estimated parameters change in relation to the precision or uncertainty of the estimate.This can also be done for the fixed parameters in a multilevel model (see page 134 for the formula) and the explanation in session8.r. H

14、ere I did program the time consuming way of this leaving-one-out procedure, whereas approximations exist, but these are not easily accessible in R (or at least I dont know how to obtain them).The random effects analogon of Cooks distance (9.6) can also not be obtained, because we dont have a (good)

15、estimate of the covariance matrix of the variance parameters.So we are left with the Cooks distance for the fixed part, which is not so bad, because there is another way to check the random part, and that is through the standardized multivariate residual for each level 2 unit, given by equation 9.9

16、on page 136. This statistic has an approximate chi-squared distribution so we can obtain a p-value for it, where a small p-value indicates an abnormally large residual. Because we obtain a statistic and a p-value for each level 2 unit, chance capitalization is a problem, so we should be careful not to overinterpret the p-values (or use a correction).Both diagnostics give us an idea of possibly outlying level 2 units. We can order the level 2 units by both criteria to investigate the worst cases. We can also plot Cooks distances by level 2 unit size, because theoretically Cooks distance is p