*CALCUATES ICCs and Variances of ICCs based on Hedges Derivations program iccvar, rclass syntax, [UNBalance] [Alpha(real 0.05)] version 11.1 tempname cs l4vc l4vcv l3vc l3vcv l2vc l2vcv l1vc tv tv4 icc2 icc3 icc4 ns m N /// se2 se3 se4 v1 v2 v3 v4 b V p r q a b c d e f l u higher lower z cv23 cv24 /// cv34 c23 icc2bar icc3bar icc4bar qp qp2 qp2tv4 qtv4 tempvar n if e(cmd) == "xtmixed" | e(cmd) == "mixed" { local levels = 0 local test1 = e(revars) foreach v in `test1' { if "`v'" == "_cons" { local ++levels } else if "`v'" == "." { display as error "must have at least one group-level random effect" exit 9 } else { display as error "only _cons can have random effects" exit 9 } } capture assert `levels' <= 3 & `levels' > 0 if _rc == 9 { display as error "this command is for models with up to 4 levels only" exit 9 } scalar `cs' = colsof(e(N_g)) assert `cs' == `levels' *Display Headers display _newline as text "Intraclass Correlation Estimates" _newline if `levels' == 1 { if "`unbalance'" == "unbalance" { display as error "there is no balance option for two level models" exit 9 } local l2var = e(ivars) local y = e(depvar) quietly : _diparm lns1_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l2vc' = r(est) scalar `l2vcv' = (r(se))^2 quietly : _diparm lnsig_e, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l1vc' = r(est) scalar `tv' = `l1vc' + `l2vc' scalar `tv4' = ((`tv')^.5)^4 scalar `icc2' = `l2vc' / `tv' scalar `v2' = (((1-`icc2')^2)*`l2vcv')/`tv4' scalar `se2' = (`v2')^.5 matrix define `b' = (`icc2') matrix define `V' = (`v2') matrix colnames `b' = `l2var' matrix rownames `b' = `y' matrix colnames `V' = `l2var' matrix rownames `V' = `l2var' } else if `levels' == 2 { preserve gen `n' = e(sample) local lvars = e(ivars) local k = 3 foreach v in `lvars' { local l`k'var "`v'" local --k } local y = e(depvar) quietly : _diparm lns1_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l3vc' = r(est) scalar `l3vcv' = (r(se))^2 quietly : _diparm lns2_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l2vc' = r(est) scalar `l2vcv' = (r(se))^2 quietly : _diparm lnsig_e, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l1vc' = r(est) scalar `tv' = `l1vc' + `l2vc' + `l3vc' scalar `tv4' = ((`tv')^.5)^4 scalar `icc2' = `l2vc' / `tv' scalar `icc3' = `l3vc' / `tv' if "`unbalance'" == "" { collapse (max) `n', by(`l2var' `l3var') collapse (sum) `n', by(`l3var') quietly : drop if `n' == 0 quietly: means `n' scalar `p' = r(mean_h) display as text "Harmonic Mean of Level 2 Units per Level 3 Unit" _col(50) "=" _col(51) as result %12.3f `p' return scalar p = `p' *variance of icc2 scalar `v2' = (((`p'*((1-`icc2')^2))+(2*`icc2'*(1-`icc2')))*`l2vcv')/(`p'*`tv4') + (((`icc2'^2)*`l3vcv')/`tv4') scalar `se2' = `v2'^.5 *variance of icc3 scalar `v3' = ((((`p'*(`icc3'^2))+(2*`icc3'*(1-`icc3')))*`l2vcv') / (`p'*`tv4'))+((((1-`icc3')^2)*`l3vcv')/`tv4') scalar `se3' = `v3'^.5 *covariance scalar `cv23' = (((((`p'*`icc3'*(1-`icc2'))+(`icc2'*`icc3')+((1-`icc2')*(1-`icc3')))*-1)*`l2vcv')/(`p'*`tv4'))-((`icc2'*(1-`icc3')*`l3vcv')/`tv4') matrix define `b' = (`icc3', `icc2') matrix define `V' = (`v3', `cv23' \ `cv23', `v2') matrix colnames `b' = `l3var' `l2var' matrix rownames `b' = `y' matrix colnames `V' = `l3var' `l2var' matrix rownames `V' = `l3var' `l2var' } else if "`unbalance'" == "unbalance" { *getting covariance between variance components collapse (sum) `n' if e(sample), by(`l2var' `l3var') quietly : gen a = `n'/((`n'*(`l2vc'))+(`l1vc')) quietly : gen b = (`n'^2)/(((`n'*(`l2vc'))+(`l1vc'))^2) quietly : gen p = 1 if `n' > 0 & `n' != . collapse (sum) a b p, by(`l3var') quietly: means p scalar `p' = r(mean_h) gen d = b/(1+(a*(`l3vc'))) gen e = (a^2)/(1+(a*(`l3vc'))) collapse (sum) d e gen c = -1*((d*`l2vcv')/e) scalar `c23' = c display as text "Harmonic Mean of Level 2 Units per Level 3 Unit" _col(50) "=" _col(51) as result %12.3f `p' return scalar p = `p' return scalar c23 = `c23' *variance of icc2 scalar `v2' = ((((1-`icc2')^2)*`l2vcv')/`tv4')+(((`icc2'^2)*`l3vcv')/`tv4')-(((2*`icc2')*(1-`icc2')*`c23')/`tv4') scalar `se2' = `v2'^.5 *variance of icc3 scalar `v3' = (((`icc3'^2)*`l2vcv')/`tv4') + ((((1-`icc3')^2)*`l3vcv')/`tv4') - (((2*`icc3')*(1-`icc3')*`c23')/`tv4') scalar `se3' = `v3'^.5 *covaraince? scalar `cv23' = . matrix define `b' = (`icc3', `icc2') matrix define `V' = (`v3', `cv23' \ `cv23', `v2') matrix colnames `b' = `l3var' `l2var' matrix rownames `b' = `y' matrix colnames `V' = `l3var' `l2var' matrix rownames `V' = `l3var' `l2var' } restore } else if `levels' == 3 { if "`unbalance'" == "unbalance" { display as error "there is no balance option for four level models" exit 9 } preserve gen `n' = e(sample) local lvars = e(ivars) local k = 4 foreach v in `lvars' { local l`k'var "`v'" local --k } local y = e(depvar) quietly : _diparm lns1_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l4vc' = r(est) scalar `l4vcv' = (r(se))^2 quietly : _diparm lns2_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l3vc' = r(est) scalar `l3vcv' = (r(se))^2 quietly : _diparm lns3_1_1, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l2vc' = r(est) scalar `l2vcv' = (r(se))^2 quietly : _diparm lnsig_e, function((exp(@))^2) derivative(2*exp(@)^2) scalar `l1vc' = r(est) scalar `tv' = `l1vc' + `l2vc' + `l3vc' + `l4vc' scalar `tv4' = ((`tv')^.5)^4 scalar `icc2' = `l2vc' / `tv' scalar `icc3' = `l3vc' / `tv' scalar `icc4' = `l4vc' / `tv' scalar `icc4bar' = 1-`icc4' scalar `icc3bar' = 1-`icc3' scalar `icc2bar' = 1-`icc2' collapse (max) `n' if e(sample), by(`l2var' `l3var' `l4var') quietly drop if `n' == 0 collapse (sum) `n', by(`l3var' `l4var') quietly: means `n' scalar `p' = r(mean_h) quietly replace `n' = 1 collapse (sum) `n', by( `l4var') quietly: means `n' scalar `q' = r(mean_h) display as text "Harmonic Mean of Level 2 Units per Level 3 Unit" _col(50) "=" _col(51) as result %12.3f `p' return scalar p = `p' display as text "Harmonic Mean of Level 3 Units per Level 4 Unit" _col(50) "=" _col(51) as result %12.3f `q' return scalar q = `q' scalar `qp' = `q'*`p' scalar `qp2' = `q'*(`p'^2) scalar `qp2tv4' = `q'*(`p'^2)*`tv4' scalar `qtv4' = `q'*`tv4' *variance of icc2 scalar `v2' = ((((`qp2'*(`icc2bar'^2))+(2*`qp'*`icc2'*`icc2bar')+(2*(`icc2'^2)))*`l2vcv')/`qp2tv4') /// + (((`q'-2)*(`icc2'^2)*`l3vcv')/`qtv4') /// + (((`icc2'^2)*`l4vcv')/`tv4') scalar `se2' = `v2'^.5 *variance of icc3 scalar `v3' = ((((`qp2'*(`icc3'^2))+(2*(`qp'-1)*`icc3'*`icc3bar'))*`l2vcv')/`qp2tv4') /// + ((((`q'*(`icc3bar'^2))+(2*`icc3'*`icc3bar'))*`l3vcv')/`qtv4') /// + (((`icc3'^2)*`l4vcv')/`tv4') scalar `se3' = `v3'^.5 *variance of icc4 scalar `v4' = ((((`qp'*(`p'-2)*`icc4'^2)-(2*`icc4'*`icc4bar'))*`l2vcv')/`qp2tv4') /// + ((((`q'*(`icc4'^2)) + (2*`icc4'*`icc4bar'))*`l3vcv')/`qtv4') /// + (((`icc4bar'^2)*`l4vcv')/`tv4') scalar `se4' = `v4'^.5 *covariance between icc2 and icc3 scalar `cv23' = ((((-`qp2'*`icc2bar'*`icc3')+(`qp'*((`icc2'*`icc3')-(`icc2bar'*`icc3bar')))-(`icc2'*`icc3bar')+(`icc2'*`icc3'))*`l2vcv')/`qp2tv4') /// - ((((`q'*`icc2'*`icc3bar')-(`icc2'*`icc3bar')+(`icc2'*`icc3'))*`l3vcv')/`qtv4') /// + ((`icc2'*`icc3'*`l4vcv')/`tv4') *covariance between icc2 and icc4 scalar `cv24' = ((((-1*(`q'*(`p'^2)*`icc2bar'*`icc4'))+(`q'*`p'*((`icc2bar'*`icc4')-(`icc2'*`icc4')))+(`icc2'*`icc4')-(`icc2'*`icc4bar'))*`l2vcv')/(`q'*(`p'^2)*`tv4')) + ((((`q'*`icc2'*`icc4')-(`icc2'*`icc4')+(`icc2'*`icc4bar'))*`l3vcv')/(`q'*`tv4')) - ((`icc2'*`icc4bar'*`l4vcv')/`tv4') *covariance between icc3 and icc4 scalar `cv34' = ((((`qp2'*`icc3'*`icc4')+(`qp'*((`icc3bar'*`icc4')-(`icc3'*`icc4')))+(`icc3'*`icc4')+(`icc3bar'*`icc4bar'))*`l2vcv')/`qp2tv4') /// - ((((`q'*`icc3bar'*`icc4')+(`icc3'*`icc4')+(`icc3bar'*`icc4bar'))*`l3vcv')/`qtv4') /// - ((`icc3'*`icc4bar'*`l4vcv')/`tv4') matrix define `b' = (`icc4', `icc3', `icc2') matrix define `V' = (`v4', `cv34' , `cv24' \ `cv34', `v3' , `cv23' \ `cv24' , `cv23' , `v2') matrix colnames `b' = `l4var' `l3var' `l2var' matrix rownames `b' = `y' matrix colnames `V' = `l4var' `l3var' `l2var' matrix rownames `V' = `l4var' `l3var' `l2var' restore } *post results return local model = e(cmdline) return scalar tv = `tv' return scalar l1vc = `l1vc' return matrix b = `b' return matrix V = `V' *display results & calculate ci scalar `z' = invnormal(1-(`alpha'/2)) local lv = `levels'+1 display _newline as text "{hline 13}" "{c TT}" "{hline 51}" display as text _col(14) "{c |}" _col(19) "ICC" _col(28) "Std. Err." /// _col(43) %12.0f "[" (1-`alpha')*100 "% Conf. Interval]" display as text "{hline 13}" "{c +}" "{hline 51}" forvalues level = `lv'(-1)2 { return scalar l`level'vc = `l`level'vc' return scalar l`level'vc_v = `l`level'vcv' scalar `lower' = `icc`level'' - (`z'*`se`level'') scalar `higher' = `icc`level'' + (`z'*`se`level'') if `lower' < 0 { scalar `lower' = 0 } if `higher' > 1 { scalar `higher' = 1 } display as text as text %12s abbrev("`l`level'var'",12) _col(14) "{c |}" /// _col(15) as result %9.5f `icc`level'' _col(27) as result %9.5f `se`level'' _col(41) /// as result %9.5f `lower' _col(54) as result %9.5f `higher' } display as text "{hline 13}" "{c BT}" "{hline 51}" } else { display as error "can't find last xtmixed or mixed estimations" exit 301 } end