[Semtools - compRelSEM] issue obtaining reliability for higher order factor

539 views
Skip to first unread message

adri...@gmail.com

unread,
Sep 25, 2022, 10:54:30 AM9/25/22
to lavaan
I am using Lavaan 0.6.12 and semTools 0.5.6.
I am using ordinal (5-point Likert) items
I am having issues retrieving reliability figure for a higher order construct.
I have specified the following measurement model:
Measurement2ndOrder.cfa <- '
# Experience
Experience      =~ Exp1 + Exp2 + Exp3 + Exp4 + Exp5 + Exp6

# Service Quality
Tangibility     =~ SQTang1 + SQTang2 + SQTang3
Responsiveness  =~ SQResp1 + SQResp2 + SQResp3
Consumables     =~ SQCons1 + SQCons2 + SQCons3
Communication   =~ SQCom1 + SQCom2+ SQCom3
SafetyPerception =~ Saf1 + Saf2 + Saf3 + Saf4 + Saf5 + Saf6 + Saf7 + Saf8

ServQual =~ Tangibility + Responsiveness + Consumables + Communication + SafetyPerception
'
This model is then fit using the following command:
Measurement2ndOrder.fit <- cfa(Measurement2ndOrder.cfa,
                              TsetNorm,
                              ordered = c(
                                "Exp1", "Exp2", "Exp3", "Exp4", "Exp5", "Exp6",
                                "SQTang1", "SQTang2", "SQTang3",
                                "SQResp1", "SQResp2", "SQResp3",
                                "SQCons1", "SQCons2", "SQCons3",
                                "SQCom1", "SQCom2",  "SQCom3",
                                "Saf1", "Saf2", "Saf3", "Saf4", "Saf5", "Saf6", "Saf7", "Saf8"
                              ),
                              std.lv=TRUE,
                              estimator = "WLSMV"
)
When I use the command:
compRelSEM(Measurement2ndOrder.fit, higher = "ServQual")

I get the following error message:
Error in diag(scales) %*% truevar : non-conformable arguments

The deprecated command reliabilityL2(Measurement2ndOrder.fit, "ServQual")
returns the following information
       omegaL1        omegaL2 partialOmegaL1
     0.8200871      0.9162002      0.9632436

Not sure whether/what I am doing wrong.

Ian McPhail

unread,
Sep 29, 2022, 6:29:02 PM9/29/22
to lavaan
Thanks for posting on the compRelSEM function. I am attempting to use this function and am getting a different error warning, and I'd like to add it to this thread in the hopes someone can point out what I might be doing wrong.

