devtools::install_github("simsem/semTools/semTools")
library(semTools)
?measEq.syntax
I can't seem to be able to install semTools from GitHub using devtools
library(devtools)
--
You received this message because you are subscribed to the Google Groups "lavaan" group.
To unsubscribe from this group and stop receiving emails from it, send an email to lavaan+un...@googlegroups.com.
To post to this group, send email to lav...@googlegroups.com.
Visit this group at https://groups.google.com/group/lavaan.
For more options, visit https://groups.google.com/d/optout.
I’m trying to use the new measEq.Syntax() function to study measurement invariance on the ISSP religion data. Here is the structure of the data frame:
> str(df)
'data.frame': 35198 obs. of 20 variables:
$ V28 : Ord.factor w/ 3 levels "I don't believe in God"<..: 3 3 2 2 2 2 2 2 2 2 ...
$ V30 : Ord.factor w/ 4 levels "No, definitely not"<..: 2 4 2 2 4 NA 1 NA 2 1 ...
$ V31 : Ord.factor w/ 4 levels "No, definitely not"<..: 2 4 2 1 NA 1 1 NA 1 3 ...
$ V32 : Ord.factor w/ 4 levels "No, definitely not"<..: 2 2 2 1 NA 1 1 1 2 3 ...
$ V33 : Ord.factor w/ 4 levels "No, definitely not"<..: 3 4 2 1 1 3 1 2 2 3 ...
$ V35 : Ord.factor w/ 5 levels "Strongly disagree"<..: 2 5 3 1 2 NA 1 3 3 NA ...
$ V37 : Ord.factor w/ 5 levels "Strongly disagree"<..: 3 2 3 1 1 1 1 3 3 NA ...
$ V46 : Ord.factor w/ 4 levels "Never"<"Yearly"<..: 4 3 3 1 2 2 3 2 1 3 ...
$ V47 : Ord.factor w/ 4 levels "Never"<"Yearly"<..: 3 NA 2 1 2 NA 1 1 1 3 ...
$ V48 : Ord.factor w/ 4 levels "Never"<"Yearly"<..: 4 4 2 1 2 4 3 1 2 3 ...
$ V49 : Ord.factor w/ 5 levels "Never"<"Yearly"<..: 4 1 2 1 1 1 1 1 1 3 ...
$ V50 : Ord.factor w/ 5 levels "Never"<"Yearly"<..: 2 1 1 2 1 1 1 1 1 2 ...
$ V51 : Ord.factor w/ 5 levels "Highly non-religious"<..: 4 1 2 1 1 4 2 2 2 2 ...
$ ATTEND : Ord.factor w/ 4 levels "Never"<"Yearly"<..: 4 1 2 1 1 2 1 1 1 3 ...
$ AGE : Ord.factor w/ 5 levels "0-24"<"25-44"<..: 4 2 2 2 2 3 3 3 3 3 ...
$ SEX : Factor w/ 2 levels "Male","Female": 2 2 1 1 1 1 2 1 1 2 ...
$ DEGREE : Ord.factor w/ 6 levels "No formal qualification"<..: 2 4 4 6 4 6 3 3 4 3 ...
$ RELIGGRP : Factor w/ 5 levels "No religion",..: 2 2 4 1 1 1 2 1 1 2 ...
$ COUNTRY.NAME: Factor w/ 26 levels "AU-Australia",..: 2 3 4 5 7 8 9 12 13 17 ...
$ YEAR : Ord.factor w/ 1 level "2008": 1 1 1 1 1 1 1 1 1 1 ...
>
>
I defined the following items as ordinal and converted to numeric, and for simplicity I used “SEX” as the group factor (just two levels):
> indicators <- c("V28", "V30", "V31", "V32", "V33",
+ "V35", "V37",
+ "V46", "V47", "V48", "V49", "V50",
+ "V51", "ATTEND")
> for (ind in indicators) {
+ df[,ind] <- as.numeric(df[,ind])
+ }
> group <- "SEX"
>
I now tried to compute a configural model, first with “cfa” and then with “measEq.syntax”, for a single factor:
> model <- "afterlife =~ V30 + V31 + V32 + V33"
> # CONFIGURAL
> config <- cfa(model, data = df, ordered = indicators, parameterization = "theta",
+ estimator = "WLSMV", group = group, group.equal = "")
> config <- measEq.syntax(model, data = df, ordered = indicators, parameterization = "theta",
+ ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016", estimator = "WLSMV",
+ group = group, group.equal = "", return.fit = TRUE)
Error in specify[[1]][RR, CC] : subscript out of bounds
>
The first runs, the second gives an error. However, if I define this 2-factor model:
model <- "afterlife =~ V30 + V31 + V32 + V33
belief.God =~ V28 + V35 + V37"
it runs OK:
> config <- cfa(model, data = df, ordered = indicators, parameterization = "theta",
+ estimator = "WLSMV", group = group, group.equal = "")
> config <- measEq.syntax(model, data = df, ordered = indicators, parameterization = "theta",
+ ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016", estimator = "WLSMV",
+ group = group, group.equal = "", return.fit = TRUE)
>
What am I doing wrong?
I proceeded to compute a model with constrained thresholds and loadings, which tests for metric invariance (scales of both latent response variables and factors fixed across groups):
> metric <- measEq.syntax(model, data = df, ordered = indicators, parameterization = "theta",
+ ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016", estimator = "WLSMV",
+ group = group, group.equal = c("thresholds","loadings"),
return.fit = TRUE)
>
and it runs OK and has good fit (cfi = 0.99936, rmsea = 0.03707745, srmr = 0.0183). I inspected the parameter table and checked that the identification conditions in eq. (19) of Wu and Estabrook are nicely in place (thank you very much for this!!!). However, when I check the estimator:
> unlist(lavInspect(metric, what = "Options"))["estimator"]
estimator
"DWLS"
>
I don’t get “WLSMV”, which is recommended for MGCFA with ordinal variables.
Finally, I have another doubt. I’m trying to run multiple imputed files with measEq.syntax(), but since this takes a long time with “mice”, I thought of saving the imputed data frames in a “mids” object and then running the analysis. However, I could not figure out how to do this with the help page of measEq.syntax().
Thank you very much once again for writing this wonderful function and for your help and attention.
Sincerely,
Carlos M. Lemos
Hi Terrence,
Thank you very much for this function.
I have a problem running my model, it is a longitudinal invariance analysis with categorical variables.
They are a total of seven measurements (DEP0 to DEP6), with nine items in each measurement. The items score from 0 to 3, however, in some follow-ups the items only reach scores from 0 to 2.
model <- "
DEP0 =~ f0V1 + f0V2 + f0V3 + f0V4 + f0V5 + f0V6 + f0V7 + f0V8 + f0V9
DEP1 =~ f1V1 + f1V2 + f1V3 + f1V4 + f1V5 + f1V6 + f1V7 + f1V8 + f1V9
DEP2 =~ f2V1 + f2V2 + f2V3 + f2V4 + f2V5 + f2V6 + f2V7 + f2V8 + f2V9
DEP3 =~ f3V1 + f3V2 + f3V3 + f3V4 + f3V5 + f3V6 + f3V7 + f3V8 + f3V9
DEP4 =~ f4V1 + f4V2 + f4V3 + f4V4 + f4V5 + f4V6 + f4V7 + f4V8 + f4V9
DEP5 =~ f5V1 + f5V2 + f5V3 + f5V4 + f5V5 + f5V6 + f5V7 + f5V8 + f5V9
DEP6 =~ f6V1 + f6V2 + f6V3 + f6V4 + f6V5 + f6V6 + f6V7 + f6V8 + f6V9
"
var0 <- c("f0V1", "f0V2", "f0V3", "f0V4", "f0V5", "f0V6", "f0V7", "f0V8", "f0V9")
var1 <- c("f1V1", "f1V2", "f1V3", "f1V4", "f1V5", "f1V6", "f1V7", "f1V8", "f1V9")
var2 <- c("f2V1", "f2V2", "f2V3", "f2V4", "f2V5", "f2V6", "f2V7", "f2V8", "f2V9")
var3 <- c("f3V1", "f3V2", "f3V3", "f3V4", "f3V5", "f3V6", "f3V7", "f3V8", "f3V9")
var4 <- c("f4V1", "f4V2", "f4V3", "f4V4", "f4V5", "f4V6", "f4V7", "f4V8", "f4V9")
var5 <- c("f5V1", "f5V2", "f5V3", "f5V4", "f5V5", "f5V6", "f5V7", "f5V8", "f5V9")
var6 <- c("f6V1", "f6V2", "f6V3", "f6V4", "f6V5", "f6V6", "f6V7", "f6V8", "f6V9")
longFacNames <- list(FU = c("DEP0","DEP1","DEP2","DEP3","DEP4","DEP5","DEP6"))
syntax.config <- measEq.syntax(configural.model = model, data = DEP_SAL, parameterization = "theta",
ID.fac = "std.lv", ID.cat = "millsap", estimator = "WLSMV",
return.fit = TRUE, longFacNames=longFacNames,
ordered = c("f0V1", "f0V2", "f0V3", "f0V4", "f0V5", "f0V6", "f0V7", "f0V8", "f0V9",
"f1V1", "f1V2", "f1V3", "f1V4", "f1V5", "f1V6", "f1V7", "f1V8", "f1V9",
"f2V1", "f2V2", "f2V3", "f2V4", "f2V5", "f2V6", "f2V7", "f2V8", "f2V9",
"f3V1", "f3V2", "f3V3", "f3V4", "f3V5", "f3V6", "f3V7", "f3V8", "f3V9",
"f4V1", "f4V2", "f4V3", "f4V4", "f4V5", "f4V6", "f4V7", "f4V8", "f4V9",
"f5V1", "f5V2", "f5V3", "f5V4", "f5V5", "f5V6", "f5V7", "f5V8", "f5V9",
"f6V1", "f6V2", "f6V3", "f6V4","f6V5", "f6V6", "f6V7", "f6V8", "f6V9"))
cat(as.character(syntax.config))
summary(syntax.config, fit.measures=TRUE, standardized=TRUE)
When I mail the analysis I get this error
Error in lav_samplestats_step1(Y = Data, ov.names = ov.names, ov.types = ov.types, :
lavaan ERROR: some categories of variable `f0V8' are empty in group 1; frequencies are [1670 45 43 0]
I defined the following items as ordinal and converted to numeric
Error in specify[[1]][RR, CC] : subscript out of bounds
when I check the estimator:
> unlist(lavInspect(metric, what = "Options"))["estimator"]
estimator
"DWLS"
>
I don’t get “WLSMV”, which is recommended for MGCFA with ordinal variables.
Finally, I have another doubt. I’m trying to run multiple imputed files with measEq.syntax(), but since this takes a long time with “mice”, I thought of saving the imputed data frames in a “mids” object and then running the analysis. However, I could not figure out how to do this with the help page of measEq.syntax().
set.seed(12345)
HSMiss <- HolzingerSwineford1939[ , paste("x", 1:9, sep = "")]
HSMiss$x5 <- ifelse(HSMiss$x1 <= quantile(HSMiss$x1, .3), NA, HSMiss$x5)
HSMiss$x9 <- ifelse(is.na(HSMiss$x5), NA, HSMiss$x9)
HSMiss$school <- HolzingerSwineford1939$school
HS.amelia <- amelia(HSMiss, m = 20, noms = "school")
imps <- HS.amelia$imputations
HS.model <- '
visual =~ x1 + a*x2 + b*x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
ab := a*b
'
mgfit1 <- cfa.mi(HS.model, data = imps, std.lv = TRUE, group = "school")
syntax.mi <- measEq.syntax(configural.model = mgfit1, group = "school",
group.equal = c("loadings","intercepts"))
cat(as.character(syntax.mi))
summary(syntax.mi)
## return fit
fit.mi <- measEq.syntax(configural.model = mgfit1, group = "school",
group.equal = c("loadings","intercepts"),
return.fit = TRUE)
summary(fit.mi)
When I mail the analysis I get this error
Error in lav_samplestats_step1(Y = Data, ov.names = ov.names, ov.types = ov.types, :
lavaan ERROR: some categories of variable `f0V8' are empty in group 1; frequencies are [1670 45 43 0]
long.partial = c("x2.1 ~ 1","x2.2 ~ 1"))And this is the warning message:Warning message: In if (long.partial == "") { : the condition has length > 1 and only the first element will be usedIf I use only long.partial="x2 ~ 1", it doesn't give a warning message, but it still equals the intercepts.
If I understood correctly, in long.partial I can not use actual observed-variables names
I have to use a name refered to in longIndNames.
Is there maybe a way to access the automated indicator names, and then use those names when testing partial invariance?
mod.cat <- ' FU1 =~ u1 + u2 + u3 + u4
FU2 =~ u5 + u6 + u7 + u8 '
## the 2 factors are actually the same factor (FU) measured twice
longFacNames <- list(FU = c("FU1","FU2"))
## configural model: no constraints across groups or repeated measures
syntax.config <- measEq.syntax(configural.model = mod.cat, data = datCat,
ordered = paste0("u", 1:8),
parameterization = "theta",
group = "g", longFacNames = longFacNames)
## extract default longIndNames
syntax.config@call$longIndNames
mi.i@external$measEq.syntax@call$longIndNames
Thank you for all the effort put in SEM tools,
It seems that the function is working OK
devtools::install_github("simsem/semTools/semTools")
Is there any function in lavaan to test the equality of the group variance-covariance matrices (general and/or multiple pairwise)?
varnames <- paste0("x", 1:6)
## specify model parameters
covstruc <- outer(varnames, varnames, function(x, y) paste(x, "~~", y))
satMod <- c(paste(varnames, "~ 1"), # mean structure
covstruc[lower.tri(covstruc, diag = TRUE)]) # covariance structure
cat(satMod, sep = "\n")
## fit model, unconstrained across groups
fit1 <- lavaan(satMod, data = HolzingerSwineford1939, group = "school")
## fit model, all parameters constrained across groups
fit0 <- lavaan(satMod, data = HolzingerSwineford1939, group = "school",
group.equal = c("intercepts","residuals","residual.covariances"))
How can I get Chi-square contributions for each group other than using summary()?
fit@test[[1]]$stat.group # naïve
fit@test[[2]]$stat.group # robust
- Some authors state that 0.05 <= RMSEA < 0.08 means acceptable fit, but others (particularly Hu & Bentler (1999)) recommend 0.06 as cutoff for acceptable fit. I'm using the more restrictive limit, particularly because a recent paper warned on the use of criteria based on ML estimation for assessing fit of models with ordinal items obtained using WLSMV:
Do you know about any recent study that recommends criteria for assessing fit in measurement invariance tests of complex models with large samples and ordinal items?
retest.model.FPP <- 'FPP_T1 =~ "PPQ7_Emp"+ "PPQ19_Emp"+ "PPQ25_Emp"+ "PPQ31_Emp"+ "PPQ37_Emp"+"PPQ2_Furchtl"+ "PPQ8_Furchtl"+ "PPQ14_Furchtl"+ "PPQ26_Furchtl"+ "PPQ38_Furchtl"+"PPQ3_Ego"+ "PPQ9_Ego"+ "PPQ21_Ego"+ "PPQ27_Ego"+ "PPQ39_Ego"+"PPQ4_Impul"+ "PPQ10_Impul"+ "PPQ16_Impul"+ "PPQ22_Impul"+ "PPQ40_Impul"+"PPQ11_Soz"+ "PPQ17_Soz"+ "PPQ23_Soz"+ "PPQ29_Soz"+ "PPQ35_Soz" +"PPQ6_Macht" + "PPQ12_Macht"+ "PPQ18_Macht"+ "PPQ24_Macht"+ "PPQ36_Macht"FPP_T2 =~ "PPQ7_A_Emp" + "PPQ19_A_Emp" + "PPQ25_A_Emp"+ "PPQ31_A_Emp"+ "PPQ37_A_Emp" +"PPQ2_A_Furchtl"+ "PPQ8_A_Furchtl"+ "PPQ14_A_Furchtl"+ "PPQ26_A_Furchtl"+ "PPQ38_A_Furchtl"+"PPQ3_A_Ego"+ "PPQ9_A_Ego"+ "PPQ21_A_Ego"+ "PPQ27_A_Ego"+ "PPQ39_A_Ego"+"PPQ4_A_Impul"+ "PPQ10_A_Impul"+ "PPQ16_A_Impul"+ "PPQ22_A_Impul"+ "PPQ40_A_Impul"+"PPQ11_A_Soz"+ "PPQ17_A_Soz" + "PPQ23_A_Soz"+ "PPQ29_A_Soz"+ "PPQ35_A_Soz"+"PPQ6_A_Macht"+ "PPQ12_A_Macht" + "PPQ18_A_Macht"+ "PPQ24_A_Macht"+ "PPQ36_A_Macht"'var1 <- c("PPQ7_Emp", "PPQ19_Emp", "PPQ25_Emp", "PPQ31_Emp", "PPQ37_Emp","PPQ2_Furchtl", "PPQ8_Furchtl", "PPQ14_Furchtl", "PPQ26_Furchtl", "PPQ38_Furchtl","PPQ3_Ego", "PPQ9_Ego", "PPQ21_Ego", "PPQ27_Ego", "PPQ39_Ego","PPQ4_Impul", "PPQ10_Impul", "PPQ16_Impul", "PPQ22_Impul", "PPQ40_Impul","PPQ11_Soz", "PPQ17_Soz", "PPQ23_Soz", "PPQ29_Soz", "PPQ35_Soz" ,"PPQ6_Macht" , "PPQ12_Macht", "PPQ18_Macht", "PPQ24_Macht", "PPQ36_Macht")var2 <- c("PPQ7_A_Emp" , "PPQ19_A_Emp" , "PPQ25_A_Emp", "PPQ31_A_Emp", "PPQ37_A_Emp" ,"PPQ2_A_Furchtl", "PPQ8_A_Furchtl", "PPQ14_A_Furchtl", "PPQ26_A_Furchtl", "PPQ38_A_Furchtl","PPQ3_A_Ego", "PPQ9_A_Ego", "PPQ21_A_Ego", "PPQ27_A_Ego", "PPQ39_A_Ego","PPQ4_A_Impul", "PPQ10_A_Impul", "PPQ16_A_Impul", "PPQ22_A_Impul", "PPQ40_A_Impul","PPQ11_A_Soz", "PPQ17_A_Soz" , "PPQ23_A_Soz", "PPQ29_A_Soz", "PPQ35_A_Soz","PPQ6_A_Macht", "PPQ12_A_Macht" , "PPQ18_A_Macht", "PPQ24_A_Macht", "PPQ36_A_Macht")
longFacNames <- list(FU = c("FPP_T1", "FPP_T2"))measEq.syntax(retest.model.FPP_lavaan, #spezifiziertes Modell über LavaanID.fac = "std.lv", #standardisierter latenter FaktorID.cat = "Wu", #default GuidelinelongFacNames = longFacNames,return.fit = TRUE)
I receive the error message:"lavaan 0.6-3 did not run (perhaps do.fit = FALSE)?** WARNING ** Estimates below are simply the starting values"
Would you help me out and tell me where I went wrong?Here is my code:retest.model.FPP <- 'FPP_T1 =~ "PPQ7_Emp"+ "PPQ19_Emp"+ "PPQ25_Emp"+ "PPQ31_Emp"+ "PPQ37_Emp"+"PPQ2_Furchtl"+ "PPQ8_Furchtl"+ "PPQ14_Furchtl"+ "PPQ26_Furchtl"+ "PPQ38_Furchtl"+"PPQ3_Ego"+ "PPQ9_Ego"+ "PPQ21_Ego"+ "PPQ27_Ego"+ "PPQ39_Ego"+"PPQ4_Impul"+ "PPQ10_Impul"+ "PPQ16_Impul"+ "PPQ22_Impul"+ "PPQ40_Impul"+"PPQ11_Soz"+ "PPQ17_Soz"+ "PPQ23_Soz"+ "PPQ29_Soz"+ "PPQ35_Soz" +"PPQ6_Macht" + "PPQ12_Macht"+ "PPQ18_Macht"+ "PPQ24_Macht"+ "PPQ36_Macht"FPP_T2 =~ "PPQ7_A_Emp" + "PPQ19_A_Emp" + "PPQ25_A_Emp"+ "PPQ31_A_Emp"+ "PPQ37_A_Emp" +"PPQ2_A_Furchtl"+ "PPQ8_A_Furchtl"+ "PPQ14_A_Furchtl"+ "PPQ26_A_Furchtl"+ "PPQ38_A_Furchtl"+"PPQ3_A_Ego"+ "PPQ9_A_Ego"+ "PPQ21_A_Ego"+ "PPQ27_A_Ego"+ "PPQ39_A_Ego"+"PPQ4_A_Impul"+ "PPQ10_A_Impul"+ "PPQ16_A_Impul"+ "PPQ22_A_Impul"+ "PPQ40_A_Impul"+"PPQ11_A_Soz"+ "PPQ17_A_Soz" + "PPQ23_A_Soz"+ "PPQ29_A_Soz"+ "PPQ35_A_Soz"+"PPQ6_A_Macht"+ "PPQ12_A_Macht" + "PPQ18_A_Macht"+ "PPQ24_A_Macht"+ "PPQ36_A_Macht"'
measEq.syntax(retest.model.FPP_lavaan, #spezifiziertes Modell über LavaanID.fac = "std.lv", #standardisierter latenter FaktorID.cat = "Wu", #default GuidelinelongFacNames = longFacNames,return.fit = TRUE)
# Latent common factor means
# Fix T1 common factor mean to zero
C1FAMo ~ 0*1
# Freely estimate T2-T4 common factor means
C2FAMo ~ 1
C3FAMo ~ 1
C4FAMo ~ 1
## LATENT MEANS/INTERCEPTS:
FS_early ~ alpha.1*1 + NA*1
FS_late ~ alpha.2*1 + NA*1
This is my correct code:library(semTools)retest.model.FPP_lavaan <- 'FPP_T1 =~ "PPQ7_Emp"+ "PPQ19_Emp"+ "PPQ25_Emp"+ "PPQ31_Emp"+ "PPQ37_Emp"+"PPQ2_Furchtl"+ "PPQ8_Furchtl"+ "PPQ14_Furchtl"+ "PPQ26_Furchtl"+ "PPQ38_Furchtl"+"PPQ3_Ego"+ "PPQ9_Ego"+ "PPQ21_Ego"+ "PPQ27_Ego"+ "PPQ39_Ego"+"PPQ4_Impul"+ "PPQ10_Impul"+ "PPQ16_Impul"+ "PPQ22_Impul"+ "PPQ40_Impul"+"PPQ11_Soz"+ "PPQ17_Soz"+ "PPQ23_Soz"+ "PPQ29_Soz"+ "PPQ35_Soz" +"PPQ6_Macht" + "PPQ12_Macht"+ "PPQ18_Macht"+ "PPQ24_Macht"+ "PPQ36_Macht"FPP_T2 =~ "PPQ7_A_Emp" + "PPQ19_A_Emp" + "PPQ25_A_Emp"+ "PPQ31_A_Emp"+ "PPQ37_A_Emp" +"PPQ2_A_Furchtl"+ "PPQ8_A_Furchtl"+ "PPQ14_A_Furchtl"+ "PPQ26_A_Furchtl"+ "PPQ38_A_Furchtl"+"PPQ3_A_Ego"+ "PPQ9_A_Ego"+ "PPQ21_A_Ego"+ "PPQ27_A_Ego"+ "PPQ39_A_Ego"+"PPQ4_A_Impul"+ "PPQ10_A_Impul"+ "PPQ16_A_Impul"+ "PPQ22_A_Impul"+ "PPQ40_A_Impul"+"PPQ11_A_Soz"+ "PPQ17_A_Soz" + "PPQ23_A_Soz"+ "PPQ29_A_Soz"+ "PPQ35_A_Soz"+"PPQ6_A_Macht"+ "PPQ12_A_Macht" + "PPQ18_A_Macht"+ "PPQ24_A_Macht"+ "PPQ36_A_Macht"'
retest.model.FPP_lavaan <- '
FPP_T1 =~ PPQ7_Emp+ PPQ19_Emp+ PPQ25_Emp+ PPQ31_Emp+ PPQ37_Emp+
PPQ2_Furchtl+ PPQ8_Furchtl+ PPQ14_Furchtl+ PPQ26_Furchtl+ PPQ38_Furchtl+
PPQ3_Ego+ PPQ9_Ego+ PPQ21_Ego+ PPQ27_Ego+ PPQ39_Ego+
PPQ4_Impul+ PPQ10_Impul+ PPQ16_Impul+ PPQ22_Impul+ PPQ40_Impul+
PPQ11_Soz+ PPQ17_Soz+ PPQ23_Soz+ PPQ29_Soz+ PPQ35_Soz +
PPQ6_Macht + PPQ12_Macht+ PPQ18_Macht+ PPQ24_Macht+ PPQ36_Macht
FPP_T2 =~ PPQ7_A_Emp + PPQ19_A_Emp + PPQ25_A_Emp+ PPQ31_A_Emp+ PPQ37_A_Emp +
PPQ2_A_Furchtl+ PPQ8_A_Furchtl+ PPQ14_A_Furchtl+ PPQ26_A_Furchtl+ PPQ38_A_Furchtl+
PPQ3_A_Ego+ PPQ9_A_Ego+ PPQ21_A_Ego+ PPQ27_A_Ego+ PPQ39_A_Ego+
PPQ4_A_Impul+ PPQ10_A_Impul+ PPQ16_A_Impul+ PPQ22_A_Impul+ PPQ40_A_Impul+
PPQ11_A_Soz+ PPQ17_A_Soz + PPQ23_A_Soz+ PPQ29_A_Soz+ PPQ35_A_Soz+
PPQ6_A_Macht+ PPQ12_A_Macht + PPQ18_A_Macht+ PPQ24_A_Macht+ PPQ36_A_Macht
'
retest.model.FPP_lavaan <- c(paste("FPP_T1 =~", var1), paste("FPP_T2 =~", var2))
retest.model.FPP_lavaan <- c(paste("FPP_T1 =~", var1), paste("FPP_T2 =~", var2),
"PPQ24_A_Macht ~~ PPQ36_A_Macht")
when I try to replicate the syntax with measEq.syntax() (as shown below) the first mean is not fixed to zero but is freely estimated:## LATENT MEANS/INTERCEPTS:
FS_early ~ alpha.1*1 + NA*1 FS_late ~ alpha.2*1 + NA*1
devtools::install_github("simsem/semTools/semTools")
Hello Dr. Terrence.
First of all, Thanks for this measEq.syntax() function as it is really helpful.
I have a doubt regarding the group.equal ( ), as stated underneath:
type5 <- 't5=~ g4q1.5+ g4q3.5+ g4q5.5 + g4q7.5 + g5q1.5+ g5q3.5+ g5q5.5 + g5q7.5
g5q1.5 ~~ g5q5.5'
After testing invariance for threshold and loadings, using group.partial ( ) for one indicator, the fit indices for latent-varaible invariance is given as:
test_th.lo.var<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016", group="group", group.equal=c("thresholds","loadings","lv.variances"), group.partial = c("type5=~g5q3.5","g5q3.5|t1"), warn=T, return.fit=T, ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5", "g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"), parameterization="theta", estimator="wlsmv", meanstructure=TRUE)
chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
#63.048 46.000 0.048 0.955 0.020 0.945
Now, I am interested in testing the invariance of latent-variable means , but it is yielding the same fit-indices as above:
test_th.lo.var.means<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016", group="group", group.equal=c("thresholds","loadings","lv.variances","means"), warn=T, return.fit=T, ordered = c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5", "g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"), parameterization = "theta", estimator="wlsmv", meanstructure=TRUE)
#chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
#63.048 46.000 0.048 0.955 0.020 0.945
I am unsure why these two fit-indices are same, though the “test_th.lo.var.means” has means set to be equal on both groups? Thanks in advance,
Pavneet Kaur
group.equal=c("thresholds","loadings","lv.variances","means")
I am unsure why these two fit-indices are same, though the “test_th.lo.var.means” has means set to be equal on both groups?
Thanks Professor Terrence for your reply back. That was useful.
I have fixed the “intercepts” to be equal in group.equal argument, but even them I am having the same issue.
I am wondering as group.equal for “means” provided a different fit index for model “t2”, but not “t5”. Even if all the indices remain same, I would expect a change in degrees of freedom. Please let me know if I am doing something wrong in the given code:
type5 <- 't5=~ g4q1.5+ g4q3.5+ g4q5.5 + g4q7.5 + g5q1.5+ g5q3.5+ g5q5.5 + g5q7.5
g5q1.5 ~~ g5q5.5'
#===== CHECKING EQUAL VARIANCES ========#
test_th.lo.var<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("loadings","thresholds", "intercepts", "lv.variances"),
warn=T, return.fit=T,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var, fit.measures=T, rsquare=T)
#chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
#62.004 46.000 0.058 0.957 0.019 0.948
#===== CHECKING EQUAL MEANS ========#
test_th.lo.var.means<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016", group="group", group.equal=c("loadings","thresholds", "intercepts", "lv.variances","means"), warn=T, return.fit=T,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var.means, fit.measures=T, rsquare=T)
#chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
# 62.004 46.000 0.058 0.957 0.019 0.948
#==============================#
type2 <- 'type2=~ g4q1.2+ g4q3.2+ g4q5.2+ g4q7.2+ g5q1.2+ g5q3.2+ g5q5.2+ g5q7.2
g4q3.2 ~~ g5q5.2
g5q1.2 ~~ g5q5.2'
#===== CHECKING EQUAL VARIANCES ========#
test_partial_th.lo2.var<-measEq.syntax(type2, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings", "intercepts","lv.variances"),
group.partial=c("type2=~g5q3.2","g5q3.2|t1"),warn=T, return.fit=T,
ordered=c("g4q1.2", "g4q3.2", "g4q5.2", "g4q7.2",
"g5q1.2", "g5q3.2", "g5q5.2", "g5q7.2"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_partial_th.lo2.var, fit.measures=T, rsquare=T)
#chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
# 48.741 43.000 0.253 0.954 0.012 0.940
#===== CHECKING EQUAL MEANS ========#
test_partial_th.lo2_var.means<-measEq.syntax(type2, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings","intercepts", "lv.variances", "means"), warn=T, return.fit=T,
ordered=c("g4q1.2", "g4q3.2", "g4q5.2", "g4q7.2",
"g5q1.2", "g5q3.2", "g5q5.2", "g5q7.2"),
parameterization="theta", estimator="wlsmv", meanstructure=TRUE)
summary(test_partial_th.lo2_var.means, fit.measures=T, rsquare=T)
#chisq.scaled df.scaled pvalue.scaled cfi.scaled rmsea.scaled tli.scaled
# 55.472 44.000 0.115 0.909 0.017 0.884
I have fixed the “intercepts” to be equal in group.equal argument, but even them I am having the same issue.
devtools::install_github("simsem/semTools/semTools")
I am wondering as group.equal for “means” provided a different fit index for model “t2”, but not “t5”.
I have fixed the “intercepts” to be equal in group.equal argument, but even them I am having the same issue.
Do you have the issue after installing the latest development version of semTools?
devtools::install_github("simsem/semTools/semTools")
Yeah, the issue is still there, even though I have installed semTools again using above code.
If so, I suggest you leave return.fit=FALSE so you can inspect the model syntax to see what is happening. For example, does adding "means" to group.equal= not lose a df because the means are already constrained in the model without "means" in group.equal=? Or because adding "means" to group.equal= does not give them the same labels to add the equality constraint, or because the means are not free parameters? This will help me figure out if I need to update the function.
So, I tried to use fit.measure=F and to clarify I am providing the syntax and results underneath:
Return.fit=FALSE
test_th.lo.var <- measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings", "lv.variances"),
warn=T, return.fit=F,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var)
output:
This lavaan model syntax specifies a CFA with 8 manifest indicators (8 of which are ordinal) of 1 common factor(s).
To identify the location and scale of each common factor, the factor means and variances were fixed to 0 and 1, respectively, unless equality constraints on measurement parameters allow them to be freed.
The location and scale of each latent item-response underlying 8 ordinal indicators were identified using the "theta" parameterization, and the identification constraints recommended by Wu & Estabrook (2016). For details, read:
https://doi.org/10.1007/s11336-016-9506-0
Pattern matrix indicating num(eric), ord(ered), and lat(ent) indicators per factor:
t5
g4q1.5 ord
g4q3.5 ord
g4q5.5 ord
g4q7.5 ord
g5q1.5 ord
g5q3.5 ord
g5q5.5 ord
g5q7.5 ord
The following types of parameter were constrained to equality across groups:
thresholds
loadings
lv.variances
test_th.lo.var.means <- measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings","lv.variances","means"),
warn=T, return.fit=F,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var.means)
OUTPUT:
This lavaan model syntax specifies a CFA with 8 manifest indicators (8 of which are ordinal) of 1 common factor(s).
To identify the location and scale of each common factor, the factor means and variances were fixed to 0 and 1, respectively, unless equality constraints on measurement parameters allow them to be freed.
The location and scale of each latent item-response underlying 8 ordinal indicators were identified using the "theta" parameterization, and the identification constraints recommended by Wu & Estabrook (2016). For details, read:
https://doi.org/10.1007/s11336-016-9506-0
Pattern matrix indicating num(eric), ord(ered), and lat(ent) indicators per factor:
t5
g4q1.5 ord
g4q3.5 ord
g4q5.5 ord
g4q7.5 ord
g5q1.5 ord
g5q3.5 ord
g5q5.5 ord
g5q7.5 ord
The following types of parameter were constrained to equality across groups:
thresholds
loadings
lv.variances
means
Syntax and output with fit.measure=T
test_th.lo.var<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings", "lv.variances"),
warn=T, return.fit=T,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var, fit.measures=T, rsquare=T)
lavaan 0.6-3 ended normally after 57 iterations
Optimization method NLMINB
Number of free parameters 42
Number of equality constraints 16
Number of observations per group
0 919
1 926
Estimator DWLS Robust
Model Fit Test Statistic 59.510 63.048
Degrees of freedom 46 46
P-value (Chi-square) 0.087 0.048
Scaling correction factor 0.943
Shift parameter for each group:
0 -0.018
1 -0.018
for simple second-order correction (Mplus variant)
Chi-square for each group:
0 34.682 36.748
1 24.827 26.301
Model test baseline model:
Minimum Function Test Statistic 466.605 432.226
Degrees of freedom 56 56
P-value 0.000 0.000
User model versus baseline model:
Comparative Fit Index (CFI) 0.967 0.955
Tucker-Lewis Index (TLI) 0.960 0.945
Robust Comparative Fit Index (CFI) NA
Robust Tucker-Lewis Index (TLI) NA
Root Mean Square Error of Approximation:
RMSEA 0.018 0.020
90 Percent Confidence Interval 0.000 0.030 0.002 0.031
P-value RMSEA <= 0.05 1.000 1.000
Robust RMSEA NA
90 Percent Confidence Interval NA NA
Standardized Root Mean Square Residual:
SRMR 0.073 0.073
Parameter Estimates:
Information Expected
Information saturated (h1) model Unstructured
Standard Errors Robust.sem
Group 1 [0]:
Latent Variables:
Estimate Std.Err z-value P(>|z|)
t5 =~
g4q1.5 (l.1_) 0.555 0.065 8.508 0.000
g4q3.5 (l.2_) 1.013 0.143 7.076 0.000
g4q5.5 (l.3_) 0.601 0.071 8.412 0.000
g4q7.5 (l.4_) 0.273 0.110 2.471 0.013
g5q1.5 (l.5_) 0.264 0.048 5.540 0.000
g5q3.5 (l.6_) -0.019 0.043 -0.452 0.652
g5q5.5 (l.7_) 0.377 0.065 5.803 0.000
g5q7.5 (l.8_) -0.149 0.044 -3.389 0.001
Covariances:
Estimate Std.Err z-value P(>|z|)
.g5q1.5 ~~
.g5q5.5 (t.7_) 0.173 0.077 2.255 0.024
Intercepts:
Estimate Std.Err z-value P(>|z|)
.g4q1.5 (n.1.) 0.000
.g4q3.5 (n.2.) 0.000
.g4q5.5 (n.3.) 0.000
.g4q7.5 (n.4.) 0.000
.g5q1.5 (n.5.) 0.000
.g5q3.5 (n.6.) 0.000
.g5q5.5 (n.7.) 0.000
.g5q7.5 (n.8.) 0.000
t5 (a.1.) 0.000
Thresholds:
Estimate Std.Err z-value P(>|z|)
g41.5|1 (g41.) 0.638 0.053 12.151 0.000
g43.5|1 (g43.) 0.581 0.073 7.932 0.000
g45.5|1 (g45.) -0.364 0.050 -7.239 0.000
g47.5|1 (g47.) 2.188 0.118 18.611 0.000
g51.5|1 (g51.) -0.132 0.043 -3.064 0.002
g53.5|1 (g53.) -0.094 0.041 -2.275 0.023
g55.5|1 (g55.) 1.349 0.066 20.578 0.000
g57.5|1 (g57.) -0.087 0.042 -2.077 0.038
Variances:
Estimate Std.Err z-value P(>|z|)
.g4q1.5 (t.1_) 1.000
.g4q3.5 (t.2_) 1.000
.g4q5.5 (t.3_) 1.000
.g4q7.5 (t.4_) 1.000
.g5q1.5 (t.5_) 1.000
.g5q3.5 (t.6_) 1.000
.g5q5.5 (t.7_) 1.000
.g5q7.5 (t.8_) 1.000
t5 (p.1_) 1.000
Scales y*:
Estimate Std.Err z-value P(>|z|)
g4q1.5 0.875
g4q3.5 0.703
g4q5.5 0.857
g4q7.5 0.965
g5q1.5 0.967
g5q3.5 1.000
g5q5.5 0.936
g5q7.5 0.989
R-Square:
Estimate
g4q1.5 0.235
g4q3.5 0.506
g4q5.5 0.265
g4q7.5 0.069
g5q1.5 0.065
g5q3.5 0.000
g5q5.5 0.125
g5q7.5 0.022
Group 2 [1]:
Latent Variables:
Estimate Std.Err z-value P(>|z|)
t5 =~
g4q1.5 (l.1_) 0.555 0.065 8.508 0.000
g4q3.5 (l.2_) 1.013 0.143 7.076 0.000
g4q5.5 (l.3_) 0.601 0.071 8.412 0.000
g4q7.5 (l.4_) 0.273 0.110 2.471 0.013
g5q1.5 (l.5_) 0.264 0.048 5.540 0.000
g5q3.5 (l.6_) -0.019 0.043 -0.452 0.652
g5q5.5 (l.7_) 0.377 0.065 5.803 0.000
g5q7.5 (l.8_) -0.149 0.044 -3.389 0.001
Covariances:
Estimate Std.Err z-value P(>|z|)
.g5q1.5 ~~
.g5q5.5 (t.7_) 0.154 0.067 2.289 0.022
Intercepts:
Estimate Std.Err z-value P(>|z|)
.g4q1.5 (n.1.) -0.044 0.071 -0.619 0.536
.g4q3.5 (n.2.) 0.085 0.086 0.997 0.319
.g4q5.5 (n.3.) 0.046 0.070 0.655 0.513
.g4q7.5 (n.4.) 0.022 0.145 0.154 0.878
.g5q1.5 (n.5.) -0.025 0.061 -0.412 0.680
.g5q3.5 (n.6.) 0.025 0.059 0.430 0.667
.g5q5.5 (n.7.) 0.328 0.080 4.114 0.000
.g5q7.5 (n.8.) -0.040 0.059 -0.684 0.494
t5 (a.1.) 0.000
Thresholds:
Estimate Std.Err z-value P(>|z|)
g41.5|1 (g41.) 0.638 0.053 12.151 0.000
g43.5|1 (g43.) 0.581 0.073 7.932 0.000
g45.5|1 (g45.) -0.364 0.050 -7.239 0.000
g47.5|1 (g47.) 2.188 0.118 18.611 0.000
g51.5|1 (g51.) -0.132 0.043 -3.064 0.002
g53.5|1 (g53.) -0.094 0.041 -2.275 0.023
g55.5|1 (g55.) 1.349 0.066 20.578 0.000
g57.5|1 (g57.) -0.087 0.042 -2.077 0.038
Variances:
Estimate Std.Err z-value P(>|z|)
.g4q1.5 (t.1_) 1.000
.g4q3.5 (t.2_) 1.000
.g4q5.5 (t.3_) 1.000
.g4q7.5 (t.4_) 1.000
.g5q1.5 (t.5_) 1.000
.g5q3.5 (t.6_) 1.000
.g5q5.5 (t.7_) 1.000
.g5q7.5 (t.8_) 1.000
t5 (p.1_) 1.000
Scales y*:
Estimate Std.Err z-value P(>|z|)
g4q1.5 0.875
g4q3.5 0.703
g4q5.5 0.857
g4q7.5 0.965
g5q1.5 0.967
g5q3.5 1.000
g5q5.5 0.936
g5q7.5 0.989
R-Square:
Estimate
g4q1.5 0.235
g4q3.5 0.506
g4q5.5 0.265
g4q7.5 0.069
g5q1.5 0.065
g5q3.5 0.000
g5q5.5 0.125
g5q7.5 0.022
test_th.lo.var.means<-measEq.syntax(type5, data=data, ID.cat="Wu.Estabrook.2016",
group="group", group.equal=c("thresholds","loadings","lv.variances","means"),
warn=T, return.fit=T,
ordered=c("g4q1.5", "g4q3.5", "g4q5.5", "g4q7.5",
"g5q1.5", "g5q3.5", "g5q5.5", "g5q7.5"),
parameterization="theta", estimator="wlsmv",
meanstructure=TRUE)
summary(test_th.lo.var.means, fit.measures=T, rsquare=T)
lavaan 0.6-3 ended normally after 57 iterations
Optimization method NLMINB
Number of free parameters 42
Number of equality constraints 16
Number of observations per group
0 919
1 926
Estimator DWLS Robust
Model Fit Test Statistic 59.510 63.048
Degrees of freedom 46 46
P-value (Chi-square) 0.087 0.048
Scaling correction factor 0.943
Shift parameter for each group:
0 -0.018
1 -0.018
for simple second-order correction (Mplus variant)
Chi-square for each group:
0 34.682 36.748
1 24.827 26.301
Model test baseline model:
Minimum Function Test Statistic