Linear Growth Model with Time-Invariant Covariates - Multilevel & SEM Implementation in R
Overview
This tutorial illustrates fitting of linear growth models with time-invariant covariates in the multilevel and SEM frameworks in R.
Example data and code are drawn from Chapter 5 of Grimm, Ram, and Estabrook (2017). Specifically, using the NLSY-CYA Dataset we examine how change in children’s mathematics achievement across grade is related to low birth weight and early (kindergarten to first grade) antisocial behaviors. Please see the book chapter for additional interpretations and insights about the analyses.
Preliminaries - Loading libraries used in this script.
library(psych) #for basic functions
library(plyr) #for data management
library(ggplot2) #for plotting
library(nlme) #for mixed effects models
library(lme4) #for mixed effects models
library(lavaan) #for SEM
library(semPlot) #for making SEM diagramsPreliminaries - Data Preparation and Description
For our examples, we use the mathematics achievement scores from the NLSY-CYA Long Data.
Load the repeated measures data
#set filepath for data file
filepath <- "https://raw.githubusercontent.com/LRI-2/Data/main/GrowthModeling/nlsy_math_long_R.dat"
#read in the text data file using the url() function
dat <- read.table(file=url(filepath),
na.strings = ".") #indicates the missing data designator
#copy data with new name
nlsy_math_long <- dat
#Add names the columns of the data set
names(nlsy_math_long) = c('id' , 'female', 'lb_wght',
'anti_k1', 'math' , 'grade' ,
'occ' , 'age' , 'men' ,
'spring' , 'anti')
#view the first few observations in the data set
head(nlsy_math_long, 10)## id female lb_wght anti_k1 math grade occ age men spring anti
## 1 201 1 0 0 38 3 2 111 0 1 0
## 2 201 1 0 0 55 5 3 135 1 1 0
## 3 303 1 0 1 26 2 2 121 0 1 2
## 4 303 1 0 1 33 5 3 145 0 1 2
## 5 2702 0 0 0 56 2 2 100 NA 1 0
## 6 2702 0 0 0 58 4 3 125 NA 1 2
## 7 2702 0 0 0 80 8 4 173 NA 1 2
## 8 4303 1 0 0 41 3 2 115 0 0 1
## 9 4303 1 0 0 58 4 3 135 0 1 2
## 10 5002 0 0 4 46 4 2 117 NA 1 4
Our specific interest is intraindividual change in the repeated
measures of math change across grade, and how
the interindividual differences in those trajectories (intercept and
linear slope) are related to lb_wght and
anti_k1.
As noted in Chapter 2 , it is important to plot the data to obtain a
better understanding of the structure and form of the observed
phenomenon. Here, we want to examine the data to see how the repeated
measures of math are structured with respect to
age.
Longitudinal Plot of Math across Grade at Testing
#intraindividual change trajetories
ggplot(data=nlsy_math_long, #data set
aes(x = grade, y = math, group = id)) + #setting variables
geom_point(size=.5) + #adding points to plot
geom_line(aes(linetype=factor(lb_wght), alpha= anti_k1)) + #adding lines to plot
theme_bw() + #changing style/background
#setting the x-axis with breaks and labels
scale_x_continuous(limits=c(2,8),
breaks = c(2,3,4,5,6,7,8),
name = "Grade at Testing") +
#setting the y-axis with limits breaks and labels
scale_y_continuous(limits=c(10,90),
breaks = c(10,30,50,70,90),
name = "PIAT Mathematics")Note that we are able to get all variables of interest into one plot. Cooll!
Multilevel: Linear Growth Model with Time-Invariant Covariates
To capture the systematic growth in the mathematics scores over
grade, we use the grade variable, so that the rate of
growth indicates (linear) school-related change in mathematics
scores.
Using the notation from the lecture, the model equations are …
\[y_{it} = b_{1i}1 + b_{2i}grade_{it} + u_{it}\]
same as we used previously. However, now we also include person-level (Level-2) predictors for the person-specific intercepts captured by \(b_{1i}\) and linear slopes captured by \(b_{2i}\). Specifically, \[b_{1i} = \beta_{01} + \beta_{11}LowBirthWeight_{i} + \beta_{21}AntiSocial_{i} + d_{1i}\] and \[b_{2i} = \beta_{02} + \beta_{12}LowBirthWeight_{i} + \beta_{22}AntiSocial_{i} + d_{2i}\]
Substitution of one equation into the other we get … \[y_{it} = (\beta_{01} + \beta_{11}LowBirthWeight_{i} + \beta_{21}AntiSocial_{i} + d_{1i}) + (\beta_{02} + \beta_{12}LowBirthWeight_{i} + \beta_{22}AntiSocial_{i} + d_{2i})grade_{it} + u_{it}\]
where \[u_{it} \tilde{} N(0,\sigma^{2}_{u})\], and \[\left[\begin{array} {r} d_{1i} \\ d_{2i} \end{array}\right] \tilde{} N(\left[\begin{array} {r} 0 \\ 0 \end{array}\right] ,\left[\begin{array} {r} \sigma^{2}_{d_1} & \sigma_{d_1d_2} \\ \sigma_{d_1d_2} & \sigma^{2}_{d_2} \end{array}\right]\].
This is then the model that gets coded into the multilevel modeling functions. This can be done in a variety of ways. As usual, we concentrate on use of the nlme:nlme() function.
Linear Growth Model with Time-Invariant Covariates using
nlme() function
This is then the model that gets coded into the multilevel modeling functions. This can be done in a variety of ways. We concentrate on use of the nlme:nlme() function, with an explicit coding style.
Note that in fitting the linear growth model, we recenter the time variable, so that the intercept is the expected score at Grade 2 (the beginning of the data).
#fitting linear growth model with tic and assigning it to an object
lg_math_tic_nlme <- nlme(math ~ (beta_01 + beta_11*lb_wght + beta_21*anti_k1 + d_1i) +
(beta_02 + beta_12*lb_wght + beta_22*anti_k1 +d_2i)*(grade-2),
data=nlsy_math_long,
fixed=beta_01+beta_11+beta_21+beta_02+beta_12+beta_22~1,
random=d_1i+d_2i~1,
group=~id,
start=c(beta_01=30, beta_11=0, beta_21=0,
beta_02=4, beta_12=0, beta_22=0),
na.action=na.exclude)
#obtaining summary of the model using the object we just created
summary(lg_math_tic_nlme)## Nonlinear mixed-effects model fit by maximum likelihood
## Model: math ~ (beta_01 + beta_11 * lb_wght + beta_21 * anti_k1 + d_1i) + (beta_02 + beta_12 * lb_wght + beta_22 * anti_k1 + d_2i) * (grade - 2)
## Data: nlsy_math_long
## AIC BIC logLik
## 15943.14 16000.2 -7961.57
##
## Random effects:
## Formula: list(d_1i ~ 1, d_2i ~ 1)
## Level: id
## Structure: General positive-definite, Log-Cholesky parametrization
## StdDev Corr
## d_1i 7.9413099 d_1i
## d_2i 0.8445068 -0.012
## Residual 6.0213575
##
## Fixed effects: beta_01 + beta_11 + beta_21 + beta_02 + beta_12 + beta_22 ~ 1
## Value Std.Error DF t-value p-value
## beta_01 36.28987 0.4969765 1284 73.02130 0.0000
## beta_11 -2.71621 1.2953408 1284 -2.09691 0.0362
## beta_21 -0.55088 0.2327789 1284 -2.36655 0.0181
## beta_02 4.31516 0.1207520 1284 35.73572 0.0000
## beta_12 0.62464 0.3335739 1284 1.87256 0.0614
## beta_22 -0.01929 0.0589361 1284 -0.32731 0.7435
## Correlation:
## bet_01 bet_11 bet_21 bet_02 bet_12
## beta_11 -0.194
## beta_21 -0.671 -0.025
## beta_02 -0.529 0.096 0.358
## beta_12 0.091 -0.532 0.026 -0.168
## beta_22 0.343 0.026 -0.529 -0.660 -0.055
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.080148313 -0.525141719 -0.008659962 0.530786044 2.534567886
##
## Number of Observations: 2221
## Number of Groups: 932
We get information about how the model fits the data, as well as some information about the relations of the trajectories with the time-invariant covariates. Specifically, we see from the model that there is evidence that differneces in low birth weight status and antisocial behavior are related to individual differences in intercept, but not in linear slope. More detailed interpetation is given in the accompanying video.
Plotting the predicted trajectories.
#obtaining predicted scores for individuals
nlsy_math_long$pred_lg <- predict(lg_math_tic_nlme)
#obtaining predicted scores for prototype
nlsy_math_long$proto_lg <- predict(lg_math_tic_nlme, level=0)
#intraindividual change trajetories
#plotting predicted trajectory
ggplot(data = nlsy_math_long, aes(x = grade, y = pred_lg, group = id)) +
geom_point(color="gray", size=.5) +
geom_line(aes(linetype=factor(lb_wght), alpha= anti_k1)) + #adding lines to plot
theme_bw() + #changing style/background
#setting the x-axis with breaks and labels
scale_x_continuous(limits=c(2,8),
breaks = c(2,3,4,5,6,7,8),
name = "Grade at Testing") +
#setting the y-axis with limits breaks and labels
scale_y_continuous(limits=c(10,90),
breaks = c(10,30,50,70,90),
name = "Predicted PIAT Mathematics")Linear Growth Model with Time-Invariant Covariates using
lmer() function (in the lme4 library)
The model can also be fit in other multilevel package, here
lme4.
#making the recentered time variable
nlsy_math_long$grade_c2 <- nlsy_math_long$grade-2
#fitting linear growth model and assigning it to an object
lg_math_tic_lmer <- lmer(math ~ 1 + lb_wght + anti_k1 +
grade_c2 + I(grade_c2*lb_wght) + I(grade_c2*anti_k1) +
(1 + grade_c2 | id),
data = nlsy_math_long,
REML = TRUE,
na.action = na.exclude)
summary(lg_math_tic_lmer)## Linear mixed model fit by REML ['lmerMod']
## Formula: math ~ 1 + lb_wght + anti_k1 + grade_c2 + I(grade_c2 * lb_wght) +
## I(grade_c2 * anti_k1) + (1 + grade_c2 | id)
## Data: nlsy_math_long
##
## REML criterion at convergence: 15930.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.07031 -0.52511 -0.00901 0.53039 2.53205
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## id (Intercept) 63.4359 7.9647
## grade_c2 0.7344 0.8569 -0.02
## Residual 36.2650 6.0220
## Number of obs: 2221, groups: id, 932
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 36.28899 0.49725 72.979
## lb_wght -2.71513 1.29602 -2.095
## anti_k1 -0.55066 0.23290 -2.364
## grade_c2 4.31611 0.12089 35.703
## I(grade_c2 * lb_wght) 0.62408 0.33389 1.869
## I(grade_c2 * anti_k1) -0.01945 0.05899 -0.330
##
## Correlation of Fixed Effects:
## (Intr) lb_wgh ant_k1 grd_c2 I(_2*l
## lb_wght -0.194
## anti_k1 -0.671 -0.025
## grade_c2 -0.529 0.097 0.358
## I(grd_2*l_) 0.091 -0.532 0.025 -0.168
## I(gr_2*a_1) 0.343 0.026 -0.529 -0.660 -0.055
SEM: Linear Growth Model with Time-Invariant Covariates
Preliminaries - Data Preparation
For the SEM implementation, we use the mathematics achievement scores and time-invariant covariates from the NLSY-CYA Wide Data.
Load the repeated measures data
#set filepath for data file
filepath <- "https://raw.githubusercontent.com/LRI-2/Data/main/GrowthModeling/nlsy_math_wide_R.dat"
#read in the text data file using the url() function
dat <- read.table(file=url(filepath),
na.strings = ".") #indicates the missing data designator
#copy data with new name
nlsy_math_wide <- dat
# Give the variable names
names(nlsy_math_wide)=c('id' ,'female' ,'lb_wght','anti_k1',
'math2' ,'math3' ,'math4' ,'math5' ,'math6' ,'math7' ,'math8' ,
'age2' ,'age3' ,'age4' ,'age5' ,'age6' ,'age7' ,'age8' ,
'men2' ,'men3' ,'men4' ,'men5' ,'men6' ,'men7' ,'men8' ,
'spring2','spring3','spring4','spring5','spring6','spring7','spring8',
'anti2' ,'anti3' ,'anti4' ,'anti5' ,'anti6' ,'anti7' ,'anti8')
#view the first few observations (and columns) in the data set
head(nlsy_math_wide[ ,1:11], 10)## id female lb_wght anti_k1 math2 math3 math4 math5 math6 math7 math8
## 1 201 1 0 0 NA 38 NA 55 NA NA NA
## 2 303 1 0 1 26 NA NA 33 NA NA NA
## 3 2702 0 0 0 56 NA 58 NA NA NA 80
## 4 4303 1 0 0 NA 41 58 NA NA NA NA
## 5 5002 0 0 4 NA NA 46 NA 54 NA 66
## 6 5005 1 0 0 35 NA 50 NA 60 NA 59
## 7 5701 0 0 2 NA 62 61 NA NA NA NA
## 8 6102 0 0 0 NA NA 55 67 NA 81 NA
## 9 6801 1 0 0 NA 54 NA 62 NA 66 NA
## 10 6802 0 0 0 NA 55 NA 66 NA 68 NA
Our specific interest is in change across the the repeated measures
of math, math2 to math8 and how those changes
are related to individual differences in lb_wght and
anti_k1.
Linear Growth Model with Time-Invariant Covariates in
lavaan
Specifying the SEM model
#writing out linear growth model with tic in full SEM way
lg_math_tic_lavaan_model <- '
#latent variable definitions
#intercept
eta1 =~ 1*math2+
1*math3+
1*math4+
1*math5+
1*math6+
1*math7+
1*math8
#linear slope
eta2 =~ 0*math2+
1*math3+
2*math4+
3*math5+
4*math6+
5*math7+
6*math8
#factor variances
eta1 ~~ eta1
eta2 ~~ eta2
#factor covariance
eta1 ~~ eta2
#manifest variances (set equal by naming theta)
math2 ~~ theta*math2
math3 ~~ theta*math3
math4 ~~ theta*math4
math5 ~~ theta*math5
math6 ~~ theta*math6
math7 ~~ theta*math7
math8 ~~ theta*math8
#latent means (freely estimated)
eta1 ~ 1
eta2 ~ 1
#manifest means (fixed to zero)
math2 ~ 0*1
math3 ~ 0*1
math4 ~ 0*1
math5 ~ 0*1
math6 ~ 0*1
math7 ~ 0*1
math8 ~ 0*1
#Time invariant covaraite
#regression of time-invariant covariate on intercept and slope factors
eta1~lb_wght+anti_k1
eta2~lb_wght+anti_k1
#variance of TIV covariates
lb_wght ~~ lb_wght
anti_k1 ~~ anti_k1
#covariance of TIV covaraites
lb_wght ~~ anti_k1
#means of TIV covariates (freely estimated)
lb_wght ~ 1
anti_k1 ~ 1
' #end of model definitionAnd fit the model to the data … using maximum likelihood estimator and full information maximum likelihood to handle missing values.
#estimating the model using sem() function
lg_math_tic_lavaan_fit <- sem(lg_math_tic_lavaan_model,
data = nlsy_math_wide,
meanstructure = TRUE,
estimator = "ML",
missing = "fiml")## Warning in lav_data_full(data = data, group = group, cluster = cluster, : lavaan
## WARNING: some observed variances are (at least) a factor 1000 times larger than
## others; use varTable(fit) to investigate
## Warning in lav_data_full(data = data, group = group, cluster = cluster, : lavaan
## WARNING: due to missing values, some pairwise combinations have less than 10%
## coverage
## Warning in lav_mvnorm_missing_h1_estimate_moments(Y = X[[g]], wt = WT[[g]], : lavaan WARNING:
## Maximum number of iterations reached when computing the sample
## moments using EM; use the em.h1.iter.max= argument to increase the
## number of iterations
We can look more specifically at the numerical summaries of the model.
#obtaining summary of the model using the object we just created
summary(lg_math_tic_lavaan_fit, fit.measures=TRUE)## lavaan 0.6-9 ended normally after 113 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 21
## Number of equality constraints 6
##
## Number of observations 933
## Number of missing patterns 61
##
## Model Test User Model:
##
## Test statistic 220.221
## Degrees of freedom 39
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 892.616
## Degrees of freedom 36
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.788
## Tucker-Lewis Index (TLI) 0.805
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -9785.085
## Loglikelihood unrestricted model (H1) -9674.975
##
## Akaike (AIC) 19600.171
## Bayesian (BIC) 19672.747
## Sample-size adjusted Bayesian (BIC) 19625.108
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.071
## 90 Percent confidence interval - lower 0.062
## 90 Percent confidence interval - upper 0.080
## P-value RMSEA <= 0.05 0.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.097
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Observed
## Observed information based on Hessian
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## eta1 =~
## math2 1.000
## math3 1.000
## math4 1.000
## math5 1.000
## math6 1.000
## math7 1.000
## math8 1.000
## eta2 =~
## math2 0.000
## math3 1.000
## math4 2.000
## math5 3.000
## math6 4.000
## math7 5.000
## math8 6.000
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## eta1 ~
## lb_wght -2.716 1.294 -2.099 0.036
## anti_k1 -0.551 0.232 -2.369 0.018
## eta2 ~
## lb_wght 0.625 0.333 1.873 0.061
## anti_k1 -0.019 0.059 -0.327 0.743
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .eta1 ~~
## .eta2 -0.078 1.145 -0.068 0.945
## lb_wght ~~
## anti_k1 0.007 0.014 0.548 0.584
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .eta1 36.290 0.497 73.052 0.000
## .eta2 4.315 0.122 35.420 0.000
## .math2 0.000
## .math3 0.000
## .math4 0.000
## .math5 0.000
## .math6 0.000
## .math7 0.000
## .math8 0.000
## lb_wght 0.080 0.009 9.031 0.000
## anti_k1 1.454 0.050 29.216 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .eta1 63.064 5.609 11.242 0.000
## .eta2 0.713 0.326 2.185 0.029
## .math2 (thet) 36.257 1.868 19.406 0.000
## .math3 (thet) 36.257 1.868 19.406 0.000
## .math4 (thet) 36.257 1.868 19.406 0.000
## .math5 (thet) 36.257 1.868 19.406 0.000
## .math6 (thet) 36.257 1.868 19.406 0.000
## .math7 (thet) 36.257 1.868 19.406 0.000
## .math8 (thet) 36.257 1.868 19.406 0.000
## lb_wght 0.074 0.003 21.599 0.000
## anti_k1 2.312 0.107 21.599 0.000
#Other summaries
#parameterEstimates(lg_math_tic_lavaan_fit)
#inspect(lg_math_tic_lavaan_fit, what="est")We get information about how the model fits the data, as well as some information about the data itself. Interpetation is covered in more detail in the accompanying video.
We can use the fitted model object to make a diagram (and also check if we specified and estimated the model as intended).
#diagram of fitted model
semPaths(lg_math_tic_lavaan_fit,what = "path", whatLabels = "par")Testing the significance of the Time-Invariant Covariates in
lavaan
To determine whether the addition of low birth weight and kindergarten antisocial behaviors is useful for model fit, we can fit the linear growth model with regression effects set = zero.
Specifying the SEM model with predictor effects at 0.
#writing out linear growth model with tic in full SEM way
lg_math_ticZERO_lavaan_model <- '
#latent variable definitions
#intercept
eta1 =~ 1*math2+
1*math3+
1*math4+
1*math5+
1*math6+
1*math7+
1*math8
#linear slope
eta2 =~ 0*math2+
1*math3+
2*math4+
3*math5+
4*math6+
5*math7+
6*math8
#factor variances
eta1 ~~ eta1
eta2 ~~ eta2
#factor covariance
eta1 ~~ eta2
#manifest variances (set equal by naming theta)
math2 ~~ theta*math2
math3 ~~ theta*math3
math4 ~~ theta*math4
math5 ~~ theta*math5
math6 ~~ theta*math6
math7 ~~ theta*math7
math8 ~~ theta*math8
#latent means (freely estimated)
eta1 ~ 1
eta2 ~ 1
#manifest means (fixed to zero)
math2 ~ 0*1
math3 ~ 0*1
math4 ~ 0*1
math5 ~ 0*1
math6 ~ 0*1
math7 ~ 0*1
math8 ~ 0*1
#Time invariant covaraite
#regression of time-invariant covariate on intercept and slope factors
#FIXED to 0
eta1~0*lb_wght+0*anti_k1
eta2~0*lb_wght+0*anti_k1
#variance of TIV covariates
lb_wght ~~ lb_wght
anti_k1 ~~ anti_k1
#covariance of TIV covaraites
lb_wght ~~ anti_k1
#means of TIV covariates (freely estimated)
lb_wght ~ 1
anti_k1 ~ 1
' #end of model definitionAnd fit the model to the data … using maximum likelihood estimator and full information maximum likelihood to handle missing values.
#estimating the model using sem() function
lg_math_ticZERO_lavaan_fit <- sem(lg_math_ticZERO_lavaan_model,
data = nlsy_math_wide,
meanstructure = TRUE,
estimator = "ML",
missing = "fiml")## Warning in lav_data_full(data = data, group = group, cluster = cluster, : lavaan
## WARNING: some observed variances are (at least) a factor 1000 times larger than
## others; use varTable(fit) to investigate
## Warning in lav_data_full(data = data, group = group, cluster = cluster, : lavaan
## WARNING: due to missing values, some pairwise combinations have less than 10%
## coverage
## Warning in lav_mvnorm_missing_h1_estimate_moments(Y = X[[g]], wt = WT[[g]], : lavaan WARNING:
## Maximum number of iterations reached when computing the sample
## moments using EM; use the em.h1.iter.max= argument to increase the
## number of iterations
We can look more specifically at the numerical summaries of the model.
#obtaining summary of the model using the object we just created
summary(lg_math_ticZERO_lavaan_fit, fit.measures=TRUE)## lavaan 0.6-9 ended normally after 90 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 17
## Number of equality constraints 6
##
## Number of observations 933
## Number of missing patterns 61
##
## Model Test User Model:
##
## Test statistic 234.467
## Degrees of freedom 43
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 892.616
## Degrees of freedom 36
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.776
## Tucker-Lewis Index (TLI) 0.813
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -9792.208
## Loglikelihood unrestricted model (H1) -9674.975
##
## Akaike (AIC) 19606.416
## Bayesian (BIC) 19659.638
## Sample-size adjusted Bayesian (BIC) 19624.703
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.069
## 90 Percent confidence interval - lower 0.061
## 90 Percent confidence interval - upper 0.078
## P-value RMSEA <= 0.05 0.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.103
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Observed
## Observed information based on Hessian
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## eta1 =~
## math2 1.000
## math3 1.000
## math4 1.000
## math5 1.000
## math6 1.000
## math7 1.000
## math8 1.000
## eta2 =~
## math2 0.000
## math3 1.000
## math4 2.000
## math5 3.000
## math6 4.000
## math7 5.000
## math8 6.000
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## eta1 ~
## lb_wght 0.000
## anti_k1 0.000
## eta2 ~
## lb_wght 0.000
## anti_k1 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .eta1 ~~
## .eta2 -0.181 1.150 -0.158 0.875
## lb_wght ~~
## anti_k1 0.007 0.014 0.548 0.584
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .eta1 35.267 0.355 99.229 0.000
## .eta2 4.339 0.088 49.136 0.000
## .math2 0.000
## .math3 0.000
## .math4 0.000
## .math5 0.000
## .math6 0.000
## .math7 0.000
## .math8 0.000
## lb_wght 0.080 0.009 9.031 0.000
## anti_k1 1.454 0.050 29.216 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .eta1 64.562 5.659 11.408 0.000
## .eta2 0.733 0.327 2.238 0.025
## .math2 (thet) 36.230 1.867 19.410 0.000
## .math3 (thet) 36.230 1.867 19.410 0.000
## .math4 (thet) 36.230 1.867 19.410 0.000
## .math5 (thet) 36.230 1.867 19.410 0.000
## .math6 (thet) 36.230 1.867 19.410 0.000
## .math7 (thet) 36.230 1.867 19.410 0.000
## .math8 (thet) 36.230 1.867 19.410 0.000
## lb_wght 0.074 0.003 21.599 0.000
## anti_k1 2.312 0.107 21.599 0.000
#Other summaries
#parameterEstimates(lg_math_tic_lavaan_fit)
#inspect(lg_math_tic_lavaan_fit, what="est")And can the new constrained model as a diagram (and also check if we specified and estimated the model as intended).
#diagram of fitted model
semPaths(lg_math_ticZERO_lavaan_fit,what = "path", whatLabels = "par")The chi-square difference test (obtained from the output or using the
anova() function) indicates that the time invariant
predictors are useful to the model. For a more in-depth discussion of
model fit and comparison, see pages 110-111 in the book.
anova(lg_math_tic_lavaan_fit,lg_math_ticZERO_lavaan_fit)## Chi-Squared Difference Test
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## lg_math_tic_lavaan_fit 39 19600 19673 220.22
## lg_math_ticZERO_lavaan_fit 43 19606 19660 234.47 14.245 4 0.006552
##
## lg_math_tic_lavaan_fit
## lg_math_ticZERO_lavaan_fit **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Conclusion
This tutorial has presented how time invariant covariates can be added to the linear growth model in both the multilevel and SEM frameworks.
Examining interindividual differences in intraindividual change! Bring it!