Brand Awareness Study
A major Canadian coffee shop chain is trying to break into the U.S. Market. They assess the following variables twice on a random sample of coffee-drinking adults. The two measurements of each variable are conducted at different times by different interviewers asking somewhat different questions, in such a way that the errors of measurement may be assumed independent. The variables are
> # Brand awareness
>
> rm(list=ls()); options(scipen=999)
> # install.packages("lavaan", dependencies = TRUE) # Only need to do this once
> library(lavaan)
This is lavaan 0.5-23.1097 lavaan is BETA software! Please report any bugs.
>
>
> coffee = read.table("http://www.utstat.toronto.edu/~brunner/openSEM/data/timmy1.data.txt")
> head(coffee)
w1 w2 w3 w4 w5 w6 v1 v2 v3 v4
1 40 23 26 21 48 38 22 22 15 15
2 45 24 29 23 49 48 26 13 8 13
3 29 21 21 13 42 37 18 12 13 13
4 38 26 18 19 47 42 20 9 12 10
5 47 31 30 18 48 52 26 16 22 16
6 31 24 18 13 39 40 20 12 16 18
>
> # Observed variables
> # w1 = Brand Awareness 1
> # w2 = Brand Awareness 2
> # w3 = Ad Awareness 1
> # w4 = Ad Awareness 2
> # w5 = Interest 1
> # w6 = Interest 2
> # v1 = Purchase Intention 1
> # v2 = Purchase Intention 2
> # v3 = Purchase Behaviour 1
> # v4 = Purchase Behaviour 2
> # Latent variables
> # BrAw = True brand awareness
> # AdAw = True advertising awareness
> # Inter = True interest in the product category
> # PI = True purchase intention
> # PBeh = True purchase behaviour
>
> torus1 =
+ '
+ # Latent variable model
+ PI ~ gamma1*BrAw + gamma2*AdAw + gamma3*Inter
+ PBeh ~ gamma4*Inter + beta*PI
+ # Measurement model (simple double measurement)
+ BrAw =~ 1*w1 + 1*w2
+ AdAw =~ 1*w3 + 1*w4
+ Inter =~ 1*w5 + 1*w6
+ PI =~ 1*v1 + 1*v2
+ PBeh =~ 1*v3 + 1*v4
+ # Variances and covariances
+ # Exogenous latent variables
+ BrAw ~~ phi11*BrAw # Var(BrAw) = phi11
+ BrAw ~~ phi12*AdAw # Cov(BrAw,AdAw) = phi12
+ BrAw ~~ phi13*Inter # Cov(BrAw,Inter) = phi13
+ AdAw ~~ phi22*AdAw # Var(AdAw) = phi22
+ AdAw ~~ phi23*Inter # Cov(AdAw,Inter) = phi23
+ Inter ~~ phi33*Inter # Var(Inter) = phi33
+ # Errors in the latent model (epsilons)
+ PI ~~ psi1*PI # Var(epsilon1) = psi1
+ PBeh ~~ psi2*PBeh # Var(epsilon2) = psi2
+ # Measurement errors
+ w1 ~~ omega1*w1 # Var(e1) = omega1
+ w2 ~~ omega2*w2 # Var(e2) = omega2
+ w3 ~~ omega3*w3 # Var(e3) = omega3
+ w4 ~~ omega4*w4 # Var(e4) = omega4
+ w5 ~~ omega5*w5 # Var(e5) = omega5
+ w6 ~~ omega6*w6 # Var(e6) = omega6
+ v1 ~~ omega7*v1 # Var(e7) = omega7
+ v2 ~~ omega8*v2 # Var(e8) = omega8
+ v3 ~~ omega9*v3 # Var(e9) = omega9
+ v4 ~~ omega10*v4 # Var(e10) = omega10
+ # Bounds (Variances are positive)
+ phi11 > 0; phi22 > 0; phi33 > 0
+ psi1 > 0; psi2 > 0
+ omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0
+ omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0
+ ' # End of model torus1
>
> fit1 = lavaan(torus1, data=coffee)
> show(fit1)
lavaan (0.5-23.1097) converged normally after 113 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 77.752 Degrees of freedom 32 P-value (Chi-square) 0.000
It didn't fit. Split the problem up into parts. Look first at the measurement model.
> torus2 = + ' + # Measurement model (still simple double measurement) + BrAw =~ 1*w1 + 1*w2 + AdAw =~ 1*w3 + 1*w4 + Inter =~ 1*w5 + 1*w6 + PI =~ 1*v1 + 1*v2 + PBeh =~ 1*v3 + 1*v4 + # Variances and covariances + # Latent variables + BrAw ~~ phi11*BrAw + BrAw ~~ phi12*AdAw + BrAw ~~ phi13*Inter + BrAw ~~ phi14*PI + BrAw ~~ phi15*PBeh + + AdAw ~~ phi22*AdAw + AdAw ~~ phi23*Inter + AdAw ~~ phi24*PI + AdAw ~~ phi25*PBeh + + Inter ~~ phi33*Inter + Inter ~~ phi34*PI + Inter ~~ phi35*PBeh + + PI ~~ phi44*PI + PI ~~ phi45*PBeh + + PBeh ~~ phi55*PBeh + + # Measurement errors + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0; phi44 > 0; phi55 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus2 > > fit2 = lavaan(torus2, data=coffee) > show(fit2)
lavaan (0.5-23.1097) converged normally after 124 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 76.380 Degrees of freedom 30 P-value (Chi-square) 0.000
The measurement model does not fit. Try a true double measurement model, allowing covariances within sets.
> torus3 = + ' + # Measurement model (double measurement with correlated measurement error) + BrAw =~ 1*w1 + 1*w2 + AdAw =~ 1*w3 + 1*w4 + Inter =~ 1*w5 + 1*w6 + PI =~ 1*v1 + 1*v2 + PBeh =~ 1*v3 + 1*v4 + # Variances and covariances + # Latent variables + BrAw ~~ phi11*BrAw + BrAw ~~ phi12*AdAw + BrAw ~~ phi13*Inter + BrAw ~~ phi14*PI + BrAw ~~ phi15*PBeh + + AdAw ~~ phi22*AdAw + AdAw ~~ phi23*Inter + AdAw ~~ phi24*PI + AdAw ~~ phi25*PBeh + + Inter ~~ phi33*Inter + Inter ~~ phi34*PI + Inter ~~ phi35*PBeh + + PI ~~ phi44*PI + PI ~~ phi45*PBeh + + PBeh ~~ phi55*PBeh + + # Measurement errors + # Variances + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Covariances within sets only (odd and even) + # Odd + w1 ~~ omega13*w3 + w1 ~~ omega15*w5 + w1 ~~ omega17*v1 + w1 ~~ omega19*v3 + w3 ~~ omega35*w5 + w3 ~~ omega37*v1 + w3 ~~ omega39*v3 + w5 ~~ omega57*v1 + w5 ~~ omega59*v3 + v1 ~~ omega79*v3 + # Even + w2 ~~ omega24*w4 + w2 ~~ omega26*w6 + w2 ~~ omega28*v2 + w2 ~~ omega1210*v4 + w4 ~~ omega46*w6 + w4 ~~ omega48*v2 + w4 ~~ omega410*v4 + w6 ~~ omega68*v2 + w6 ~~ omega610*v4 + v2 ~~ omega810*v4 + + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0; phi44 > 0; phi55 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus3 > > fit3 = lavaan(torus3, data=coffee)
Warning message:
In lav_object_post_check(object) :
lavaan WARNING: the covariance matrix of the residuals of the observed
variables (theta) is not positive definite;
use inspect(fit,"theta") to investigate.
Likely one of the measurement error covariance matrices left the parameter space, and its estimate is not positive definite.
> show(fit3)
lavaan (0.5-23.1097) converged normally after 186 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 33.310 Degrees of freedom 10 P-value (Chi-square) 0.000
Still does not fit.
Testing the covariances between error terms with an invalid likelihood ratio test:
> anova(fit2,fit3)
Chi Square Difference Test
Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
fit3 10 10995 11144 33.31
fit2 30 10998 11081 76.38 43.07 20 0.002001 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
This shows why you should not believe a LR test based on a full model that does not fit. Later we will conclude that those covariances could easily be zero.
Consider that the two measurements of each latent variable are different. Maybe they're not really equivalent. Perhaps one in each set (say number 2) should have a coefficient not equal to one. Have to drop the correlated measurement errors for identifiability.
> torus4 = + ' + # Measurement model + BrAw =~ 1*w1 + lambda2*w2 + AdAw =~ 1*w3 + lambda4*w4 + Inter =~ 1*w5 + lambda6*w6 + PI =~ 1*v1 + lambda8*v2 + PBeh =~ 1*v3 + lambda10*v4 + # Variances and covariances + # Latent variables + BrAw ~~ phi11*BrAw + BrAw ~~ phi12*AdAw + BrAw ~~ phi13*Inter + BrAw ~~ phi14*PI + BrAw ~~ phi15*PBeh + + AdAw ~~ phi22*AdAw + AdAw ~~ phi23*Inter + AdAw ~~ phi24*PI + AdAw ~~ phi25*PBeh + + Inter ~~ phi33*Inter + Inter ~~ phi34*PI + Inter ~~ phi35*PBeh + + PI ~~ phi44*PI + PI ~~ phi45*PBeh + + PBeh ~~ phi55*PBeh + + # Measurement errors + # Variances + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0; phi44 > 0; phi55 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus4 > > fit4 = lavaan(torus4, data=coffee) > show(fit4)
lavaan (0.5-23.1097) converged normally after 192 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 17.837 Degrees of freedom 25 P-value (Chi-square) 0.849
The measurement model fits! Now combine it with the latent variable model.
> torus5 = + ' + # Latent variable model + PI ~ gamma1*BrAw + gamma2*AdAw + gamma3*Inter + PBeh ~ gamma4*Inter + beta*PI + # Measurement model + BrAw =~ 1*w1 + lambda2*w2 + AdAw =~ 1*w3 + lambda4*w4 + Inter =~ 1*w5 + lambda6*w6 + PI =~ 1*v1 + lambda8*v2 + PBeh =~ 1*v3 + lambda10*v4 + # Variances and covariances + # Exogenous latent variables + BrAw ~~ phi11*BrAw # Var(BrAw) = phi11 + BrAw ~~ phi12*AdAw # Cov(BrAw,AdAw) = phi12 + BrAw ~~ phi13*Inter # Cov(BrAw,Inter) = phi13 + AdAw ~~ phi22*AdAw # Var(AdAw) = phi22 + AdAw ~~ phi23*Inter # Cov(AdAw,Inter) = phi23 + Inter ~~ phi33*Inter # Var(Inter) = phi33 + # Errors in the latent model (epsilons) + PI ~~ psi1*PI # Var(epsilon1) = psi1 + PBeh ~~ psi2*PBeh # Var(epsilon2) = psi2 + # Measurement errors + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0 + psi1 > 0; psi2 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus5 > > fit5 = lavaan(torus5, data=coffee)
Warning messages:
1: In lav_model_vcov(lavmodel = lavmodel, lavsamplestats = lavsamplestats, :
lavaan WARNING: could not compute standard errors!
lavaan NOTE: this may be a symptom that the model is not identified.
2: In lav_object_post_check(object) :
lavaan WARNING: covariance matrix of latent variables
is not positive definite;
use inspect(fit,"cov.lv") to investigate.
> # Trouble, and there should not be. Follow the suggestion. > inspect(fit5,"cov.lv") # Estimated covariances of latent variables -- tells me nothing.
BrAw AdAw Inter PI PBeh
BrAw 13.442
AdAw 12.525 15.453
Inter 13.805 11.240 14.857
PI 16.127 14.960 15.601 21.214
PBeh 7.043 6.696 6.657 9.530 17.020
> parTable(fit5)
id lhs op rhs user block group free ustart exo label plabel start est se 1 1 PI ~ BrAw 1 1 1 1 NA 0 gamma1 .p1. 0.000 -86.370 NA 2 2 PI ~ AdAw 1 1 1 2 NA 0 gamma2 .p2. 0.000 26.318 NA 3 3 PI ~ Inter 1 1 1 3 NA 0 gamma3 .p3. 0.000 61.392 NA 4 4 PBeh ~ Inter 1 1 1 4 NA 0 gamma4 .p4. 0.000 -0.104 NA 5 5 PBeh ~ PI 1 1 1 5 NA 0 beta .p5. 0.000 0.526 NA 6 6 BrAw =~ w1 1 1 1 0 1 0 .p6. 1.000 1.000 0 7 7 BrAw =~ w2 1 1 1 6 NA 0 lambda2 .p7. 0.418 0.548 NA 8 8 AdAw =~ w3 1 1 1 0 1 0 .p8. 1.000 1.000 0 9 9 AdAw =~ w4 1 1 1 7 NA 0 lambda4 .p9. 0.377 0.548 NA 10 10 Inter =~ w5 1 1 1 0 1 0 .p10. 1.000 1.000 0 11 11 Inter =~ w6 1 1 1 8 NA 0 lambda6 .p11. 0.771 1.088 NA 12 12 PI =~ v1 1 1 1 0 1 0 .p12. 1.000 1.000 0 13 13 PI =~ v2 1 1 1 9 NA 0 lambda8 .p13. 0.508 0.706 NA 14 14 PBeh =~ v3 1 1 1 0 1 0 .p14. 1.000 1.000 0 15 15 PBeh =~ v4 1 1 1 10 NA 0 lambda10 .p15. 0.818 1.037 NA 16 16 BrAw ~~ BrAw 1 1 1 11 NA 0 phi11 .p16. 0.050 13.442 NA 17 17 BrAw ~~ AdAw 1 1 1 12 NA 0 phi12 .p17. 0.000 12.525 NA 18 18 BrAw ~~ Inter 1 1 1 13 NA 0 phi13 .p18. 0.000 13.805 NA 19 19 AdAw ~~ AdAw 1 1 1 14 NA 0 phi22 .p19. 0.050 15.453 NA 20 20 AdAw ~~ Inter 1 1 1 15 NA 0 phi23 .p20. 0.000 11.240 NA 21 21 Inter ~~ Inter 1 1 1 16 NA 0 phi33 .p21. 0.050 14.857 NA 22 22 PI ~~ PI 1 1 1 17 NA 0 psi1 .p22. 0.050 62.631 NA 23 23 PBeh ~~ PBeh 1 1 1 18 NA 0 psi2 .p23. 0.050 12.702 NA 24 24 w1 ~~ w1 1 1 1 19 NA 0 omega1 .p24. 12.120 10.800 NA 25 25 w2 ~~ w2 1 1 1 20 NA 0 omega2 .p25. 9.162 14.293 NA 26 26 w3 ~~ w3 1 1 1 21 NA 0 omega3 .p26. 11.474 7.494 NA 27 27 w4 ~~ w4 1 1 1 22 NA 0 omega4 .p27. 9.046 13.445 NA 28 28 w5 ~~ w5 1 1 1 23 NA 0 omega5 .p28. 10.593 6.330 NA 29 29 w6 ~~ w6 1 1 1 24 NA 0 omega6 .p29. 11.965 6.357 NA 30 30 v1 ~~ v1 1 1 1 25 NA 0 omega7 .p30. 14.725 8.225 NA 31 31 v2 ~~ v2 1 1 1 26 NA 0 omega8 .p31. 10.439 10.300 NA 32 32 v3 ~~ v3 1 1 1 27 NA 0 omega9 .p32. 10.797 4.572 NA 33 33 v4 ~~ v4 1 1 1 28 NA 0 omega10 .p33. 11.085 3.852 NA 34 34 phi11 > 0 1 0 0 0 NA 0 0.000 13.442 0 35 35 phi22 > 0 1 0 0 0 NA 0 0.000 15.453 0 36 36 phi33 > 0 1 0 0 0 NA 0 0.000 14.857 0 37 37 psi1 > 0 1 0 0 0 NA 0 0.000 62.631 0 38 38 psi2 > 0 1 0 0 0 NA 0 0.000 12.702 0 39 39 omega1 > 0 1 0 0 0 NA 0 0.000 10.800 0 40 40 omega2 > 0 1 0 0 0 NA 0 0.000 14.293 0 41 41 omega3 > 0 1 0 0 0 NA 0 0.000 7.494 0 42 42 omega4 > 0 1 0 0 0 NA 0 0.000 13.445 0 43 43 omega5 > 0 1 0 0 0 NA 0 0.000 6.330 0 44 44 omega6 > 0 1 0 0 0 NA 0 0.000 6.357 0 45 45 omega7 > 0 1 0 0 0 NA 0 0.000 8.225 0 46 46 omega8 > 0 1 0 0 0 NA 0 0.000 10.300 0 47 47 omega9 > 0 1 0 0 0 NA 0 0.000 4.572 0 48 48 omega10 > 0 1 0 0 0 NA 0 0.000 3.852 0
It looks like gamma1, 2 3 are too large. Other parameter estimates are not too bad.
Model 4 fit, with no structure in the latent model. Use estimated covariances of latent variables to estimate the regression coefficients. See path diagram above.
> S = inspect(fit4,"cov.lv"); S # Estimated covariance matrix of latent variables
BrAw AdAw Inter PI PBeh
BrAw 19.135
AdAw 12.297 15.914
Inter 13.502 11.306 14.980
PI 16.248 15.070 15.564 21.128
PBeh 7.883 6.144 6.533 9.619 17.155
> S[4,1:3] # Covariances of BrAw, AdAw and Inter with Purchase Intention.
BrAw AdAw Inter
16.24755 15.07001 15.56431
> covxy = as.matrix(S[4,1:3]); covxy
[,1]
BrAw 16.24755
AdAw 15.07001
Inter 15.56431
> Phi = S[1:3,1:3]; Phi
BrAw AdAw Inter
BrAw 19.13517 12.29664 13.50212
AdAw 12.29664 15.91386 11.30579
Inter 13.50212 11.30579 14.98028
> gammahat = solve(Phi) %*% covxy; gammahat
[,1]
BrAw 0.1996463
AdAw 0.3932808
Inter 0.5622264
> torus6 = + ' + # Latent variable model + PI ~ gamma1*BrAw + gamma2*AdAw + gamma3*Inter + + start(0.1996463)*BrAw + start(0.3932808)*AdAw + start(0.5622264)*Inter + PBeh ~ gamma4*Inter + beta*PI + # Measurement model + BrAw =~ 1*w1 + lambda2*w2 + AdAw =~ 1*w3 + lambda4*w4 + Inter =~ 1*w5 + lambda6*w6 + PI =~ 1*v1 + lambda8*v2 + PBeh =~ 1*v3 + lambda10*v4 + # Variances and covariances + # Exogenous latent variables + BrAw ~~ phi11*BrAw # Var(BrAw) = phi11 + BrAw ~~ phi12*AdAw # Cov(BrAw,AdAw) = phi12 + BrAw ~~ phi13*Inter # Cov(BrAw,Inter) = phi13 + AdAw ~~ phi22*AdAw # Var(AdAw) = phi22 + AdAw ~~ phi23*Inter # Cov(AdAw,Inter) = phi23 + Inter ~~ phi33*Inter # Var(Inter) = phi33 + # Errors in the latent model (epsilons) + PI ~~ psi1*PI # Var(epsilon1) = psi1 + PBeh ~~ psi2*PBeh # Var(epsilon2) = psi2 + # Measurement errors + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0 + psi1 > 0; psi2 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus6 > > fit6 = lavaan(torus6, data=coffee) # Still numerical problems.
Warning messages:
1: In lav_model_vcov(lavmodel = lavmodel, lavsamplestats = lavsamplestats, :
lavaan WARNING: could not compute standard errors!
lavaan NOTE: this may be a symptom that the model is not identified.
2: In lav_object_post_check(object) :
lavaan WARNING: covariance matrix of latent variables
is not positive definite;
use inspect(fit,"cov.lv") to investigate.
> parTable(fit6)
id lhs op rhs user block group free ustart exo label plabel start est se 1 1 PI ~ BrAw 1 1 1 1 0.200 0 gamma1 .p1. 0.200 19.023 NA 2 2 PI ~ AdAw 1 1 1 2 0.393 0 gamma2 .p2. 0.393 -61.221 NA 3 3 PI ~ Inter 1 1 1 3 0.562 0 gamma3 .p3. 0.562 31.771 NA 4 4 PBeh ~ Inter 1 1 1 4 NA 0 gamma4 .p4. 0.000 -0.157 NA 5 5 PBeh ~ PI 1 1 1 5 NA 0 beta .p5. 0.000 0.571 NA 6 6 BrAw =~ w1 1 1 1 0 1.000 0 .p6. 1.000 1.000 0 7 7 BrAw =~ w2 1 1 1 6 NA 0 lambda2 .p7. 0.418 0.535 NA 8 8 AdAw =~ w3 1 1 1 0 1.000 0 .p8. 1.000 1.000 0 9 9 AdAw =~ w4 1 1 1 7 NA 0 lambda4 .p9. 0.377 0.552 NA 10 10 Inter =~ w5 1 1 1 0 1.000 0 .p10. 1.000 1.000 0 11 11 Inter =~ w6 1 1 1 8 NA 0 lambda6 .p11. 0.771 1.094 NA 12 12 PI =~ v1 1 1 1 0 1.000 0 .p12. 1.000 1.000 0 13 13 PI =~ v2 1 1 1 9 NA 0 lambda8 .p13. 0.508 0.708 NA 14 14 PBeh =~ v3 1 1 1 0 1.000 0 .p14. 1.000 1.000 0 15 15 PBeh =~ v4 1 1 1 10 NA 0 lambda10 .p15. 0.818 1.034 NA 16 16 BrAw ~~ BrAw 1 1 1 11 NA 0 phi11 .p16. 0.050 18.712 NA 17 17 BrAw ~~ AdAw 1 1 1 12 NA 0 phi12 .p17. 0.000 12.505 NA 18 18 BrAw ~~ Inter 1 1 1 13 NA 0 phi13 .p18. 0.000 13.409 NA 19 19 AdAw ~~ AdAw 1 1 1 14 NA 0 phi22 .p19. 0.050 9.680 NA 20 20 AdAw ~~ Inter 1 1 1 15 NA 0 phi23 .p20. 0.000 11.621 NA 21 21 Inter ~~ Inter 1 1 1 16 NA 0 phi33 .p21. 0.050 14.855 NA 22 22 PI ~~ PI 1 1 1 17 NA 0 psi1 .p22. 0.050 103.141 NA 23 23 PBeh ~~ PBeh 1 1 1 18 NA 0 psi2 .p23. 0.050 12.621 NA 24 24 w1 ~~ w1 1 1 1 19 NA 0 omega1 .p24. 12.120 5.531 NA 25 25 w2 ~~ w2 1 1 1 20 NA 0 omega2 .p25. 9.162 12.959 NA 26 26 w3 ~~ w3 1 1 1 21 NA 0 omega3 .p26. 11.474 13.273 NA 27 27 w4 ~~ w4 1 1 1 22 NA 0 omega4 .p27. 9.046 15.143 NA 28 28 w5 ~~ w5 1 1 1 23 NA 0 omega5 .p28. 10.593 6.334 NA 29 29 w6 ~~ w6 1 1 1 24 NA 0 omega6 .p29. 11.965 6.160 NA 30 30 v1 ~~ v1 1 1 1 25 NA 0 omega7 .p30. 14.725 8.341 NA 31 31 v2 ~~ v2 1 1 1 26 NA 0 omega8 .p31. 10.439 10.303 NA 32 32 v3 ~~ v3 1 1 1 27 NA 0 omega9 .p32. 10.797 4.523 NA 33 33 v4 ~~ v4 1 1 1 28 NA 0 omega10 .p33. 11.085 3.905 NA 34 34 phi11 > 0 1 0 0 0 NA 0 0.000 18.712 0 35 35 phi22 > 0 1 0 0 0 NA 0 0.000 9.680 0 36 36 phi33 > 0 1 0 0 0 NA 0 0.000 14.855 0 37 37 psi1 > 0 1 0 0 0 NA 0 0.000 103.141 0 38 38 psi2 > 0 1 0 0 0 NA 0 0.000 12.621 0 39 39 omega1 > 0 1 0 0 0 NA 0 0.000 5.531 0 40 40 omega2 > 0 1 0 0 0 NA 0 0.000 12.959 0 41 41 omega3 > 0 1 0 0 0 NA 0 0.000 13.273 0 42 42 omega4 > 0 1 0 0 0 NA 0 0.000 15.143 0 43 43 omega5 > 0 1 0 0 0 NA 0 0.000 6.334 0 44 44 omega6 > 0 1 0 0 0 NA 0 0.000 6.160 0 45 45 omega7 > 0 1 0 0 0 NA 0 0.000 8.341 0 46 46 omega8 > 0 1 0 0 0 NA 0 0.000 10.303 0 47 47 omega9 > 0 1 0 0 0 NA 0 0.000 4.523 0 48 48 omega10 > 0 1 0 0 0 NA 0 0.000 3.905 0
> show(fit6)This time, hold those parameters fixed at gamma-hat and search for the rest.
> torus7 = + ' + # Latent variable model + PI ~ 0.1996463*BrAw + 0.3932808*AdAw + 0.5622264*Inter + PBeh ~ gamma4*Inter + beta*PI + # Measurement model + BrAw =~ 1*w1 + lambda2*w2 + AdAw =~ 1*w3 + lambda4*w4 + Inter =~ 1*w5 + lambda6*w6 + PI =~ 1*v1 + lambda8*v2 + PBeh =~ 1*v3 + lambda10*v4 + # Variances and covariances + # Exogenous latent variables + BrAw ~~ phi11*BrAw # Var(BrAw) = phi11 + BrAw ~~ phi12*AdAw # Cov(BrAw,AdAw) = phi12 + BrAw ~~ phi13*Inter # Cov(BrAw,Inter) = phi13 + AdAw ~~ phi22*AdAw # Var(AdAw) = phi22 + AdAw ~~ phi23*Inter # Cov(AdAw,Inter) = phi23 + Inter ~~ phi33*Inter # Var(Inter) = phi33 + # Errors in the latent model (epsilons) + PI ~~ psi1*PI # Var(epsilon1) = psi1 + PBeh ~~ psi2*PBeh # Var(epsilon2) = psi2 + # Measurement errors + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0 + psi1 > 0; psi2 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + ' # End of model torus7 > > fit7 = lavaan(torus7, data=coffee) > show(fit7)
lavaan (0.5-23.1097) converged normally after 139 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 19.014 Degrees of freedom 30 P-value (Chi-square) 0.940
Fits nicely! But it's not quite the MLE.
Start a search at this end point. It's not really too bad.
> coef(fit7); gammahat
gamma4 beta lambda2 lambda4 lambda6 lambda8 lambda10 phi11 phi12 phi13 phi22 phi23
-0.124 0.541 0.525 0.547 1.091 0.707 1.040 19.368 12.319 13.506 15.717 11.295
phi33 psi1 psi2 omega1 omega2 omega3 omega4 omega5 omega6 omega7 omega8 omega9
14.961 3.275 12.638 4.887 12.988 7.218 13.399 6.227 6.112 8.265 10.294 4.613
omega10
3.808
[,1]
BrAw 0.1996463
AdAw 0.3932808
Inter 0.5622264
> torus8 = + ' + # Latent variable model + PI ~ gamma1*BrAw + gamma2*AdAw + gamma3*Inter + + start(0.1996463)*BrAw + start(0.3932808)*AdAw + start(0.5622264)*Inter + PBeh ~ gamma4*Inter + beta*PI + + start(-0.124)*Inter + start(0.541)*PI + # Measurement model + BrAw =~ 1*w1 + lambda2*w2 + start(0.525)*w2 + AdAw =~ 1*w3 + lambda4*w4 + start(0.547)*w4 + Inter =~ 1*w5 + lambda6*w6 + start(1.091)*w6 + PI =~ 1*v1 + lambda8*v2 + start(0.707)*v2 + PBeh =~ 1*v3 + lambda10*v4 + start(1.040)*v4 + # Variances and covariances + # Exogenous latent variables + BrAw ~~ phi11*BrAw + start(19.368)*BrAw # Var(BrAw) = phi11 + BrAw ~~ phi12*AdAw + start(12.319)*AdAw # Cov(BrAw,AdAw) = phi12 + BrAw ~~ phi13*Inter + start(13.506)*Inter # Cov(BrAw,Inter) = phi13 + AdAw ~~ phi22*AdAw + start(15.717)*AdAw # Var(AdAw) = phi22 + AdAw ~~ phi23*Inter + start(11.295)*Inter # Cov(AdAw,Inter) = phi23 + Inter ~~ phi33*Inter + start(14.961)*Inter # Var(Inter) = phi33 + # No more starting values needed. + # Errors in the latent model (epsilons) + PI ~~ psi1*PI # Var(epsilon1) = psi1 + PBeh ~~ psi2*PBeh # Var(epsilon2) = psi2 + # Measurement errors + w1 ~~ omega1*w1 # Var(e1) = omega1 + w2 ~~ omega2*w2 # Var(e2) = omega2 + w3 ~~ omega3*w3 # Var(e3) = omega3 + w4 ~~ omega4*w4 # Var(e4) = omega4 + w5 ~~ omega5*w5 # Var(e5) = omega5 + w6 ~~ omega6*w6 # Var(e6) = omega6 + v1 ~~ omega7*v1 # Var(e7) = omega7 + v2 ~~ omega8*v2 # Var(e8) = omega8 + v3 ~~ omega9*v3 # Var(e9) = omega9 + v4 ~~ omega10*v4 # Var(e10) = omega10 + # Bounds (Variances are positive) + phi11 > 0; phi22 > 0; phi33 > 0 + psi1 > 0; psi2 > 0 + omega1 > 0; omega2 > 0; omega3 > 0; omega4 > 0; omega5 > 0 + omega6 > 0; omega7 > 0; omega8 > 0; omega9 > 0; omega10 > 0 + # The lambda coefficients linking latent variables to their measurements are + # different from zero. They had better be! But which of them are different + # from one? Test "non-linear" constraints. + diff2 := lambda2 - 1 + diff4 := lambda4 - 1 + diff6 := lambda6 - 1 + diff8 := lambda8 - 1 + diff10 := lambda10 - 1 + ' # End of model torus8 > > fit8 = lavaan(torus8, data=coffee) > show(fit8)
lavaan (0.5-23.1097) converged normally after 118 iterations Number of observations 200 Estimator ML Minimum Function Test Statistic 18.962 Degrees of freedom 27 P-value (Chi-square) 0.871
> summary(fit8)
lavaan (0.5-23.1097) converged normally after 118 iterations
Number of observations 200
Estimator ML
Minimum Function Test Statistic 18.962
Degrees of freedom 27
P-value (Chi-square) 0.871
Parameter Estimates:
Information Expected
Standard Errors Standard
Latent Variables:
Estimate Std.Err z-value P(>|z|)
BrAw =~
w1 1.000
w2 (lmb2) 0.528 0.077 6.861 0.000
AdAw =~
w3 1.000
w4 (lmb4) 0.543 0.090 6.013 0.000
Inter =~
w5 1.000
w6 (lmb6) 1.092 0.081 13.528 0.000
PI =~
v1 1.000
v2 (lmb8) 0.707 0.066 10.745 0.000
PBeh =~
v3 1.000
v4 (lm10) 1.040 0.110 9.457 0.000
Regressions:
Estimate Std.Err z-value P(>|z|)
PI ~
BrAw (gmm1) 0.229 0.145 1.581 0.114
AdAw (gmm2) 0.369 0.161 2.285 0.022
Inter (gmm3) 0.553 0.170 3.253 0.001
PBeh ~
Inter (gmm4) -0.129 0.257 -0.502 0.615
PI (beta) 0.546 0.224 2.438 0.015
Covariances:
Estimate Std.Err z-value P(>|z|)
BrAw ~~
AdAw (ph12) 12.301 1.864 6.598 0.000
Inter (ph13) 13.480 1.831 7.360 0.000
AdAw ~~
Inter (ph23) 11.312 1.694 6.679 0.000
Variances:
Estimate Std.Err z-value P(>|z|)
BrAw (ph11) 19.200 3.110 6.174 0.000
AdAw (ph22) 15.910 3.033 5.246 0.000
Inter (ph33) 14.961 2.153 6.949 0.000
.PI (psi1) 3.301 1.340 2.463 0.014
.PBeh (psi2) 12.620 2.097 6.019 0.000
.w1 (omg1) 5.041 2.075 2.430 0.015
.w2 (omg2) 12.974 1.413 9.179 0.000
.w3 (omg3) 7.038 2.218 3.172 0.002
.w4 (omg4) 13.400 1.477 9.074 0.000
.w5 (omg5) 6.224 0.960 6.484 0.000
.w6 (omg6) 6.098 1.063 5.735 0.000
.v1 (omg7) 8.280 1.479 5.598 0.000
.v2 (omg8) 10.299 1.215 8.477 0.000
.v3 (omg9) 4.612 1.682 2.742 0.006
.v4 (om10) 3.810 1.789 2.129 0.033
Defined Parameters:
Estimate Std.Err z-value P(>|z|)
diff2 -0.472 0.077 -6.137 0.000
diff4 -0.457 0.090 -5.058 0.000
diff6 0.092 0.081 1.137 0.256
diff8 -0.293 0.066 -4.454 0.000
diff10 0.040 0.110 0.361 0.718
Constraints:
|Slack|
phi11 - 0 19.200
phi22 - 0 15.910
phi33 - 0 14.961
psi1 - 0 3.301
psi2 - 0 12.620
omega1 - 0 5.041
omega2 - 0 12.974
omega3 - 0 7.038
omega4 - 0 13.400
omega5 - 0 6.224
omega6 - 0 6.098
omega7 - 0 8.280
omega8 - 0 10.299
omega9 - 0 4.612
omega10 - 0 3.810