UPDATE: This problem was (embarrassingly) solved by specifying the intercept in the regression equation as shown below:
lcs ~ 1 + 0*Y1_bl_ctr #gamma is set to 0 for equivilance with the t-test
This question is distinct enough from (Change Score Model in lavaan) that I am migrating it here, but I am including the link for reference purposes.
At the aforementioned post, I was attempting to calculate a latent change score for two waves of observation. As indicated there, I believe the model I specified was over-identified.
Having thought about this for a day and reading a little more, I think I may have been over-complicating the problem.
In a recent paper by Colman and colleagues (http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3794455/), it was shown how under certain conditions a latent change score model is equivalent to a simple t-test.
The authors specify a path model in the paper as follows:
As noted in the path diagram, the authors state that the model is equivalent to a paired sample t-test when $\gamma=0$.
The authors additionally make the following constraints on the model:
- No $Y_2$ residual error;
- Intercept of $Y_2$ is set to zero;
- The auto-regressive path is set to 1;
- The $LCS$ to $Y_2$ path is set to 1.
NOTE: Within the body of the paper, the authors also specify that they are baseline-centering the $Y_1$ and $Y_2$ values.
The authors supply the data for their paper at (http://trippcenter.uchc.edu/modeling). For convenience, I provide the dput
output of the data in the following code chunk and assign it to an object named dat
.
dat <- structure(list(Y1 = c(7.6, 8.4, 8.4, 8.7, 9.7, 10.9, 9.7, 7.4,
7.9, 9.3, 10.3, 11.3, 14, 8.6, 13.3, 10.9, 7.9, 9.5, 9.9, 11.5,
11.8, 9.9, 8.8, 10.6, 11.3, 7.6, 8, 8.1, 8.6, 8.2, 7.4, 8.5,
9.3, 9.4, 10.6, 10.3, 10.8, 8.2, 10.1, 7.8, 7.2, 7.8, 8.1, 8.1,
8.2, 9.1, 8.7, 8.6, 7.6, 8.8, 10.3, 11.3, 7.7, 7.8, 8.5, 8.4,
8.3, 9.2, 9.9, 9.8, 10.9, 7.8, 10.1, 9.5, 8, 8.3, 8.2, 7.8, 8.8,
9.4, 8.7, 8.8, 10.8, 11.5, 7.6, 7.6, 7.9, 8.2, 9.7, 10, 8.8,
10, 12.7, 7.5, 13.4, 8.3, 13.8, 14, 8.4, 14, 10, 9.5, 6.2, 11.7,
9.7, 14, 7.3), Y2 = c(7.9, 7.8, 5.9, 7.5, 7.5, 9.1, 8, 7.9, 6.7,
8.1, 8.5, 12, 14, 6.9, 10, 10.9, 7.2, 8.9, 10.8, 7.1, 7.4, 9.7,
10.3, 9.3, 11.6, 8, 7.1, 6.7, 8.2, 9.3, 7.6, 9.9, 8.9, 8.8, 7.2,
10.1, 6.7, 6.2, 8.9, 7.3, 7.6, 7.5, 7.3, 9.6, 8.1, 7.8, 8.7,
8.4, 11.4, 9, 10.2, 12.5, 7, 8.7, 8, 7.2, 8.9, 10.4, 9.4, 10.8,
9.9, 6.3, 5.7, 10.1, 7.8, 8.2, 7.4, 7.7, 11.8, 7.1, 6.8, 8.1,
9.2, 10.2, 8.4, 7.1, 9, 6.9, 8.7, 8.8, 9.3, 8.6, 8.5, 7.7, 13.8,
8.7, 10.5, 14, 10.1, 14, 14, 9, 14, 14, 10.3, 14, 8.4)), .Names = c("Y1",
"Y2"), class = "data.frame", row.names = c(NA, -97L))
In laavan
I am trying to implement the model as follows:
# baseline mean center Y1 and Y2
dat$Y1_bl_ctr = dat$Y1 - mean(dat$Y1)
dat$Y2_bl_ctr = dat$Y2 - mean(dat$Y1)
test <- '
# measurement model
lcs =~ 1*Y2_bl_ctr #4 - the LCS to Y2 path is set to 1
# regressions
lcs ~ 0*Y1_bl_ctr #gamma is set to 0 for equivilance with the t-test
Y2_bl_ctr ~ 1*Y1_bl_ctr #3 - The auto-regressive path is set to 1
# residual error
lcs ~~ 1*lcs
Y2_bl_ctr ~~ 0*Y2_bl_ctr #1 - No Y2 residual error
'
summary(test <- lavaan(test
,data=dat
,int.lv.free = TRUE #intercepts of LCS is to be estimated
,int.ov.free = FALSE #2- Intercept of Y2 is set to zero;
)
)
As is probably obvious, this overly-restricted specification results in an unestimated model.
#Error in lav_syntax_parse_rhs(rhs = rhs.formula[[2L]], op = op) :
# lavaan ERROR: I'm confused parsing this line: offsetY2_bl_ctr
If I respecify the measurement model as lcs =~ Y2_bl_ctr
, I get the following output from laavan
.
#lavaan (0.5-20) converged normally after 7 iterations
#
# Number of observations 97
#
# Estimator ML
# Minimum Function Test Statistic 11.543
# Degrees of freedom 1
# P-value (Chi-square) 0.001
#
#Parameter Estimates:
#
# Information Expected
# Standard Errors Standard
#
#Latent Variables:
# Estimate Std.Err Z-value P(>|z|)
# lcs =~
# Y2_bl_ctr 1.780 0.128 13.928 0.000
#
#Regressions:
# Estimate Std.Err Z-value P(>|z|)
# lcs ~
# Y1_bl_ctr 0.000
# Y2_bl_ctr ~
# Y1_bl_ctr 1.000
#
#Variances:
# Estimate Std.Err Z-value P(>|z|)
# lcs 1.000
# Y2_bl_ctr 0.000
However, I still don't appear to get an intercept estimated - just the $Y_2$ to $LCS$ path coefficient.
t.test(dat$Y1,dat$Y2,paired=TRUE)
# Paired t-test
#
#data: dat$Y1 and dat$Y2
#t = 2.1734, df = 96, p-value = 0.03221
#alternative hypothesis: true difference in means is not equal to 0
#95 percent confidence interval:
# 0.03423662 0.75545410
#sample estimates:
#mean of the differences
# 0.3948454
The 0.395 value is the correct value based on the paper.
Any thoughts on how this model can be specified in laavan
in order to produce equivalent results to the t-test?