3716 views

Skip to first unread message

Aug 26, 2018, 6:02:43 PM8/26/18

to lavaan

For the past few years, there have been 3 separate functions for automatically testing the omnibus null hypotheses of measurement invariance, each under very constrained conditions (e.g., simple structure, first-order measurement).

- measurementInvariance(): only for numeric indicators, only tests across multiple groups
- measurementInvarianceCat(): only for ordered indicators, only tests across multiple groups
- longInvariance(): only for numeric indicators of 1 (repeatedly measured) factor, only tests across repeated measures

In the next version of semTools, these functions will still available, but they have been deprecated.

The new function measEq.syntax() addresses all of the above limitations by writing lavaan model syntax and (optionally) fitting the model to data. By specifying a configural model with minimal syntax (i.e., the pattern of factor loadings), users can set arguments to specify models with invariance constraints imposed on *any* measurement and/or structural parameters across *any* combination of multiple groups and/or repeated measures. The configural model may have complex structure (i.e., cross-loadings and residual covariances are allowed), higher-order factors, and *any* combination of numeric and ordered indicators. Different repeatedly measured factors do *not* need to all be measured on the same number of occasions, nor do the same indicators even need to be measured on all occasions (e.g., if some test items are replaced when measuring abilities at later ages). Beyond the choice between std.lv = TRUE or FALSE in the previous functions, users can also choose effects-coding to identify common-factor distributions (except in models with cross-loadings, such as bifactor models). Users can also choose between 4 different methods for identifying the latent item-response distributions underlying any ordered indicators, no longer being confined to use the Millsap & Tein (2004) recommendations implemented in measurementInvarianceCat(). The measEq.syntax() function (and now the compareFit() function) can also be applied to lavaan.mi objects, if users fit lavaan models to multiple imputations. If a user still has a more complex situation than can be automated using the measEq.syntax() function, it is still possible to specify a model with similar needs, then print the lavaan model syntax and edit it manually to fit the desired model.

I encourage anyone who has had difficulty using the previous invariance-testing functions to try using this, and please let me know if there are any difficulties. Users can install the development version of semTools and find examples on the help page to get started:

`devtools::install_github("simsem/semTools/semTools")`

library(semTools)

?measEq.syntax

The help page examples show how the user can mimic the Console output of the previous functions by looping over a sequence of constraints, saving the fitted models in a list, and comparing their fit using the compareFit() function. Eventually, I will write an interactive wrapper function to do this, integrating follow-up tests rather than only providing omnibus tests that assume all previous levels of invariance are tenable. But that will only be possible after updating the existing partialInvariance() function to work in this new framework, and I probably will not have time to work on that soon.

Thanks for any feedback (please post here so others can benefit),

Terrence D. Jorgensen

Postdoctoral Researcher, Methods and Statistics

Research Institute for Child Development and Education, the University of Amsterdam

UvA web page: http://www.uva.nl/profile/t.d.jorgensen

Aug 28, 2018, 5:51:56 AM8/28/18

to lavaan

Thank you for the hard work Terrence, I'll try this out soon.

Giacomo

Aug 28, 2018, 11:24:25 AM8/28/18

to lavaan

I can't seem to be able to install semTools from GitHub using devtools. Is that possible?

Best,

G

On Sunday, August 26, 2018 at 11:02:43 PM UTC+1, Terrence Jorgensen wrote:

Aug 29, 2018, 8:35:36 AM8/29/18

to lavaan

I can't seem to be able to install semTools from GitHub using devtools

My colleague was able to install it. Do you have the latest R (3.5.1) and development version of lavaan (0.6-3.1295 installed?

And is devtools installed?

`library(devtools)`

Aug 30, 2018, 1:17:16 PM8/30/18

to lav...@googlegroups.com

I can't thank you enough for this, Terrence!

I had a hard time implementing the Wu & Estabrook procedure following your advice, and I can imagine the amount of expertise and work required to write this function.

I will certainly try it and give feedback.

Best regards,

Carlos

--

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.

Sep 10, 2018, 7:06:32 AM9/10/18

to lav...@googlegroups.com

Hi again Terrence,

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

Em seg, 27 de ago de 2018 às 00:02, Terrence Jorgensen <tjorge...@gmail.com> escreveu:

Sep 12, 2018, 3:58:58 PM9/12/18

to lavaan

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]

