The data set used in this chapter is aggression.dat. See http://bear.soe.berkeley.edu/EIRM for more data for this book (aggression.dat is the same as data verbal aggression poly.txt).
The do-file is aggression.do.
The programs we use are gllamm and gllapred. You can download these programs by issuing the command ssc describe gllamm and following instructions; see here for more information on installing the GLLAMM programs. For more information on these commands and a manual, see http://www.gllamm.org.
Read data
insheet using aggression.dat, clear list person gender anger y i1-i12 in 1/12, clean person gender anger y i1 i2 i3 i4 i5 i6 i7 i8 i9 i10 i11 i12 1. 1 1 20 0 1 0 0 0 0 0 0 0 0 0 0 0 2. 1 1 20 0 0 1 0 0 0 0 0 0 0 0 0 0 3. 1 1 20 0 0 0 1 0 0 0 0 0 0 0 0 0 4. 1 1 20 0 0 0 0 1 0 0 0 0 0 0 0 0 5. 1 1 20 0 0 0 0 0 1 0 0 0 0 0 0 0 6. 1 1 20 0 0 0 0 0 0 1 0 0 0 0 0 0 7. 1 1 20 0 0 0 0 0 0 0 1 0 0 0 0 0 8. 1 1 20 0 0 0 0 0 0 0 0 1 0 0 0 0 9. 1 1 20 1 0 0 0 0 0 0 0 0 1 0 0 0 10. 1 1 20 2 0 0 0 0 0 0 0 0 0 1 0 0 11. 1 1 20 0 0 0 0 0 0 0 0 0 0 0 1 0 12. 1 1 20 0 0 0 0 0 0 0 0 0 0 0 0 1 list person gender anger y i13-i24 in 13/24, clean person gender anger y i13 i14 i15 i16 i17 i18 i19 i20 i21 i22 i23 i24 13. 1 1 20 1 1 0 0 0 0 0 0 0 0 0 0 0 14. 1 1 20 0 0 1 0 0 0 0 0 0 0 0 0 0 15. 1 1 20 1 0 0 1 0 0 0 0 0 0 0 0 0 16. 1 1 20 1 0 0 0 1 0 0 0 0 0 0 0 0 17. 1 1 20 0 0 0 0 0 1 0 0 0 0 0 0 0 18. 1 1 20 0 0 0 0 0 0 1 0 0 0 0 0 0 19. 1 1 20 1 0 0 0 0 0 0 1 0 0 0 0 0 20. 1 1 20 0 0 0 0 0 0 0 0 1 0 0 0 0 21. 1 1 20 0 0 0 0 0 0 0 0 0 1 0 0 0 22. 1 1 20 2 0 0 0 0 0 0 0 0 0 1 0 0 23. 1 1 20 2 0 0 0 0 0 0 0 0 0 0 1 0 24. 1 1 20 2 0 0 0 0 0 0 0 0 0 0 0 1
Recode the ordinal response variable into a dichotomous response.
recode y 2=1
We first estimate the Rasch Model as described in Section 2.3.2. Table 12.3 on page 867 gives estimates using 20-point ordinary quadrature. We use adaptive quadrature here because it requires fewer quadrature points for a given precision and is therefore faster. Specifically, we use 5-point adaptive quadrature although this may not be quite sufficient.
gllamm y i1-i24, nocons i(person) link(logit) family(binom) adapt nip(5) Running adaptive quadrature Iteration 0: log likelihood = -4250.4922 Iteration 1: log likelihood = -4048.0262 Iteration 2: log likelihood = -4037.1085 Iteration 3: log likelihood = -4036.9365 Iteration 4: log likelihood = -4036.9304 Iteration 5: log likelihood = -4036.9304 Adaptive quadrature has converged, running Newton-Raphson Iteration 0: log likelihood = -4036.9304 Iteration 1: log likelihood = -4036.9304 (backed up) Iteration 2: log likelihood = -4036.9256 Iteration 3: log likelihood = -4036.9256 number of level 1 units = 7584 number of level 2 units = 316 Condition Number = 6.0244357 gllamm model log likelihood = -4036.9256 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- i1 | 1.22025 .1631531 7.48 0.000 .9004757 1.540024 i2 | .5642088 .1547349 3.65 0.000 .260934 .8674836 i3 | .0796954 .1526948 0.52 0.602 -.2195809 .3789717 i4 | 1.747735 .1754059 9.96 0.000 1.403946 2.091524 i5 | .7071053 .1559911 4.53 0.000 .4013685 1.012842 i6 | .011342 .1526763 0.07 0.941 -.287898 .310582 i7 | .5289166 .1544716 3.42 0.001 .2261579 .8316753 i8 | -.6865499 .1563173 -4.39 0.000 -.9929263 -.3801736 i9 | -1.527167 .1708641 -8.94 0.000 -1.862054 -1.192279 i10 | 1.081277 .1607845 6.73 0.000 .7661452 1.396409 i11 | -.3496521 .1536777 -2.28 0.023 -.6508548 -.0484493 i12 | -1.044171 .1610314 -6.48 0.000 -1.359786 -.7285549 i13 | 1.22025 .1631531 7.48 0.000 .9004757 1.540024 i14 | .3891213 .1536078 2.53 0.011 .0880555 .6901871 i15 | -.8713759 .158498 -5.50 0.000 -1.182026 -.5607256 i16 | .8719614 .1578279 5.52 0.000 .5626244 1.181298 i17 | -.0570312 .1527238 -0.37 0.709 -.3563644 .2423019 i18 | -1.482074 .169763 -8.73 0.000 -1.814804 -1.149345 i19 | -.2113447 .1530745 -1.38 0.167 -.5113653 .0886759 i20 | -1.504531 .1703064 -8.83 0.000 -1.838325 -1.170737 i21 | -2.975819 .2321202 -12.82 0.000 -3.430766 -2.520871 i22 | .7071053 .1559911 4.53 0.000 .4013685 1.012842 i23 | -.3844778 .153873 -2.50 0.012 -.6860634 -.0828922 i24 | -1.999901 .184953 -10.81 0.000 -2.362403 -1.6374 ------------------------------------------------------------------------------ Variances and covariances of random effects ------------------------------------------------------------------------------ ***level 2 (person) var(1): 1.920391 (.19507899) ------------------------------------------------------------------------------
Note that the item parameter estimates have the opposite sign as in Table 12.3 because we used the parameterization logit[Pr(yip=1|θp)] = θp + βi instead of logit[Pr(yip=1|θp)] = θp - βi. The estimated person variance (1.92) is somewhat lower than in Table 12.3 but agrees with the variance and standard error obtained using SAS NLMIXED with adaptive quadrature (see page 368).
On page 57, there is a brief discussion of person measurement and the associated precision. We can use gllapred to obtain expected a posteriori (EAP) or empirical Bayes predictions. These are just means of the posterior distributions of the person parameters given the item responses, with parameter estimates pugged in. The corresponding posterior standard deviations are often interpreted as standard errors in item response theory (IRT).
gllapred score, u list person scorem1 scores1 if i1==1&person<21, clean person scorem1 scores1 1. 1 -.47821717 .44810137 25. 2 -2.620004 .66216859 49. 3 -.27943152 .44387067 73. 4 .50677069 .4481273 97. 5 -.27943152 .44387067 121. 6 -.08340467 .44187609 145. 7 .50677069 .4481273 169. 8 -1.6093829 .52012249 193. 9 -.27943152 .44387067 217. 10 1.8899446 .53219586 241. 11 .71026818 .45437453 265. 12 1.1399934 .47437231 289. 13 -.27943152 .44387067 313. 14 1.6199691 .50790979 337. 15 .30790344 .44402597 361. 16 .11177496 .4419504 385. 17 .71026818 .45437453 409. 18 -.47821717 .44810137 433. 19 -3.1158026 .75530537 457. 20 .11177496 .4419504
scorem1 are the EAP scores and scores1 the corresponding standard errors
We now estimate the latent regression Rasch model as described in Section 2.4.2, using estimates from the Rasch model as starting values (zero will be substituted for the new parameters).
matrix a=e(b) eq f1: anger gender gllamm y i1-i24, nocons i(person) link(logit) family(binom) geqs(f1) from(a) adapt nip(5) Running adaptive quadrature Iteration 0: log likelihood = -4036.9254 Iteration 1: log likelihood = -4030.1465 Iteration 2: log likelihood = -4029.9865 Iteration 3: log likelihood = -4029.987 Adaptive quadrature has converged, running Newton-Raphson Iteration 0: log likelihood = -4029.987 Iteration 1: log likelihood = -4029.987 Iteration 2: log likelihood = -4029.981 Iteration 3: log likelihood = -4029.981 number of level 1 units = 7584 number of level 2 units = 316 Condition Number = 436.55046 gllamm model log likelihood = -4029.981 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- i1 | -.0120518 .3796101 -0.03 0.975 -.756074 .7319704 i2 | -.6683203 .3767924 -1.77 0.076 -1.40682 .0701792 i3 | -1.152842 .3764523 -3.06 0.002 -1.890675 -.4150085 i4 | .5158246 .3844048 1.34 0.180 -.2375949 1.269244 i5 | -.525397 .377159 -1.39 0.164 -1.264615 .2138211 i6 | -1.221187 .3765124 -3.24 0.001 -1.959138 -.483236 i7 | -.7036172 .3767212 -1.87 0.062 -1.441977 .0347427 i8 | -1.918899 .3786809 -5.07 0.000 -2.6611 -1.176698 i9 | -2.759257 .3857852 -7.15 0.000 -3.515382 -2.003131 i10 | -.1510957 .3787571 -0.40 0.690 -.893446 .5912545 i11 | -1.582106 .3772714 -4.19 0.000 -2.321545 -.842668 i12 | -2.276395 .3810068 -5.97 0.000 -3.023154 -1.529635 i13 | -.0120518 .3796101 -0.03 0.975 -.756074 .7319704 i14 | -.8434248 .3765123 -2.24 0.025 -1.581375 -.1054743 i15 | -2.10366 .3797684 -5.54 0.000 -2.847993 -1.359328 i16 | -.3604948 .3777444 -0.95 0.340 -1.10086 .3798707 i17 | -1.28955 .3765989 -3.42 0.001 -2.02767 -.5514295 i18 | -2.714173 .385247 -7.05 0.000 -3.469243 -1.959102 i19 | -1.443833 .3768919 -3.83 0.000 -2.182528 -.7051387 i20 | -2.736625 .3855125 -7.10 0.000 -3.492215 -1.981034 i21 | -4.208565 .4186407 -10.05 0.000 -5.029086 -3.388045 i22 | -.525397 .377159 -1.39 0.164 -1.264615 .2138211 i23 | -1.616922 .3773847 -4.28 0.000 -2.356583 -.877262 i24 | -3.231983 .3928153 -8.23 0.000 -4.001887 -2.462079 ------------------------------------------------------------------------------ Variances and covariances of random effects ------------------------------------------------------------------------------ ***level 2 (person) var(1): 1.8304066 (.18663415) Regressions of latent variables on covariates ------------------------------------------------------------------------------ random effect 1 has 2 covariates: anger: .05786975 (.01702726) gender: .32265841 (.19398198) ------------------------------------------------------------------------------
The person property effect estimates reported on page 60 (0.057 (SE=0.016) for anger and 0.29 (SE=0.20) for gender) differ a little from the estimates here.
We can again obtain EAP scores. Using the fac option gives predictions based on the two the explanatory variables in the latent regression as well as the item responses, whereas the u option gives predictions of the residuals in the latent regression.
gllapred a, fac gllapred b, u sort person list person am1 bm1 if i1==1&person<21, clean person am1 bm1 1. 1 .78366929 -.69638403 25. 2 -1.4222138 -2.3814394 49. 3 .9277491 -.05603657 73. 4 1.7346317 .51936699 97. 5 .9277491 -.05603657 121. 6 1.147595 -.06766963 145. 7 1.8490626 -.40785743 169. 8 -.36809639 -1.583361 193. 9 .97115 -.41772389 217. 10 3.061921 2.1360051 241. 11 1.8980956 1.0300494 265. 12 2.3424642 1.3008088 289. 13 1.0796236 -1.3263457 313. 14 2.8925603 1.2967675 337. 15 1.5058998 .57998391 361. 16 1.3233475 .28169205 385. 17 1.9498918 .61888764 409. 18 .72317424 -.20274168 433. 19 -1.8426932 -3.0579578 457. 20 1.366363 -.08038062
We now estimate the linear logistic test model (LLTM) described in Section 2.5.2. First we have to define covariates for the item properties (see Table 1.1, page 9, and Figure 2.8, page 63).
bysort person (i1-i24): gen item=25-_n gen do_want = cond(item>12,1,0) gen other_self = cond(item<7|(item>12&item<19),1,0) gen blame = cond(mod(item,3)==0,-1,0.5) gen express = cond(mod(item,3)==2,-1,0.5) gllamm y do_want other_self blame express, i(person) link(logit) family(binom) adapt nip(5) Running adaptive quadrature Iteration 0: log likelihood = -4317.1946 Iteration 1: log likelihood = -4125.9568 Iteration 2: log likelihood = -4116.8186 Iteration 3: log likelihood = -4116.6333 Iteration 4: log likelihood = -4116.6333 Adaptive quadrature has converged, running Newton-Raphson Iteration 0: log likelihood = -4116.6333 Iteration 1: log likelihood = -4116.6333 (backed up) Iteration 2: log likelihood = -4116.6333 number of level 1 units = 7584 number of level 2 units = 316 Condition Number = 2.905931 gllamm model log likelihood = -4116.6333 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- do_want | -.6714401 .0571109 -11.76 0.000 -.7833754 -.5595049 other_self | 1.027672 .0580049 17.72 0.000 .9139847 1.14136 blame | 1.361158 .0499605 27.24 0.000 1.263237 1.459079 express | .7032311 .0462105 15.22 0.000 .6126601 .7938021 _cons | -.3169545 .0904736 -3.50 0.000 -.4942795 -.1396295 ------------------------------------------------------------------------------ Variances and covariances of random effects ------------------------------------------------------------------------------ ***level 2 (person) var(1): 1.811924 (.18426981) ------------------------------------------------------------------------------
The person variance is lower than for the Rasch model (and also lower than that reported on page 64, possibly due to using adaptive quadrature). The estimates of the item property effects and their standard errors are the same to two decimal places as those given on pages 64 to 65, with the signs of the coefficients again reversed.
Finally, we combine the LLTM model with a latent regression as described in Section 2.6.2, using previous estimates as starting values.
matrix a=e(b) eq f1: anger gender gllamm y do_want other_self blame express, i(person) link(logit) family(binom) geqs(f1) from(a) adapt nip(5) Running adaptive quadrature Iteration 0: log likelihood = -4116.6334 Iteration 1: log likelihood = -4110.5943 Iteration 2: log likelihood = -4109.6688 Iteration 3: log likelihood = -4109.6698 Adaptive quadrature has converged, running Newton-Raphson Iteration 0: log likelihood = -4109.6698 Iteration 1: log likelihood = -4109.6698 (backed up) Iteration 2: log likelihood = -4109.6681 Iteration 3: log likelihood = -4109.6681 number of level 1 units = 7584 number of level 2 units = 316 Condition Number = 90.392126 gllamm model log likelihood = -4109.6681 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- do_want | -.6714664 .0571117 -11.76 0.000 -.7834033 -.5595294 other_self | 1.027699 .0580048 17.72 0.000 .9140115 1.141386 blame | 1.36115 .0499588 27.25 0.000 1.263233 1.459068 express | .703301 .0462128 15.22 0.000 .6127257 .7938764 _cons | -1.516946 .3473162 -4.37 0.000 -2.197673 -.8362187 ------------------------------------------------------------------------------ Variances and covariances of random effects ------------------------------------------------------------------------------ ***level 2 (person) var(1): 1.7267964 (.17625438) Regressions of latent variables on covariates ------------------------------------------------------------------------------ random effect 1 has 2 covariates: anger: .05632505 (.01657122) gender: .31642333 (.18875432) ------------------------------------------------------------------------------
The person variance is estimated as 1.73 (same as on page 67).