* Prevent "more" messages from appearing
set more off
* Control line length
set linesize 150

********************************************************************************
******* BEGIN DATA MANIPULATION FOR CHAPTER 12 REPEATED MEASURES EXAMPLE *******
*******               CHANGE "filesave" to your directory                *******
********************************************************************************

* Defining global variable for file location to be replaced in code below
global filesave "C:\Dropbox\PilesOfVariance\Chapter12\STATA"

* Import chapter 12 stacked data and create variables for analysis
use "$filesave\STATA_Chapter12.dta", clear
* Log of response time
gen logrt = log(rt)
* Make copy of 'older' and dummy codes to be used in heterogeneous variance models
gen agegroup = older
gen isyounger=0
replace isyounger=1 if (older==0)
gen isolder=1
replace isolder=1 if (older==1)
* Piecewise effect of age
gen yrs65=0
replace yrs65=age-65 if (older==1)
* Item predictors
gen rel3 = relevance - 3
gen sal3 = salience - 3
* Item salience random slope for heterogeneous random effects by age
gen isysal = isyounger*sal3
gen isosal = isolder*sal3
label variable logrt "logrt: Natural Log RT in Seconds"
label variable agegroup "agegroup: Younger=0, Older=1"
label variable isyounger "isyounger: 0=Older, 1=Younger"
label variable isolder "isolder: 0=Younger, 1=Older"
label variable yrs65 "yrs65: Subject Years of Age Over 65 (0=65)"
label variable rel3  "rel3: Item Relevance (0=3)"
label variable sal3  "sal3: Item Salience (0=3)"
label variable isysal "isysal: Younger Random Salience Slope"
label variable isosal "isosal: Older Random Salience Slope"

* Subset sample to complete cases for all predictors
egen nummiss = rowmiss(older yrs65 relevance salience logrt)
drop if nummiss>0

********************************************************************************
*******          BEGIN CHAPTER 12 CROSSED RANDOM EFFECTS MODELS          *******
*******          NOTE: MODELS WITH RANDOM SLOPES WILL NOT CONVERGE       *******
*******           BUT ARE PROVIDED BELOW FOR ILLUSTRATION                *******
********************************************************************************

* Save results to separate file
log using $filesave\STATA_Chapter12_Output, replace name(STATA_Chapter12)

display as result "Chapter 12: Descriptive Statistics for Subject-Level Variables"
preserve 
collapse older age, by(subjectid)
tabulate older, summarize(age)
restore 

display as result "Chapter 12: Descriptive Statistics for Item-Level Variables"
preserve 
collapse relevance salience, by(itemid)
summarize relevance salience
preserve 

display as result "Chapter 12: Descriptive Statistics for Trial-Level Variables"
summarize rt logrt

display as result "Eq 12.8: Empty Means, Single-Level Model for Log RT"
mixed logrt , ///
              || subjectid: , noconstant variance reml covariance(unstructured),
      estat ic, n(1),
      estimates store FitEmpty1,

display as result "Eq 12.9: Add Subject Random Intercept Variance"
mixed logrt , ///
              || subjectid: , variance reml covariance(unstructured),
      estat ic, n(1),
      estat icc,
      estimates store FitEmpty2,

display as result "Eq 12.10: Add Item Random Intercept Variance"
mixed logrt , ///
              || _all: R.subjectid , variance reml ///
              || _all: R.itemid ,
      estat ic, n(1),
      estimates store FitEmpty3,
      lrtest FitEmpty3 FitEmpty2,