Sep 12, 2018, 4:29:30 PM9/12/18

to lavaan

I defined the following items as ordinal and converted to numeric

That is not necessary. If the ordered variables are already declared as ordered in the data.frame, then you don't even need to use the ordered= argument.

Error in specify[[1]][RR, CC] : subscript out of bounds

I'm glad you tried it with both a 1-factor and 2-factor model. That helped me locate the source of the error and add a check so that syntax is only written for factor covariances if there are any. If you install the development version of semTools again, it should work

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.

WLSMV is not an estimator. That is just the shortcut for setting estimator = "DWLS", test = "scaled.shifted", and se = "robust.sem". Nothing to do with this function; check the result for config when fitted with cfa().

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().

Yes, you should impute your data first so you are not wasting time repeating the imputation process before each analysis. If you fit the configural model to the imputed data with runMI(), you can submit that object to measEq.syntax(). Here is an example.

`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)

Terrence D. Jorgensen

Postdoctoral Researcher, Methods and Statistics

Research Institute for Child Development and Education, the University of Amsterdam

Sep 12, 2018, 4:37:53 PM9/12/18

to lavaan

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]

The source of the error is clearly noted in the error message, and it has nothing to do with the syntax-writing function. If you fit the configural model yourself with cfa(), you will notice the same error occurs.

Terrence D. Jorgensen

Postdoctoral Researcher, Methods and Statistics

Research Institute for Child Development and Education, the University of Amsterdam

Sep 13, 2018, 11:08:24 AM9/13/18

to lav...@googlegroups.com

Hi again Terrence,

I'm glad that my feedback was somehow useful to you. I will reinstall semTools and retry.

Thank you very much for also clarifying my confusion about WLSMV and for sending the example on combining cfa.mi() with measEq.syntax(). I will make sure to set std.lv = TRUE in the call to cfa.mi (I missed that in my example) and then ID.fac = "std.lv" and ID.cat = "Wu.Estabrook.2016" in the call to measEq.syntax().

Best regards,

Carlos

Oct 12, 2018, 11:56:23 AM10/12/18

to lavaan

Thank you for all the work, and making our work simpler. I am now trying the function to test longitudinal variance, and I'm having problems with the long.partial argument.

This is the syntax:

mi.i<-measEq.syntax(ep.m, data= data,estimator="MLM", meanstructure=TRUE, int.lv.free=TRUE, auto.fix.first=T , orthogonal=T,

ID.fac="auto.fix.first", return.fit = TRUE,

longFacNames = list(X = c("X1","X2")), auto=T,

long.equal=c("loadings","intercepts"), long.partial = c("x2.1 ~ 1","x2.2 ~ 1"))

ID.fac="auto.fix.first", return.fit = TRUE,

longFacNames = list(X = c("X1","X2")), auto=T,

long.equal=c("loadings","intercepts"), 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 used

If I use only long.partial="x2 ~ 1", it doesn't give a warning message, but it still equals the intercepts.

Hope this feedback helps...

Kind regards,

Una

Oct 18, 2018, 7:04:39 AM10/18/18

to lavaan

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.

I fixed the warning by checking the length too, but it is not something that would have actually led to a problem for anyone.

The problem (I am guessing) is that you are using actual observed-variable names in the vector: Are x2.1 and x2.2 the same variable measured at times 1 and 2, respectively?

When you use long.partial="x2 ~ 1", that seems to refer to a name that you give to the observed variables x2.1 and x2.2 that are repeatedly measured (i.e., x2 itself is not a variable in your data set, right?). In order to do that, you need to tell the function that "x2" refers to c("x2.1","x2.2") using the longIndNames= argument, which works the same way longFacNames= does. As described in the "Repeated Measures" part of the Details section in the help page, if you don't make your own longIndNames, the default will be to assign automated labels based on the order, using the factor name. So, since your longitudinal factor is named "X", the automated indicator name would probably be "._X_ind.2" (assuming it is the second indicator specified after =~ both times). But of course I would recommend assigning your own longIndNames so you can use sensible names in long.partial=.

Oct 19, 2018, 5:37:58 AM10/19/18

to lavaan

Thank you for the feedback, I did judge wrongly how the automated
labels would look like. If I understood correctly, in long.partial I can
not use actual observed-variables names (which were in my case
x2.1
and
x2.2, the same variable measured two times, as you supposed) and I have
to use a name refered to in longIndNames. But if I only define one
indicator name in longIndNames, it seems to me it does not test the
invariance of other indicators (those not defined in longIndNames). Some
of my models are quite small, so this is not a problem to name them one
by one. But for the bigger CFA models I was quite thankful for the
automated indicator label formation, but it seems it is not compatible
with partial invariance. I tried a few options of forming the label to
use in long.partial based on ""._factor_ind.1", but did not find the
right one (the one that would free the parameter). Is there maybe a way
to access the automated indicator names, and then use those names when
testing partial invariance?

Thank you for all the effort put in SEM tools,

Oct 20, 2018, 5:17:41 PM10/20/18

to lavaan

If I understood correctly, in long.partial I can not use actual observed-variables names

Correct

I have to use a name refered to in longIndNames.

Correct

Is there maybe a way to access the automated indicator names, and then use those names when testing partial invariance?

Yes, and I should probably add this to the help page. I store it in the @call slot, even if it was not part of your actual call. From the help-page example:

`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

