library(lavaan)

setwd(mypath)
dat <- read.table("MPLUS_reliability_exercise_data.csv", sep = ",", header = F)
  
model <- '

TGO =~ NA*V1 + L1*V1 + L2*V2 + L3*V3 + L4*V4 + L5*V5

TGO ~~ 1*TGO

V1 ~~ TH1*V1
V2 ~~ TH2*V2
V3 ~~ TH3*V3
V4 ~~ TH4*V4
V5 ~~ TH5*V5

# model constraints

L1 > 0
L2 > 0
L3 > 0
L4 > 0
L5 > 0


sumL := L1 + L2 + L3 + L4 + L5 
sumTH := TH1 + TH2 + TH3 + TH4 + TH5
omega := (sumL^2)/(sumL^2 + sumTH)

'
mod.fit <- sem(model, dat, se = "bootstrap", bootstrap = 5000)
parameterEstimates(mod.fit, ci = T)
##      lhs op                     rhs label   est    se
## 1    TGO =~                      V1    L1 1.358 0.095
## 2    TGO =~                      V2    L2 1.175 0.092
## 3    TGO =~                      V3    L3 1.116 0.098
## 4    TGO =~                      V4    L4 0.892 0.119
## 5    TGO =~                      V5    L5 1.107 0.125
## 6    TGO ~~                     TGO       1.000 0.000
## 7     V1 ~~                      V1   TH1 0.863 0.131
## 8     V2 ~~                      V2   TH2 1.326 0.160
## 9     V3 ~~                      V3   TH3 1.186 0.125
## 10    V4 ~~                      V4   TH4 2.008 0.259
## 11    V5 ~~                      V5   TH5 2.210 0.221
## 17  sumL :=          L1+L2+L3+L4+L5  sumL 5.649 0.336
## 18 sumTH :=     TH1+TH2+TH3+TH4+TH5 sumTH 7.594 0.391
## 19 omega := (sumL^2)/(sumL^2+sumTH) omega 0.808 0.022
##         z pvalue ci.lower ci.upper
## 1  14.362      0    1.164    1.539
## 2  12.704      0    0.987    1.349
## 3  11.387      0    0.913    1.300
## 4   7.475      0    0.643    1.116
## 5   8.845      0    0.855    1.350
## 6      NA     NA    1.000    1.000
## 7   6.577      0    0.593    1.114
## 8   8.265      0    1.014    1.635
## 9   9.503      0    0.936    1.422
## 10  7.757      0    1.516    2.534
## 11  9.994      0    1.775    2.648
## 17 16.821      0    4.966    6.308
## 18 19.420      0    6.765    8.301
## 19 37.114      0    0.760    0.845

© Copyright 2025 @Yi Feng and @Gregory R. Hancock.