** Three-level Normal Model, math achievement data **
set mem 400m

use http://www.stata-press.com/data/mlmus/achievement, clear


** explore our data a little **

gen yr=year+3.5
xtdes, i(child) t(yr) patterns(30)

xtsum math, i(school)
** sd between is quite different from sd within
xtsum math, i(child)
** sd between and within is approximately the same

** these don't capture fact that children are nested in schools **


******************************************************************
** start off using two-level variance components model ***********
** with a random intercept for school, not accounting for child **
******************************************************************
xtreg math year, i(school) nolog mle
xtmixed math year || school:, nolog mle
estimates store modelS
** above two specifications are equivalent **
gllamm math year, i(school) nolog

** if the gllamm estimation of the random effect variances is a little
** different from those estimated by xtreg and xtmixed, we can use 
** nip() and adapt
** THIS TAKES A LITTLE WHILE TO RUN !!!  (so I commented it out)

*gllamm math year, i(school) nip(15) adapt

* now compare
** sd of random intercept
display sqrt(.20724444 )
** sd within-school
display sqrt(.9978509) 

** this model assumes that math scores are a linear function of 
** time and, conditional on a school, the math scores within 
** this school are independent...not very reasonable!

******************************************************************
** start off using two-level variance components model ***********
** with a random intercept for child, not accounting for school **
******************************************************************
xtmixed math year || child:, nolog mle
estimates store modelC

** this model assumes that math scores are a linear function of 
** time and, conditional on a child, the repeated math scores 
** are independent. This is an okay model, but doesn't 
** take into account clustering according to school


*********************************************************************************
** Three-level variance components model, random intercept for school and child**
xtmixed math year|| school: || child:, nolog mle
estimates store modelSC

** test null hypothesis that variance components for child is zero **
lrtest modelS modelSC
** we need to include the random effect for child

** test null hypothesis that variance components for school is zero **
lrtest modelC modelSC
** we need to include the random effect for school

***************************************************
** intraclass correlations for three-level model **
***************************************************


------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
school: Identity             |
                   sd(_cons) |   .4280823   .0462896      .3463257    .5291391
-----------------------------+------------------------------------------------
child: Identity              |
                   sd(_cons) |   .8184857   .0160566      .7876127    .8505689
-----------------------------+------------------------------------------------
                sd(Residual) |   .5890159   .0056111      .5781204    .6001168

** rho(school), ICC between measurements from same school different child
display .43^2 /(.43^2 + .82^2 + .59^2)

** rho(school, child), ICC between measurements from same school same child
display (.43^2 + .82^2 )/(.43^2 + .82^2 + .59^2)

** note that rho(school, child) always > rho(school)


*********************************************************************************
** Three-level variance components model, random intercept for school and child**
** adding in covariates as fixed effects ****************************************

** add in some child-level covariates
xtmixed math year female hispanic black || school: || child:, nolog mle

** note that the sd of the random intercept for school changes, by controlling
** for these student-level characteristics, so the schools must be very different
** according to the make-up of their student body

** drop female

** add in some school-level covariates
xtmixed math year hispanic black lowinc size mobility|| school: || child:, mle

** the sd of the random intercept for school decreases even more, so we 
** have removed more of the variability between schools by controlling for these
** school-level covariates
** If we added in enough school-level covariates so that, having controlled for 
** these school-level covariates, the sd of the random intercept for school was 
** zero, we wouldn't need to include the random intercept for school
** The random intercept for school in effect "mops up" unexplained variability between schools

** do we need to include any of the covariates that control for SES?
test hispanic black lowinc mobility

** drop size

** add in a random slope on year at child level 

** at child-level the random parts are independent **
xtmixed math year hispanic black lowinc mobility|| school: || child: year, nolog mle
estat ic
** at child-level the random parts have unstructured corr matrix (allowing for correlation) **
xtmixed math year hispanic black lowinc mobility|| school: || child: year, cov(unstructured)nolog mle
estat ic

** stick with the unstructured since it minimizes AIC and BIC

** this correlation of .44 means that for children who tend to have a higher value of the random 
** intercept (a higher baseline math score), they also tend to have a higher random slope on year
** (they improve math scores at a greater rate)