version 12.1
set more off
// Open National Longitudinal Survey (Young women 14-26 years of age in 1968)
webuse nlswork, clear
// Level 1: annual measurement
// Level 2: young women
// Null model
xtmixed ln_w || id: , var
xtmrho // VPC in random intercept case w/out covariates equals ICC
// Prepare level-1 covariate
center tenure, replace
label var tenure "Job tenure (years)"
// Add level-1 covariate
xtmixed ln_w tenure || id: , var
xtmrho // VPC after adjusting for a covariate
// 1) Constant variance
// Create variables: constant variance function at level 2
gen lev2var = e(var_u1)
label var lev2var "Variance at level 2"
// Create variables: constant variance function at level 1
gen lev1var = e(var_e)
label var lev1var "Variance at level 1"
// Create variable: predicted values
predict fit, xb
// Create variables: constant bounds for between-neighborhood variation
gen lev2cbhi = fit + 1.96*(lev2var^.5)
gen lev2cblo = fit - 1.96*(lev2var^.5)
label var fit "Fixed relationship"
label var lev2cbhi "Upper bound between-person variation"
label var lev2cblo "Lower bound between-person variation"
// Graphs: Heterogeneity as a constant function
qui twoway (line lev2var tenure, sort) ///
(line lev1var tenure, sort) ///
, legend(ring(0) pos(5)) ///
ylabel(.085 (.01) .12,format(%6.3f)) ///
ytitle("Variance") ///
name(variancefunction1, replace)
qui twoway (line fit lev2cbhi lev2cblo tenure, sort) ///
, legend(ring(0)) ///
ylabel(1 (.5) 3,format(%6.1f)) ///
ytitle("Predicted log wage") ///
name(betweenperson1, replace)
graph combine variancefunction1 between1
, col(1) ysize(8) xcommon ///
title("Random intercept model") ///
name(constantvariance, replace)
// Clean up
drop lev2var lev1var fit lev2cbhi lev2cblo
// Add random slope to model
xtmixed ln_w tenure || id: tenure, var cov(uns)
xtmrho // VPC after adjusting for a covariate and random slope
// 2) Quadratic variance
// Getting at the covariance parameter
// and slope variation as
// -xtmrho- hasn't stored it in a better format
local cov = tanh([atr1_1_1_2]_b[_cons])* exp([lns1_1_1]_b[_cons])* exp([lns1_1_2]_b[_cons])
di `cov'
local var_tenure = exp([lns1_1_1]_b[_cons])^2
// In order to get the other variance components
// without -xtmrho-, one needs:
*local var_cons = exp([lns1_1_2]_b[_cons])^2 // Intercept variance
*local var_res = exp([lnsig_e]_b[_cons])^2 // Residual variance
// Create variables: quadratic variance function at level 2
gen lev2var = e(var_u1) ///
+ 2 * `cov' * tenure ///
+ `var_tenure' * tenure * tenure
label var lev2var "Variance function at level 2"
// Create variable: predicted values
predict fit, xb
// Create variables: quadratic bound bounds for between-neighborhood variation
gen lev2cbhi = fit + 1.96*(lev2var^.5)
gen lev2cblo = fit - 1.96*(lev2var^.5)
label var fit "Fixed relationship"
label var lev2cbhi "Upper bound between-person variation"
label var lev2cblo "Lower bound between-person variation"
// Graphs: Heterogeneity as a quadratic function
qui twoway (line lev2var tenure, sort) ///
, legend(ring(0)) ///
ytitle("Variance") ///
ylabel(, format(%6.3f)) ///
legend(ring(0) pos(11)) ///
name(variancefunction2, replace)
qui twoway (line fit lev2cbhi lev2cblo tenure, sort) ///
, legend(ring(0) pos(11)) ///
name(betweenperson2, replace)
qui graph combine variancefunction2 betweenperson2, ///
, col(1) xcommon name(quadraticvariance, replace) ///
title("Random coefficient model")
graph combine constantvariance quadraticvariance, row(1) xsize(12) ysize(8)
// Clean up
drop lev2var fit lev2cbhi lev2cblo