The error reads as: "Error in PHI[[block.label[b]]] <- Reduce("+", lapply(phiList, "[[", i = b))/m :
  attempt to select less than one element in OneIndex"

The syntax I am using is: (note that I am using multiple imputed datasets and the "cfa.mi" function)

cssoc.pre1.sf.rtc <-
                   '
                   F1 =~ BUMBYMQ3_PRE_r + BUMBYMQ7_PRE_r + BUMBYMQ8_PRE_r + BUMBYMQ12_PRE_r + BUMBYMQ15_PRE_r
                   F2 =~ BUMBYMQ24_PRE_r + BUMBYMQ26_PRE_r + BUMBYMQ28_PRE_r + BUMBYMQ32_PRE_r + BUMBYMQ33_PRE_r
                   '

fit.cssoc.pre1.sf.rtc <- cfa.mi(
  cssoc.pre1.sf.rtc,
  data = rtc.mi.list,
  m = 50,
  ordered = c("BUMBYMQ3_PRE_r", "BUMBYMQ7_PRE_r", "BUMBYMQ8_PRE_r", "BUMBYMQ12_PRE_r", "BUMBYMQ15_PRE_r",  
              "BUMBYMQ24_PRE_r", "BUMBYMQ26_PRE_r", "BUMBYMQ28_PRE_r", "BUMBYMQ32_PRE_r", "BUMBYMQ33_PRE_r"),
  miPackage = "mice",
  std.lv = TRUE,
  parameterization = "theta",
  estimator = "WLSMV",
  seed = 1234)

cssoc.pre1.sf.rtc.rel <- compRelSEM(fit.cssoc.pre1.sf.rtc, ord.scale = TRUE, return.total = TRUE)

Terrence Jorgensen

unread,
Oct 8, 2022, 5:37:08 PM10/8/22
to lavaan
Error in diag(scales) %*% truevar : non-conformable arguments

This error was due to a bug when observed indicators were ordered= and requesting reliability for a higher-order factor.  Now fixed in the development version:

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

Terrence D. Jorgensen
Assistant Professor, Methods and Statistics
Research Institute for Child Development and Education, the University of Amsterdam

Terrence Jorgensen

unread,
Oct 8, 2022, 5:39:25 PM10/8/22
to lavaan
"Error in PHI[[block.label[b]]] <- Reduce("+", lapply(phiList, "[[", i = b))/m :
  attempt to select less than one element in OneIndex"

 Strange, why would there be less than 1 element?  I don't have any "OneIndex" in my source code, but it sounds like you might have no elements in phiList, which would only happen if the model did not converge on any imputations.  Without a reprex (including data), I would not be able to track down the error.

Adriaan Van Liempt

unread,
Oct 8, 2022, 5:54:31 PM10/8/22
to lav...@googlegroups.com
Thank you. Will update.

Would it be possible to Bootstrap the confidence intervals of the reliability measure?

--
You received this message because you are subscribed to a topic in the Google Groups "lavaan" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/lavaan/NNiGnwY-pwQ/unsubscribe.
To unsubscribe from this group and all its topics, send an email to lavaan+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/lavaan/23eb25ea-8fe6-45f1-ab7b-92027cf17567n%40googlegroups.com.

Terrence Jorgensen

unread,
Oct 8, 2022, 6:10:23 PM10/8/22
to lavaan
Would it be possible to Bootstrap the confidence intervals of the reliability measure?

Sure, go ahead.  You can use the bootstrapLavaan() function, with custom function for the FUN= argument 

highRel <- function(x) {
  library(semTools)
  compRelSEM(x, higher = "ServQual")

Christian Arnold

unread,
Oct 8, 2022, 6:20:28 PM10/8/22
to lav...@googlegroups.com
Of course Terrence is right, as always. However, the problem from a practical point of view is that bootstrapLavaan doesn't give you a CI. Maybe x.boot will help you (see forum). I am happy to receive feedback.

Best

Christian 

From: lav...@googlegroups.com <lav...@googlegroups.com> on behalf of Terrence Jorgensen <tjorge...@gmail.com>
Sent: Sunday, October 9, 2022 12:10:23 AM
To: lavaan <lav...@googlegroups.com>
Subject: Re: [Semtools - compRelSEM] issue obtaining reliability for higher order factor
 
--
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 view this discussion on the web visit https://groups.google.com/d/msgid/lavaan/bd4282fa-ad78-4e69-9185-d9b8e4ac6427n%40googlegroups.com.

Adriaan Van Liempt

unread,
Oct 10, 2022, 2:02:45 PM10/10/22
to lav...@googlegroups.com
With Terrence's advice, when set to verbose, I get a list of outcomes from all successful draws. After ordering this list, could I mention the 2.5th and the 97.5th percentile as the 95% confidence interval?

@Christian
I have tried to work with x.boot (from the following post, correct? https://groups.google.com/g/lavaan/c/0RSsh4M6zQg ), but I get the following error:
> HS.boot <- x.boot(Measurement2ndOrder.fit, R = 80, verbose = TRUE)
Error in FUN(object, ...) : unused argument (warn = -1)
Called from: lav_bootstrap_internal(object = object, lavdata. = NULL, lavmodel. = NULL,
    lavsamplestats. = NULL, lavoptions. = lavoptions., lavpartable. = NULL,
    R = R, type = type., verbose = verbose, FUN = FUN, return.boot = return.boot,
    h0.rmsea = h0.rmsea, ...)

Not entirely certain what I am doing wrong as I try to follow your most basic example from the forumlink I pasted above. I get the same error adding a FUN (tried both FUN = highRel and FUN = "highRel" ) to the command.



You received this message because you are subscribed to a topic in the Google Groups "lavaan" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/lavaan/NNiGnwY-pwQ/unsubscribe.
To unsubscribe from this group and all its topics, send an email to lavaan+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/lavaan/AM0PR04MB5955FD5A31C6F4BE799BE71FFC5E9%40AM0PR04MB5955.eurprd04.prod.outlook.com.

Christian Arnold

unread,
Oct 11, 2022, 12:04:27 AM10/11/22
to lav...@googlegroups.com
Hi Adriaan,

thanks for your feedback. This one works for me:

library(lavaan)
library(semTools)
source("x.boot.r")

HS.model <- "
visual  =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed   =~ x7 + x8 + x9
higher  =~ visual + textual + speed
"
set.seed(1)
fit <- sem(HS.model, HolzingerSwineford1939)
HS.boot <- x.boot(fit, R = 400, verbose = TRUE)

myHighRel <- function(object, higher) {
  compRelSEM(object, higher = higher)
}

summary(
  x.boot.out(HS.boot, myFUN = myHighRel, myArgs = list(higher = "higher"), ci.type = "bc")
)

Best

Christian


Von: lav...@googlegroups.com <lav...@googlegroups.com> im Auftrag von Adriaan Van Liempt <adri...@gmail.com>
Gesendet: Montag, 10. Oktober 2022 20:02
An: lav...@googlegroups.com <lav...@googlegroups.com>
Betreff: Re: [Semtools - compRelSEM] issue obtaining reliability for higher order factor
 

Shu Fai Cheung (張樹輝)

unread,
Oct 11, 2022, 12:24:00 AM10/11/22
to lavaan
I also got the error on unused argument. I used the script from the Google Group:


Maybe it is related to the version of lavaan or semTools?

> library(lavaan)
This is lavaan 0.6-12
lavaan is FREE software! Please report any bugs.
> library(semTools)
 
###############################################################################
This is semTools 0.5-6
All users of R (or SEM) are invited to submit functions or ideas for functions.
###############################################################################

> source("x.boot.r")
>
> HS.model <- "
+ visual  =~ x1 + x2 + x3
+ textual =~ x4 + x5 + x6
+ speed   =~ x7 + x8 + x9
+ higher  =~ visual + textual + speed
+ "

> set.seed(1)
> fit <- sem(HS.model, HolzingerSwineford1939)
> HS.boot <- x.boot(fit, R = 400, verbose = TRUE)
Error in FUN(object, ...) : unused argument (warn = -1)
> sessionInfo()
R version 4.2.1 (2022-06-23 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19044)

Matrix products: default

locale:
[1] LC_COLLATE=Chinese (Traditional)_Hong Kong SAR.utf8
[2] LC_CTYPE=Chinese (Traditional)_Hong Kong SAR.utf8
[3] LC_MONETARY=Chinese (Traditional)_Hong Kong SAR.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=Chinese (Traditional)_Hong Kong SAR.utf8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] semTools_0.5-6 lavaan_0.6-12

loaded via a namespace (and not attached):
 [1] codetools_0.2-18   lattice_0.20-45    mvtnorm_1.1-3      zoo_1.8-11
 [5] MASS_7.3-58.1      grid_4.2.1         xtable_1.8-4       jsonlite_1.8.2
 [9] stats4_4.2.1       coda_0.19-4        rlang_1.0.6        estimability_1.4.1
[13] cli_3.4.1          multcomp_1.4-20    Matrix_1.5-1       pbivnorm_0.6.0
[17] sandwich_3.0-2     TH.data_1.1-1      splines_4.2.1      emmeans_1.8.1-1
[21] parallel_4.2.1     survival_3.4-0     compiler_4.2.1     mnormt_2.1.1
>

Hope the session info can help finding the source of the error.

-- Shu Fai

Christian Arnold

unread,
Oct 11, 2022, 2:15:59 AM10/11/22
to lavaan
Hi Shu Fai,

Thank you! Super helpful! In fact I was still using lavaan 0.6-11. With lavaan 0.6-12 I get the error too. I have uploaded a revised version: https://groups.google.com/g/lavaan/c/0RSsh4M6zQg/m/V2zKKV7FAwAJ
I hope this fixes the problem :-)

Best

Christian

adri...@gmail.com

unread,
Oct 11, 2022, 4:52:40 AM10/11/22
to lavaan
Hi Christian,

Thank you very much for providing the example.

I have updated the x.boot.r and I have ran your code and I was to able successfully run your code using the HolzingerSwineford1939 data.

I then proceeded with my data, which is ordinal, and regrettably it did not seem to work:
> summary(
+   x.boot.out(HS.boot, myFUN = highRel, myArgs = list(higher = "higher"), ci.type = "bc")
+ )
Error in h(simpleError(msg, call)) :
  error in evaluating the argument 'object' in selecting a method for function 'summary': unused argument (object = new("lavaan", version = "0.6.12", call = lavaan(slotOptions = lavoptions, slotParTable = lavpartable, slotSampleStats = bootSampleStats, slotData = lavdata, slotModel = model.boot), timing = list(0, 0, 0, 0, 0, 0.0200000000186265, 0, 0, 0, 0, 0, 0, 0.42000000004191, 0, 0, 0.0200000000186265, 0.0399999999208376, 0, 0.0500000000465661, 0, 0.550000000046566), Options = list("cfa", "lavaan", TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, "", FALSE, "delta", FALSE, TRUE,
    TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, "geomin", "bordered", list(FALSE, "none", TRUE, 0.001, 1, 0, 0, numeric(0), numeric(0), 100, "gpa", TRUE, "index", 1e-05, 1e-08, FALSE, FALSE, TRUE, 10000), FALSE, "listwise", "total", FALSE, TRUE, FALSE, "default", character(0), character(0), FALSE, NULL, "DWLS", list(), "default", "default", "LISREL", TRUE, "none", "robust.sem", "scaled.shifted", c
Called from: h(simpleError(msg, call))

This is my model:
Measurement2ndOrder.cfa <- '
# Experience
Experience      =~ Exp1 + Exp2 + Exp3 + Exp4 + Exp5 + Exp6

# Service Quality
Tangibility     =~ SQTang1 + SQTang2 + SQTang3
Responsiveness  =~ SQResp1 + SQResp2 + SQResp3
Consumables     =~ SQCons1 + SQCons2 + SQCons3
Communication   =~ SQCom1 + SQCom2+ SQCom3

# Perceptions Safety Precautions

SafetyPerception =~ Saf1 + Saf2 + Saf3 + Saf4 + Saf5 + Saf6 + Saf7 + Saf8

higher =~ Tangibility + Responsiveness + Consumables + Communication + SafetyPerception
'

This is how I have fitted the model:
Measurement2ndOrder.fit <- cfa(Measurement2ndOrder.cfa,
                              TsetNorm,
                              ordered = c(
                                "Exp1", "Exp2", "Exp3", "Exp4", "Exp5", "Exp6",
                                "SQTang1", "SQTang2", "SQTang3",
                                "SQResp1", "SQResp2", "SQResp3",
                                "SQCons1", "SQCons2", "SQCons3",
                                "SQCom1", "SQCom2",  "SQCom3",
                                "Saf1", "Saf2", "Saf3", "Saf4", "Saf5", "Saf6", "Saf7", "Saf8"
                              ),
                              std.lv=TRUE,
                              estimator = "WLSMV"
)

I use this following code:
HS.boot <- x.boot(Measurement2ndOrder.fit, R = 80, verbose = TRUE) # low amount of draws for testing
#code works, no errors
summary(
  x.boot.out(HS.boot, myFUN = highRel, myArgs = list(higher = "higher"), ci.type = "bc")
)
# get error

Hope this information is helpful to you.

Best,
Adriaan

Christian Arnold

unread,
Oct 12, 2022, 9:14:45 AM10/12/22
to lav...@googlegroups.com
Hm, this seems to be an interesting problem. x.boot doesn't actually have any problems with categorical variables:


library(lavaan)
library(semTools)
source("x.boot.r")

HS.model <- "
visual  =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed   =~ x7 + x8 + x9
higher  =~ visual + textual + speed
"

data <- matrix(NA, nrow = nrow(HolzingerSwineford1939), ncol = 9)
for(i in 1 : 9) {
  x <- c(HolzingerSwineford1939[paste0("x", i)] - min(HolzingerSwineford1939[paste0("x", i)]))[[paste0("x", i)]]
  data[, i] <- cut((x / max(x)), breaks = c(-Inf, 0.25, 0.5, 0.75, Inf), labels = c(1, 2, 3, 4))
}
colnames(data) <- paste0("x", rep(1 : 9))

     
set.seed(1)
fit <- sem(HS.model, data, estimator = "WLSMV", ordered = TRUE)
HS.boot <- x.boot(fit, R = 80, verbose = TRUE)


myHighRel <- function(object, higher) {
  compRelSEM(object, higher = higher)
}

summary(
  x.boot.out(HS.boot, myFUN = myHighRel, myArgs = list(higher = "higher"), ci.type = "bc")
)


I can only speculate, unfortunately. Perhaps compRelSEM produces an error when it evaluates a particular bootstrap draw. Maybe adding an error handler to myHighRel will help.



Gesendet: Dienstag, 11. Oktober 2022 10:52
An: lavaan <lav...@googlegroups.com>

Adriaan Van Liempt

unread,
Oct 12, 2022, 10:20:27 AM10/12/22
to lav...@googlegroups.com
Hi Christian, Could it have something to do with unsuccessful bootstrap runs? I've tried including an error handler, to no avail.

Christian Arnold

unread,
Oct 12, 2022, 10:32:01 AM10/12/22
to lav...@googlegroups.com
Non-converged bootstrap draws are automatically excluded. Maybe it is due to a bootstrap draw that is not admissible (see lavInspect "post.check"). x.boot will provide some information if you inspect the summary of the bootstrap draws. However, these are all just pure guesses. If you send me a snippet of your data, I might be able to reproduce the problem. Gladly by mail.

Best

Christian
From: lav...@googlegroups.com <lav...@googlegroups.com> on behalf of Adriaan Van Liempt <adri...@gmail.com>
Sent: Wednesday, October 12, 2022 4:20:11 PM
To: lav...@googlegroups.com <lav...@googlegroups.com>

Adriaan Van Liempt

unread,
Oct 12, 2022, 5:51:41 PM10/12/22
to lav...@googlegroups.com
Thank you Christian,

The problem wasn't because of Terrence's function, it was a convergence issue.

Christian Arnold

unread,
Oct 12, 2022, 5:58:49 PM10/12/22
to lav...@googlegroups.com
If I program a new x.boot version, I will take this into account and implement warning messages. Thanks for your feedback. Good luck with your project!

Best

Christian 

From: lav...@googlegroups.com <lav...@googlegroups.com> on behalf of Adriaan Van Liempt <adri...@gmail.com>
Sent: Wednesday, October 12, 2022 11:51:24 PM
Reply all
Reply to author
Forward
0 new messages