display as result "Eq 12.11: Add Fixed Effects of Item Predictors"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3, ///
              || _all: R.subjectid , variance reml  ///
              || _all: R.itemid ,
      estat ic, n(1),
      * Multivariate Test of 3 Item Predictor Effects
      test (c.rel3=0) (c.sal3=0) (c.rel3#c.sal3=0)
      estimates store FitItem,

display as result "Ch 12: Remove Item Random Intercept Variance"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3, ///
              || _all: R.subjectid , variance reml,
      estat ic, n(1),
      estimates store FitNoRandItem,
      lrtest FitItem FitNoRandItem,

display as result "Eq 12.12: Add Fixed Main Effects of Subject Predictors"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65, ///
              || _all: R.subjectid , variance reml ///
              || _all: R.itemid ,
      estat ic, n(1),
      * Multivariate Test of 2 Subject Predictor Effects
      test (c.older=0) (c.yrs65=0)
      estimates store FitSubject,

display as result "Ch 12: Remove Subject Random Intercept Variance"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65, ///
              || _all: R.itemid , variance reml,
      estat ic, n(1),
      estimates store FitNoRandSubject,
      lrtest FitSubject FitNoRandSubject,

****** END CHAPTER 12 MODELS ******

* Close log
log close STATA_Chapter12
* Convert log to html using custom downloaded package
log2html $filesave\STATA_Chapter12_Output, replace


****** CHAPTER 12 MODELS WITH RANDOM SLOPES (WILL NOT CONVERGE) ******

display as result "Eq 12.13: Add Subject Random Salience Slope Variance and Covariance"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65, ///
              || subjectid: sal3, variance reml covariance(unstructured) ///
              || _all: R.itemid ,
      estat ic, n(1),
      estimates store FitRandSal,
      lrtest FitRandSal FitSubject,

display as result "Ch 12: Add Subject Random Relevance Slope Variance and Covariances"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65, ///
              || subjectid: sal3 rel3, variance reml covariance(unstructured) ///
              || _all: R.itemid ,
      estat ic, n(1),
      estimates store FitRandMean,
      lrtest FitRandMean FitRandSal,

display as result "Eq 12.14: Add All Possible Fixed Effect Interactions"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65 ///
            c.rel3#c.older c.sal3#c.older c.rel3#c.sal3#c.older ///
            c.rel3#c.yrs65 c.sal3#c.yrs65 c.rel3#c.sal3#c.yrs65, /// 
              || subjectid: sal3, variance reml covariance(unstructured) ///
              || _all: R.itemid ,
      estat ic, n(1),

display as result "Eq 12.15: Keep Significant Fixed Effect Interactions"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65 ///
            c.rel3#c.older c.sal3#c.older c.rel3#c.sal3#c.older, ///
              || subjectid: sal3, variance reml covariance(unstructured) ///
              || _all: R.itemid ,
      estat ic, n(1),
      estimates store FitSigInt,

display as result "Ch 12: Remove Subject Random Slope Variance and Covariance"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65 ///
            c.rel3#c.older c.sal3#c.older c.rel3#c.sal3#c.older, ///
              || subjectid: , variance reml covariance(unstructured) ///
              || _all: R.itemid ,
      estat ic, n(1),
      estimates store FitNoRandSal,
      lrtest FitSigInt FitNoRandSal,

display as result "Ch 12: Add Separate Residual Variances by Age Group"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65 ///
            c.rel3#c.older c.sal3#c.older c.rel3#c.sal3#c.older, ///
              || subjectid: sal3, variance reml covariance(unstructured) ///
              || _all: R.itemid ,
              residuals(independent,by(agegroup)),
      estat ic, n(1),
      estimates store FitHetResVar,
      lrtest FitHetResVar FitSigInt,

display as result "Ch 12: Add Separate Subject Random Effects Variances and Covariance by Age Group"
mixed logrt c.rel3 c.sal3 c.rel3#c.sal3 c.older c.yrs65 ///
            c.rel3#c.older c.sal3#c.older c.rel3#c.sal3#c.older, ///
              || subjectid: isyounger isysal, noconstant variance reml covariance(unstructured) ///
              || subjectid: isolder isosal, noconstant  covariance(unstructured)
              || _all: R.itemid ,
              residuals(independent,by(agegroup)),
      estat ic, n(1),
      * Relevance Effect at Salience=3 for Younger
      lincom c.rel3*1 + c.rel3#c.older*0 + c.rel3#c.sal3*0 + c.rel3#c.sal3#c.older*0
      * Relevance Effect at Salience=3 for Older
      lincom c.rel3*1 + c.rel3#c.older*1 + c.rel3#c.sal3*0 + c.rel3#c.sal3#c.older*0
      * Salience Effect at Relevance=3 for Younger
      lincom c.sal3*1 + c.sal3#c.older*0 + c.rel3#c.sal3*0 + c.rel3#c.sal3#c.older*0
      * Salience Effect at Relevance=3 for Older
      lincom c.sal3*1 + c.sal3#c.older*1 + c.rel3#c.sal3*0 + c.rel3#c.sal3#c.older*0
      * Relevance by Salience for Younger
      lincom c.rel3#c.sal3*1 + c.rel3#c.sal3#c.older*0
      * Relevance by Salience for Older
      lincom c.rel3#c.sal3*1 + c.rel3#c.sal3#c.older*1
      estimates store FitHetRandVar,
      lrtest FitHetRandVar FitHetResVar,
      predict PredFinal, xb,
      margins, at (c.older=(0(1)1) c.yrs65=(0(15)15) c.older=(-2(1)2) c.older=(-2(3)1)) vsquish,
corr logrt PredFinal