Along the same purpose with tutorials of previous chapters. This tutorial walks through the fitting of multivariate growth modeling in several different frameworks (e.g., multilevel modeling framework, structural equation modeling framework), and demonstrates these models using different R packages (knowing how to fit the models in different packages can be helpful when trying to fit more complex models as each packages as its own advantages and limitations). In this tutorial, we will be using a sample data set that includes repeated measures of individuals’ mathematics achievement through elementary and middle school. Since we want to demonstrate multivariate growth modeling, we also use hyperactivity (repeated measures) a the second depedent variable in the model. Seasonality (spring) will be included as a dynamic predictor to predict math achievement when we move onto time-varying predictors.
The data, code and example provided in this tutorial are from Chapter 8 of Grimm, Ram, and Estabrook (2016), with a few additions in code and commentary. Please refer to the chapter for further interpretations and insights about the analyses.
This tutorial provides line-by-line code to
Fit a bivariate model
Fit a univariate model with a time-varying (dynamic) predictors
Obtain predicted and residual estimates
library(nlme)
library(lavaan)
library(dplyr)
library(tidyr)For this tutorial, we will access the relevant data through the QuantDev website by first setting the file path to the data file, and then using the read.csv() to load these data into R.
# Set filepath
filepath <- "https://quantdev.ssri.psu.edu/sites/qdev/files/nlsy_math_hyp_long.csv"
# Read in the .csv file using the url() function
nlsy_math_hyp_long <- read.csv(file = url(filepath), header=TRUE)
summary(nlsy_math_hyp_long)## id female lb_wght anti_k1
## Min. : 201 Min. :0.0000 Min. :0.00000 Min. :0.000
## 1st Qu.: 260702 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.000
## Median : 497403 Median :0.0000 Median :0.00000 Median :1.000
## Mean : 528449 Mean :0.4935 Mean :0.07969 Mean :1.417
## 3rd Qu.: 792704 3rd Qu.:1.0000 3rd Qu.:0.00000 3rd Qu.:2.000
## Max. :1256601 Max. :1.0000 Max. :1.00000 Max. :8.000
##
## math comp rec bpi
## Min. :12.00 Min. :15.00 Min. :15.00 Min. : 0.00
## 1st Qu.:38.00 1st Qu.:35.00 1st Qu.:38.00 1st Qu.: 20.00
## Median :46.00 Median :45.00 Median :47.50 Median : 60.00
## Mean :46.12 Mean :43.75 Mean :48.87 Mean : 71.09
## 3rd Qu.:54.00 3rd Qu.:51.00 3rd Qu.:59.00 3rd Qu.:110.00
## Max. :81.00 Max. :81.00 Max. :84.00 Max. :280.00
## NA's :22 NA's :7 NA's :105
## as anx hd hyp
## Min. :0.000 Min. :0.000 Min. :0.000 Min. :0.000
## 1st Qu.:0.000 1st Qu.:0.000 1st Qu.:1.000 1st Qu.:0.000
## Median :1.000 Median :1.000 Median :2.000 Median :1.000
## Mean :1.317 Mean :1.611 Mean :2.289 Mean :1.757
## 3rd Qu.:2.000 3rd Qu.:3.000 3rd Qu.:4.000 3rd Qu.:3.000
## Max. :6.000 Max. :5.000 Max. :5.000 Max. :5.000
## NA's :104 NA's :56 NA's :56 NA's :51
## dp wd grade occ
## Min. :0.000 Min. :0.0000 Min. :2.000 Min. :2.000
## 1st Qu.:0.000 1st Qu.:0.0000 1st Qu.:3.000 1st Qu.:2.000
## Median :1.000 Median :0.0000 Median :4.000 Median :3.000
## Mean :1.055 Mean :0.4682 Mean :4.512 Mean :2.841
## 3rd Qu.:2.000 3rd Qu.:1.0000 3rd Qu.:6.000 3rd Qu.:3.000
## Max. :4.000 Max. :3.0000 Max. :8.000 Max. :5.000
## NA's :47 NA's :53
## age men spring anti
## Min. : 82.0 Min. :0.0000 Min. :0.0000 Min. :0.000
## 1st Qu.:109.0 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.000
## Median :126.0 Median :0.0000 Median :1.0000 Median :1.000
## Mean :126.9 Mean :0.1946 Mean :0.6497 Mean :1.577
## 3rd Qu.:143.0 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:2.000
## Max. :175.0 Max. :1.0000 Max. :1.0000 Max. :8.000
## NA's :1147 NA's :51Here, we created a multivariate form of the data in three steps.
First, create a dataset that has id, math score and grade in it, plus 2 dummy variables called d_math and d_hyp, with value of 1 and 0 respectively, and a group variable called grp with value “math”
math <- data.frame(id = nlsy_math_hyp_long$id,
var = nlsy_math_hyp_long$math,
grade = nlsy_math_hyp_long$grade,
d_math = 1,
d_hyp = 0,
grp = 'math')Second, create a dataset that has id, hyp score and grade in it, plus 2 dummy variables called d_math and d_hyp, with value of 0 and 1 respectively, and a group variable called grp with value “hyp”
hyp <- data.frame(id = nlsy_math_hyp_long$id,
var = nlsy_math_hyp_long$hyp,
grade = nlsy_math_hyp_long$grade,
d_math = 0,
d_hyp = 1,
grp = 'hyp')Third, bind the two datasets by rows.
multivariate <- rbind(math, hyp)
# give summary statistics of the dataset
summary(multivariate)## id var grade d_math
## Min. : 201 Min. : 0.0 Min. :2.000 Min. :0.0
## 1st Qu.: 260702 1st Qu.: 2.0 1st Qu.:3.000 1st Qu.:0.0
## Median : 497403 Median :17.0 Median :4.000 Median :0.5
## Mean : 528449 Mean :24.2 Mean :4.512 Mean :0.5
## 3rd Qu.: 792704 3rd Qu.:46.0 3rd Qu.:6.000 3rd Qu.:1.0
## Max. :1256601 Max. :81.0 Max. :8.000 Max. :1.0
## NA's :51
## d_hyp grp
## Min. :0.0 math:2221
## 1st Qu.:0.0 hyp :2221
## Median :0.5
## Mean :0.5
## 3rd Qu.:1.0
## Max. :1.0
## First, subsetting a dataset that has only id, grade and math score, in order to adjust the format into wide (or short later)
nlsy_math_long_math <- data.frame(cbind(nlsy_math_hyp_long$id,
nlsy_math_hyp_long$grade,
nlsy_math_hyp_long$math))
names(nlsy_math_long_math) <- c("id",
"grade",
"math")Similarly, subsetting a dataset that has only id, grade and hyp score, in order to adjust the format into wide (or short later)
nlsy_math_long_hyp <- data.frame(cbind(nlsy_math_hyp_long$id,
nlsy_math_hyp_long$grade,
nlsy_math_hyp_long$hyp))
names(nlsy_math_long_hyp) <- c("id",
"grade",
"hyp")Second, use dplyr to change math and hyp to wide format
math_wide <-nlsy_math_long_math %>% spread(grade, math)
# give name to the variables with the grade
names(math_wide) <- c("id", "math.2","math.3","math.4","math.5","math.6","math.7","math.8")
hyp_wide <-nlsy_math_long_hyp %>% spread(grade, hyp)
# give name to the variables with the grade
names(hyp_wide) <- c("id", "hyp.2","hyp.3","hyp.4","hyp.5","hyp.6","hyp.7","hyp.8")Last, merge the two datasets by id, so each row has the math and hyp scores of the same id.
short <- merge(math_wide, hyp_wide, by ="id")First, subsetting a dataset that has only id, grade and spring, in order to adjust the format into wide (or short later)
nlsy_math_long_spring <- data.frame(cbind(nlsy_math_hyp_long$id,
nlsy_math_hyp_long$grade,
nlsy_math_hyp_long$spring))
names(nlsy_math_long_spring) <- c("id",
"grade",
"spring")Use dplyr to change it to wide format
spring_wide <-nlsy_math_long_spring %>% spread(grade, spring)
names(spring_wide) <- c("id","spring.2","spring.3","spring.4","spring.5","spring.6","spring.7","spring.8")Merge the covariate dataset with the math dataset
math.spring <- merge(math_wide, spring_wide, by ="id")This is to replicate the code from Page 172 (Script 8.2).
lmeCtlList <- lmeControl(maxIter = 500, msMaxIter = 200, tolerance = 1e-4, niter = 100,
msTol = 1e-5, nlmStepMax = 500,
msVerbose = FALSE,
returnObject = TRUE)# bivariate growth model
math.hyp.nlme<-nlme(var~d_math*(b_1i+b_2i*(grade-2))+
d_hyp *(h_1i+h_2i*(grade-2)),
data = multivariate,
fixed = b_1i+b_2i+h_1i+h_2i~1,
random = b_1i+b_2i+h_1i+h_2i~1,
group = ~id,
start = c(35, 4, 1, -1),
weights = varIdent(c(hyp=.3), form = ~1|grp),
na.action = na.omit,
control=lmeCtlList # this was added by Xiao
)
summary(math.hyp.nlme) ## Nonlinear mixed-effects model fit by maximum likelihood
## Model: var ~ d_math * (b_1i + b_2i * (grade - 2)) + d_hyp * (h_1i + h_2i * (grade - 2))
## Data: multivariate
## AIC BIC logLik
## 23612.95 23715.15 -11790.48
##
## Random effects:
## Formula: list(b_1i ~ 1, b_2i ~ 1, h_1i ~ 1, h_2i ~ 1)
## Level: id
## Structure: General positive-definite, Log-Cholesky parametrization
## StdDev Corr
## b_1i 8.04444256 b_1i b_2i h_1i
## b_2i 0.86712080 -0.029
## h_1i 1.24183601 -0.299 0.093
## h_2i 0.06806648 0.201 -0.698 -0.235
## Residual 6.01033637
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | grp
## Parameter estimates:
## math hyp
## 1.000000 0.174825
## Fixed effects: b_1i + b_2i + h_1i + h_2i ~ 1
## Value Std.Error DF t-value p-value
## b_1i 35.25850 0.3551066 3456 99.28990 0e+00
## b_2i 4.34308 0.0873772 3456 49.70495 0e+00
## h_1i 1.90349 0.0581507 3456 32.73373 0e+00
## h_2i -0.05663 0.0142227 3456 -3.98159 1e-04
## Correlation:
## b_1i b_2i h_1i
## b_2i -0.531
## h_1i -0.166 0.037
## h_2i 0.040 -0.067 -0.608
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.32891924 -0.53297652 -0.06387833 0.55392613 3.04894123
##
## Number of Observations: 4391
## Number of Groups: 932## Alternative Approach, use grp variable (cateogoical variable, two values "math" or "hyp") as the predictor
math.hyp.nlme.alt <- nlme(var~(beta_1+d_1i)+(beta_2+d_2i)*(grade-2),
data=multivariate,
fixed=beta_1+beta_2 ~ grp-1,
random=d_1i+d_2i ~ grp-1,
groups = ~id,
start=c(35, 1, 4, -1),
weights = varIdent(c(hyp=.3),form =~ 1|grp),
na.action=na.omit)
summary(math.hyp.nlme.alt)## Nonlinear mixed-effects model fit by maximum likelihood
## Model: var ~ (beta_1 + d_1i) + (beta_2 + d_2i) * (grade - 2)
## Data: multivariate
## AIC BIC logLik
## 23612.95 23715.15 -11790.48
##
## Random effects:
## Formula: list(d_1i ~ grp - 1, d_2i ~ grp - 1)
## Level: id
## Structure: General positive-definite, Log-Cholesky parametrization
## StdDev Corr
## d_1i.grpmath 8.04456494 d_1.grpm d_1.grph d_2.gr
## d_1i.grphyp 1.24206676 -0.299
## d_2i.grpmath 0.86723889 -0.029 0.093
## d_2i.grphyp 0.06857628 0.199 -0.236 -0.691
## Residual 6.01026418
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | grp
## Parameter estimates:
## math hyp
## 1.0000000 0.1748042
## Fixed effects: beta_1 + beta_2 ~ grp - 1
## Value Std.Error DF t-value p-value
## beta_1.grpmath 35.25848 0.3551087 3456 99.28928 0e+00
## beta_1.grphyp 1.90349 0.0581545 3456 32.73155 0e+00
## beta_2.grpmath 4.34308 0.0873783 3456 49.70433 0e+00
## beta_2.grphyp -0.05663 0.0142254 3456 -3.98062 1e-04
## Correlation:
## bt_1.grpm bt_1.grph bt_2.grpm
## beta_1.grphyp -0.166
## beta_2.grpmath -0.531 0.037
## beta_2.grphyp 0.040 -0.608 -0.067
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.32660669 -0.53304174 -0.06379288 0.55389267 3.04921418
##
## Number of Observations: 4391
## Number of Groups: 932This replicated the result on Page 175.
For the diagram for a bivariate linear growth model, see Figure 8.3 (Page 181).
Note: 1. We have slightly more measurements in this model (7 in total), unlike the one showing in the diagram (4 measurements). 2. 1 to 4 in the diagram are named as b0m, b1m, b0h, b1h in the following code, respectively.
model.syntax <- "
# factor loadings
b0m =~ 1* math.2 # this line means b0m is indicated by 1 * math.2, factor loading is 1
b0m =~ 1* math.3
b0m =~ 1* math.4
b0m =~ 1* math.5
b0m =~ 1* math.6
b0m =~ 1* math.7
b0m =~ 1* math.8
b1m =~ 0* math.2
b1m =~ 1* math.3
b1m =~ 2* math.4
b1m =~ 3* math.5
b1m =~ 4* math.6
b1m =~ 5* math.7
b1m =~ 6* math.8 # this line means b1m is indicated by 6 * math.8, factor loading is 6
b0h =~ 1* hyp.2
b0h =~ 1* hyp.3
b0h =~ 1* hyp.4
b0h =~ 1* hyp.5
b0h =~ 1* hyp.6
b0h =~ 1* hyp.7
b0h =~ 1* hyp.8
b1h =~ 0* hyp.2
b1h =~ 1* hyp.3
b1h =~ 2* hyp.4
b1h =~ 3* hyp.5
b1h =~ 4* hyp.6
b1h =~ 5* hyp.7
b1h =~ 6* hyp.8
# factor means
b0m ~ start(35)*1
b1m ~ start(4)*1
b0h ~ 1
b1h ~ 1
# manifest means (fixed at zero, all the variances in the observed variables of math score are explained by the intercept and slope of math, which are b0m and b1m)
math.2 ~ 0
math.3 ~ 0
math.4 ~ 0
math.5 ~ 0
math.6 ~ 0
math.7 ~ 0
math.8 ~ 0
hyp.2 ~ 0
hyp.3 ~ 0
hyp.4 ~ 0
hyp.5 ~ 0
hyp.6 ~ 0
hyp.7 ~ 0
hyp.8 ~ 0
# factor variances (random intercept and random slope for math and hyp scores)
b0m ~~ b0m
b1m ~~ b1m
b0h ~~ b0h
b1h ~~ b1h
# covariances among factors (the following code means covariances between intercepts and slopes for math and hyp scores are freely estimated)
b0m ~~ b1m #variance between intercept and slope of math scores
b0m ~~ b0h #variance between math intercept and hyp intercept
b0m ~~ b1h #variance between math intercept and hyp slope
b1m ~~ b0h #variance between math slope and hyp intercept
b1m ~~ b1h #variance between math slope and hyp slope
b0h ~~ b1h #variance between intercept and slope of hyp scores
# variances of observed variables at each measurement (made equivalent by naming theta1)
math.2 ~~ theta1*math.2
math.3 ~~ theta1*math.3
math.4 ~~ theta1*math.4
math.5 ~~ theta1*math.5
math.6 ~~ theta1*math.6
math.7 ~~ theta1*math.7
math.8 ~~ theta1*math.8
hyp.2 ~~ theta2*hyp.2
hyp.3 ~~ theta2*hyp.3
hyp.4 ~~ theta2*hyp.4
hyp.5 ~~ theta2*hyp.5
hyp.6 ~~ theta2*hyp.6
hyp.7 ~~ theta2*hyp.7
hyp.8 ~~ theta2*hyp.8
# covariances between observed variables, math and hyp, at each measurement
math.2 ~~ theta3*hyp.2
math.3 ~~ theta3*hyp.3
math.4 ~~ theta3*hyp.4
math.5 ~~ theta3*hyp.5
math.6 ~~ theta3*hyp.6
math.7 ~~ theta3*hyp.7
math.8 ~~ theta3*hyp.8
"
fit.bivariate <- sem(model.syntax,
data=short,
missing = "fiml",
estimator = "ml"
)
summary(fit.bivariate)## lavaan (0.5-23.1097) converged normally after 89 iterations
##
## Number of observations 932
##
## Number of missing patterns 96
##
## Estimator ML
## Minimum Function Test Statistic 318.885
## Degrees of freedom 102
## P-value (Chi-square) 0.000
##
## Parameter Estimates:
##
## Information Observed
## Standard Errors Standard
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## b0m =~
## math.2 1.000
## math.3 1.000
## math.4 1.000
## math.5 1.000
## math.6 1.000
## math.7 1.000
## math.8 1.000
## b1m =~
## math.2 0.000
## math.3 1.000
## math.4 2.000
## math.5 3.000
## math.6 4.000
## math.7 5.000
## math.8 6.000
## b0h =~
## hyp.2 1.000
## hyp.3 1.000
## hyp.4 1.000
## hyp.5 1.000
## hyp.6 1.000
## hyp.7 1.000
## hyp.8 1.000
## b1h =~
## hyp.2 0.000
## hyp.3 1.000
## hyp.4 2.000
## hyp.5 3.000
## hyp.6 4.000
## hyp.7 5.000
## hyp.8 6.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## b0m ~~
## b1m -0.204 1.154 -0.176 0.860
## b0h -2.979 0.673 -4.426 0.000
## b1h 0.107 0.164 0.654 0.513
## b1m ~~
## b0h 0.098 0.161 0.608 0.543
## b1h -0.040 0.038 -1.061 0.289
## b0h ~~
## b1h -0.020 0.031 -0.644 0.520
## .math.2 ~~
## .hyp.2 (tht3) -0.011 0.233 -0.046 0.963
## .math.3 ~~
## .hyp.3 (tht3) -0.011 0.233 -0.046 0.963
## .math.4 ~~
## .hyp.4 (tht3) -0.011 0.233 -0.046 0.963
## .math.5 ~~
## .hyp.5 (tht3) -0.011 0.233 -0.046 0.963
## .math.6 ~~
## .hyp.6 (tht3) -0.011 0.233 -0.046 0.963
## .math.7 ~~
## .hyp.7 (tht3) -0.011 0.233 -0.046 0.963
## .math.8 ~~
## .hyp.8 (tht3) -0.011 0.233 -0.046 0.963
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## b0m 35.259 0.356 99.178 0.000
## b1m 4.343 0.088 49.139 0.000
## b0h 1.903 0.058 32.702 0.000
## b1h -0.057 0.014 -3.950 0.000
## .math.2 0.000
## .math.3 0.000
## .math.4 0.000
## .math.5 0.000
## .math.6 0.000
## .math.7 0.000
## .math.8 0.000
## .hyp.2 0.000
## .hyp.3 0.000
## .hyp.4 0.000
## .hyp.5 0.000
## .hyp.6 0.000
## .hyp.7 0.000
## .hyp.8 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## b0m 64.711 5.664 11.424 0.000
## b1m 0.752 0.329 2.282 0.023
## b0h 1.542 0.156 9.908 0.000
## b1h 0.005 0.009 0.539 0.590
## .math.2 (tht1) 36.126 1.863 19.390 0.000
## .math.3 (tht1) 36.126 1.863 19.390 0.000
## .math.4 (tht1) 36.126 1.863 19.390 0.000
## .math.5 (tht1) 36.126 1.863 19.390 0.000
## .math.6 (tht1) 36.126 1.863 19.390 0.000
## .math.7 (tht1) 36.126 1.863 19.390 0.000
## .math.8 (tht1) 36.126 1.863 19.390 0.000
## .hyp.2 (tht2) 1.104 0.057 19.404 0.000
## .hyp.3 (tht2) 1.104 0.057 19.404 0.000
## .hyp.4 (tht2) 1.104 0.057 19.404 0.000
## .hyp.5 (tht2) 1.104 0.057 19.404 0.000
## .hyp.6 (tht2) 1.104 0.057 19.404 0.000
## .hyp.7 (tht2) 1.104 0.057 19.404 0.000
## .hyp.8 (tht2) 1.104 0.057 19.404 0.000Same results with the nlme result.
math.tvc.nlme<-nlme(math~(beta_1+d_1i)+(beta_2+d_2i)*(grade-2)+b3*spring,
data = nlsy_math_hyp_long,
fixed = beta_1+beta_2+b3~1,
random = d_1i+d_2i~1,
group = ~id,
start = c(35, 4, 0),
na.action = na.omit)
summary(math.tvc.nlme)## Nonlinear mixed-effects model fit by maximum likelihood
## Model: math ~ (beta_1 + d_1i) + (beta_2 + d_2i) * (grade - 2) + b3 * spring
## Data: nlsy_math_hyp_long
## AIC BIC logLik
## 15814.31 15854.25 -7900.155
##
## Random effects:
## Formula: list(d_1i ~ 1, d_2i ~ 1)
## Level: id
## Structure: General positive-definite, Log-Cholesky parametrization
## StdDev Corr
## d_1i 7.4670761 d_1i
## d_2i 0.8301936 0.075
## Residual 5.8924537
##
## Fixed effects: beta_1 + beta_2 + b3 ~ 1
## Value Std.Error DF t-value p-value
## beta_1 32.64789 0.4021616 1287 81.18102 0
## beta_2 4.15866 0.0868214 1287 47.89896 0
## b3 4.69404 0.3931911 1287 11.93830 0
## Correlation:
## beta_1 beta_2
## beta_2 -0.334
## b3 -0.543 -0.179
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.265693646 -0.522918388 -0.003805591 0.534813023 2.679014672
##
## Number of Observations: 2221
## Number of Groups: 932This replicated the result on page 178.
model.syntax <- "
# factor loadings
b0m =~ 1* math.2
b0m =~ 1* math.3
b0m =~ 1* math.4
b0m =~ 1* math.5
b0m =~ 1* math.6
b0m =~ 1* math.7
b0m =~ 1* math.8
b1m =~ 0* math.2
b1m =~ 1* math.3
b1m =~ 2* math.4
b1m =~ 3* math.5
b1m =~ 4* math.6
b1m =~ 5* math.7
b1m =~ 6* math.8
# factor means
b0m ~ start(35)*1 # with a starting value at 35
b1m ~ start(4)*1 # with a starting value at 4
# means of the observed variables (fixed at zero)
math.2 ~ 0
math.3 ~ 0
math.4 ~ 0
math.5 ~ 0
math.6 ~ 0
math.7 ~ 0
math.8 ~ 0
# observed variables regressed on the time-varying covariate
math.2 ~ a * spring.2
math.3 ~ a * spring.3
math.4 ~ a * spring.4
math.5 ~ a * spring.5
math.6 ~ a * spring.6
math.7 ~ a * spring.7
math.8 ~ a * spring.8
spring.2 ~ 1
spring.3 ~ 1
spring.4 ~ 1
spring.5 ~ 1
spring.6 ~ 1
spring.7 ~ 1
spring.8 ~ 1
# factor variances
b0m ~~ b0m
b1m ~~ b1m
# covariances between factors
b0m ~~ b1m
# variances of observed variables (made equivalent by naming theta1)
math.2 ~~ theta1*math.2
math.3 ~~ theta1*math.3
math.4 ~~ theta1*math.4
math.5 ~~ theta1*math.5
math.6 ~~ theta1*math.6
math.7 ~~ theta1*math.7
math.8 ~~ theta1*math.8
spring.2 ~~ theta2 * spring.2
spring.3 ~~ theta2 * spring.3
spring.4 ~~ theta2 * spring.4
spring.5 ~~ theta2 * spring.5
spring.6 ~~ theta2 * spring.6
spring.7 ~~ theta2 * spring.7
spring.8 ~~ theta2 * spring.8
"
fit.time.varying <- sem(model.syntax,
data=math.spring,
missing = "fiml",
estimator = "ml"
)
summary(fit.time.varying)## lavaan (0.5-23.1097) converged normally after 81 iterations
##
## Number of observations 932
##
## Number of missing patterns 60
##
## Estimator ML
## Minimum Function Test Statistic 1227.519
## Degrees of freedom 104
## P-value (Chi-square) 0.000
##
## Parameter Estimates:
##
## Information Observed
## Standard Errors Standard
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## b0m =~
## math.2 1.000
## math.3 1.000
## math.4 1.000
## math.5 1.000
## math.6 1.000
## math.7 1.000
## math.8 1.000
## b1m =~
## math.2 0.000
## math.3 1.000
## math.4 2.000
## math.5 3.000
## math.6 4.000
## math.7 5.000
## math.8 6.000
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## math.2 ~
## spring.2 (a) 4.694 0.396 11.864 0.000
## math.3 ~
## spring.3 (a) 4.694 0.396 11.864 0.000
## math.4 ~
## spring.4 (a) 4.694 0.396 11.864 0.000
## math.5 ~
## spring.5 (a) 4.694 0.396 11.864 0.000
## math.6 ~
## spring.6 (a) 4.694 0.396 11.864 0.000
## math.7 ~
## spring.7 (a) 4.694 0.396 11.864 0.000
## math.8 ~
## spring.8 (a) 4.694 0.396 11.864 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## b0m ~~
## b1m 0.468 1.075 0.435 0.663
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## b0m 32.648 0.404 80.888 0.000
## b1m 4.159 0.088 47.461 0.000
## .math.2 0.000
## .math.3 0.000
## .math.4 0.000
## .math.5 0.000
## .math.6 0.000
## .math.7 0.000
## .math.8 0.000
## spring.2 0.609 0.026 23.610 0.000
## spring.3 0.545 0.023 23.978 0.000
## spring.4 0.646 0.024 26.585 0.000
## spring.5 0.661 0.024 27.018 0.000
## spring.6 0.708 0.024 29.605 0.000
## spring.7 0.786 0.036 21.903 0.000
## spring.8 0.718 0.040 18.132 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## b0m 55.757 5.132 10.864 0.000
## b1m 0.689 0.319 2.157 0.031
## .math.2 (tht1) 34.721 1.801 19.282 0.000
## .math.3 (tht1) 34.721 1.801 19.282 0.000
## .math.4 (tht1) 34.721 1.801 19.282 0.000
## .math.5 (tht1) 34.721 1.801 19.282 0.000
## .math.6 (tht1) 34.721 1.801 19.282 0.000
## .math.7 (tht1) 34.721 1.801 19.282 0.000
## .math.8 (tht1) 34.721 1.801 19.282 0.000
## sprng.2 (tht2) 0.223 0.007 33.324 0.000
## sprng.3 (tht2) 0.223 0.007 33.324 0.000
## sprng.4 (tht2) 0.223 0.007 33.324 0.000
## sprng.5 (tht2) 0.223 0.007 33.324 0.000
## sprng.6 (tht2) 0.223 0.007 33.324 0.000
## sprng.7 (tht2) 0.223 0.007 33.324 0.000
## sprng.8 (tht2) 0.223 0.007 33.324 0.000Same results with the nlme result.