Even when you fit the model, I store the measEq.syntax object in lavaan's @external slot. So from your posted example, here is how you can find the names:

`mi.i@external$measEq.syntax@call$longIndNames`

Thank you for all the effort put in SEM tools,

Happy to help :-)

Oct 25, 2018, 11:51:40 AM10/25/18

to lav...@googlegroups.com

Hi Terrence,

I got a warning with group.partial that looks similar to that reported by Una. I'm testing a 4-factor model with 10 ordinal indicators, and I'm trying to establish partial scalar invariance across countries (group) by freeing one indicator in each factor (Wu & Estabrook parameterization, no marker indicators) and I get the following:

> # SCALAR.PARTIAL

> scalar.partial <- measEq.syntax(model, data = df, ordered = indicators, parameterization = "theta",

> scalar.partial <- 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","intercepts"),

+ group.partial = c("V33~1","V37~1","V48~1","ATTEND~1"), return.fit = TRUE)

Warning message:

In if (group.partial == "") { :

+ group.partial = c("V33~1","V37~1","V48~1","ATTEND~1"), return.fit = TRUE)

Warning message:

In if (group.partial == "") { :

the condition has length > 1 and only the first element will be used

It seems that the function is working OK because using lavInspect(scalar.partial, what = "free") I get (only two countries shown):

...

$`NO-Norway`$nu

intrcp

V30 0

V31 0

V32 0

V33 1579

V28 0

V35 0

V37 1580

V46 0

V47 0

V48 1581

V49 0

V51 0

ATTEND 1582

...

$`IT-Italy`$nu

intrcp

V30 0

V31 0

V32 0

V33 1751

V28 0

V35 0

V37 1752

V46 0

V47 0

V48 1753

V49 0

V51 0

ATTEND 1754

$`NO-Norway`$nu

intrcp

V30 0

V31 0

V32 0

V33 1579

V28 0

V35 0

V37 1580

V46 0

V47 0

V48 1581

V49 0

V51 0

ATTEND 1582

...

$`IT-Italy`$nu

intrcp

V30 0

V31 0

V32 0

V33 1751

V28 0

V35 0

V37 1752

V46 0

V47 0

V48 1753

V49 0

V51 0

ATTEND 1754

i.e. the intercepts I want are free and different for each country. I don't know whether I'm doing some error and/or the warning is serious.

I also have some questions not directly related to measEq.syntax(). If you can answer them without spending much time I would be very grateful:

- Is there any function in lavaan to test the equality of the group variance-covariance matrices (general and/or multiple pairwise)?

- How can I get Chi-square contributions for each group other than using summary()?

- 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?

Thank you once more for your attention, and for
creating measEq.syntax()and making it available!

Carlos M. Lemos

Oct 26, 2018, 6:20:15 AM10/26/18

to lavaan

It seems that the function is working OK

Correct. As I told Una above, that warning is not something that will cause a problem for anyone, so you can ignore it. But in the development version, the warning is gone.

`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)?

No, I've always just written my own saturated-model syntax using outer() with paste().

`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"))

Note that "intercepts" constrains all observed-variable intercepts, even when the intercepts happen to be means because there are no predictors. Likewise, "residuals" and "residual.covariances" constrain all observed-variable *total *(co)variances when there happen to be no predictors.

Because the saturated model has *df *= 0, the test of equality is merely the chi-squared from fit0.

How can I get Chi-square contributions for each group other than using summary()?

Careful, it isn't really valid to compare those across multigroup models once any parameters are constrained across groups. But they are stored in the object's @test slot:

`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:

So you are using a fixed cutoff even though the article said it's not a good idea?

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?

Yes (warning: shameless self-promotion)

You can use the semTools function permuteMeasEq() to do this. This paper is specifically about getting better control of Type I errors with categorical indicators, but see our first paper showing why fixed cutoffs give inconsistent results *even under ideal circumstances*.

It also shows that permutation provides an actual test of configural invariance, even if your model does not fit perfectly. If your configural model does not fit well, consider this paper too:

Terrence D. Jorgensen

Assistant Professor, Methods and Statistics

Oct 27, 2018, 5:25:19 PM10/27/18

to lav...@googlegroups.com

Hi Terrence,

I'm once again extremely grateful for your help and guidance!

Congratulations also for your papers - they are great, and likely to be useful not only to me but to many people in this forum. I hope I get our paper published, to be able to cite yours'!

The problem of the cutoff criteria is very important, because in the model for the countries choosing the cutoff 0.06 instead of 0.08 for RMSEA makes an enormous difference - the RMSEA for the scalar-invariant model is 0.077 (0.075,0.078) [90%].

In fact, two of my questions are related. In a test comparing a metric-invariant model for religious groups (equal thresholds and loadings across groups) with a model with invariant thresholds, loadings, and factor variance-covariance matrices (= lv.variances and = lv.covariances), I get RMSEA = 0.072, so I'm rejecting the hypothesis of invariant \Phi across groups. I thought that one indirect way of confirming the rejection criteria would be to test the \Phi_{g} matrices for equality in the scalar invariant model (metric+intercepts fixed, but \Phi again free across groups), which holds OK (good fit indices, small enough delta(CFI), etc.).

So, I will have to study your papers and your code, which will take me some time...

Best regards, and many, many thanks,

Carlos M. Lemos

Oct 31, 2018, 10:42:59 AM10/31/18

to lavaan

Dear Terrence,

I tried to use the measEq.syntax() function but 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"'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)

Thank you very much in advance!

Johanna

Nov 2, 2018, 9:48:03 AM11/2/18

to lavaan

I receive the error message:"lavaan 0.6-3 did not run (perhaps do.fit = FALSE)?** WARNING ** Estimates below are simply the starting values"

That is a warning, not an error.

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)

It can't be all of your code, because you are not submitting the syntax object to measEq.syntax(). Did you fit the model first? Is that when the warning message was printed?

Terrence D. Jorgensen

Assistant Professor, Methods and Statistics

Research Institute for Child Development and Education, the University of Amsterdam

Nov 4, 2018, 10:59:26 AM11/4/18

to lav...@googlegroups.com

Thanks for the answer. 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"+