Reliability questions on CFA

3,326 views
Skip to first unread message

robin ghertner

unread,
Aug 2, 2013, 12:01:29 AM8/2/13
to lav...@googlegroups.com
I'm running a single-factor CFA and am exploring the different options to calculate reliability. All 22 indicator variables are dichotomous (0/1), and I'm using the ordered statement to use WLSV for estimation. The fit stats support a single factor.

First, I used Cronbach's alpha (equivalent to KR-20, for dichotomous variables), from the psych package:
cronbach.alpha(dt)

This gave me an estimate of 0.778
I then used the reliability function in semTools:
reliability(dt_fit)
and got alpha=0.8894284

From the help file for semTools, alpha is Cronbach's alpha. The notes in the output said it used polychoric correlation, not Pearson's, that's the only difference I can find. Would that account for a 0.1 difference, which seems like quite a bit? To make sure something wasn't wrong with cronbach.alpha, I ran the same data through Stata and got the same result. Any thoughts on this much appreciated.

A second type of reliability I ran was composite reliability (Fornell and Larcker, 1981), which is basically the square of the sum of the loadings, divided by the square of sum of loadings plus sum of error variance:
l<-inspect(dt_fit,"coef")$lambda
v<-diag(inspect(dt_fit,"coef")$theta)
cr<-sum(l)^2/(sum(l)^2+sum(v))

This results in 0.96.

Again, this is a big wow for me. I would expect Cronbach's alpha to underestimate reliability, since my loadings do not allow for essential tau-equivalence. But I did not expect a bump from 0.78 to 0.96.

Does this seem like a reasonable finding? Am I missing something here?

Any thoughts much appreciated.

Robin


Terrence Jorgensen

unread,
Aug 2, 2013, 7:01:26 AM8/2/13
to lav...@googlegroups.com
From the source code, I don't see an obvious difference in the way alpha is calculated in either package's function:

from {psych}
alpha.raw <- (1 - tr(C)/sum(C)) * (n/(n - 1))
from {semTools}
alpha[j] <- length(index)/(length(index) - 1) * (1 - sum(diag(sigma)) / sum(sigma))
-- replace C with sigma and n with length(index), and the tr(x) function is equivalent to sum(diag(x)).

But I wouldn't use either one of them.  Like you said, your indicators are not tau-equivalent (a rarely tenable assumption).  Does the value of omega differ between the two packages as well?

The jump in reliability between alpha and CR looks large, but not necessarily surprising.  There is another similar reliability measure calculated from factor loadings called Average Variance Extracted, in which you simply square your factor loadings before taking the sum (a sum of squares, not a squared sum):

CR <- sum(l)^2 / (sum(l)^2 + sum(v))
AVE <- sum(l^2) / (sum(l^2) + sum(v))

Notice that these values also differ, but one is not more right than the other.  These reliability estimates are all measures of the same concept, but they are not estimates of a single common value you would find in the population (i.e., AVE and CR would usually differ in the population, and both would differ from Omega).  Excluding alpha, if all these reliability measures (omega, CR, and AVE) indicate high reliability, then that's a consistent indication of high reliability using different quantifications of a common concept, and that's all I would expect that you need to report.

Terry

Sunthud Pornprasertmanit

unread,
Aug 4, 2013, 2:01:24 AM8/4/13
to lav...@googlegroups.com
When data consist of categorical variables, Pearson's correlations are the correlation of the observed scores. Polychoric correlations are the correlation of the latent multivariate normal distribution underlying categorical variables. Polychoric correlations always have better values. This is why when we run continuous CFA and CFA for categorical variables on categorical data, the factor loadings in CFA with categorical data are higher. If you believe that your items have underlying latent multivariate normal distribution, polychoric correlation and categorical CFA is recommended.

I agree with Terry that if your scales are not tau-equivalent, coefficient alpha is not recommended. Coefficient alpha is the lower bound of population reliability. The estimate of population reliability is the coefficient omega. 

Best,
Sunthud




--
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 http://groups.google.com/group/lavaan.
For more options, visit https://groups.google.com/groups/opt_out.
 
 

robin ghertner

unread,
Aug 4, 2013, 11:31:12 PM8/4/13
to lav...@googlegroups.com
Thanks Terry and Sunthud -

Those were my thoughts precisely. One question - I've run semTools omega (using the reliability function) and the omega function from the psych package. They do not come out the same.

semTools gives me
Omega=0.9011911, Omega2=0.9011911, Omega3=0.9529614
where Omega is from Rakov (2001), Omega2 and 3 from McDonald (1999).

psych gives me (computed using: omega(dt, nfactors=1)):
Omega Hierarchical:    0.78
Omega H asymptotic:    0.99
Omega Total            0.79

I understand fully the semTools documentation on how Omega is calculated in that package, but the psych documentation is a little dense for me, though I understand more or less what it is doing. My first question is basic: are these functions supposed to calculate the same thing? If so, why am I getting different results?

Thanks for any guidance you may be able to provide, I know this is a Lavaan forum and I don't want to get super off-topic.

Robin

robin ghertner

unread,
Aug 4, 2013, 11:35:39 PM8/4/13
to lav...@googlegroups.com
One thing I realized - semTools is using polychoric correlation - I don't know if psych is making that adjustment (my data are binary). That may be a factor. I don't know if it would explain why the alpha numbers are different  - I always understood that with binary data alpha became KR-20, and never really thought about whether the covariance matrix used in alpha needed to be a polychoric/tetrachoric matrix. If that's the case, then perhaps that would go part of the way in explaining the difference.

Sunthud Pornprasertmanit

unread,
Aug 5, 2013, 12:59:43 AM8/5/13
to lav...@googlegroups.com
Could you try estimating the reliability using the semTools package but not specifying your items as categorical (i.e., treat it as continuous)? The number should match.

FYI, I found this article interesting for reliability calculation for binary items. The authors provided Mplus script, which I think it can be easily translated to lavaan. I am thinking of adding this approach in the reliability function.

http://www.tandfonline.com.www2.lib.ku.edu:2048/doi/abs/10.1080/10705511003659417#.Uf8tPW1lV8E

Sunthud

robin ghertner

unread,
Aug 6, 2013, 12:25:27 AM8/6/13
to lav...@googlegroups.com
The key was the ordered argument - semTools reads that and understands that I'm dealing with ordinal/binary data. The other functions (both from psych and from ltm) must assume the data are continuous. The same goes for the Omega function.

However, what is confusing me (in terms of interpretation and which output to use) is that running KR20 (in Stata) produces practically the same estimate as alpha from ltm and psych and semTools w/o ordered statement. KR20 assumes binary data, and is the standard from my limited experience when using binary data (like alpha for continuous). Since semTools is using the polychoric/tetrachoric matrix as the input, I see that it produces a higher estimate of alpha than KR20. Any thoughts about which is more appropriate? This is the first time I've seen a solution to calculating reliability for binary outcomes based off of the tetrachoric correlation, so I'm not really sure where to take this.

BTW I couldn't access the paper you linked to, it sent me to a U of Kansas login screen.

Robin

kma...@aol.com

unread,
Aug 6, 2013, 8:50:34 AM8/6/13
to lav...@googlegroups.com
I could be missing something, but, for what it is worth, it seems to me
that reliability estimates based on polychoric correlations require
extreme care in interpretation. Normally, one wishes to estimate the
reliability of the observed scores. It seems to me that reliability
estimates based on polychorics estimate the reliability one would
obtain were one able to observe the underlying continuous scores. In
many contexts, that may not be a counterfactual that has much practical
import.

Keith
------------------------
Keith A. Markus
John Jay College of Criminal Justice, CUNY
http://jjcweb.jjay.cuny.edu/kmarkus
Frontiers of Test Validity Theory:
http://www.routledge.com/books/details/9781841692203/


Sunthud Pornprasertmanit

unread,
Aug 6, 2013, 9:14:27 AM8/6/13
to lav...@googlegroups.com
I agree that users need to be careful. On one hand, the current practice of calculating alpha (or KR-20) is equivalent to the tau-equivalent model treating dichotomous variables as continuous variables. The normality assumption is violated. On the other hand, reliability is defined as the proportion of true score variance over observed score variance. Using polychoric correlation does not calculate totally based on "observed" score variance but "latent" score variance. Thus, I think both ways are not ideal solutions. If I have no choice, I may go conservative by using KR-20 (the current practice). Some guidelines are available for interpretation of the magnitude. I think the following article may provide the solution.

Raykov, T., Dimitrov, D. M., & Asparouhov, T. (2010). Evaluation of scale reliability with binary measures using latent variable modeling. Structural Equation Modeling, 17(2), 265-279.

Sunthud


--
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+unsubscribe@googlegroups.com.

robin ghertner

unread,
Mar 7, 2014, 11:01:13 PM3/7/14
to lav...@googlegroups.com
Wanted to close the loop on this, as it just dawned on me that Fornell & Larcker's composite reliability and McDonald's Omega are the "same" (well, not really, Tenko Raykov pointed it out to me a while back and I never got around to thinking about it until now). Sorry if this is too basic for everyone - sometimes it takes me a while to realize simple things.

See the formula's laid out here:

F&L, equation 10: http://faculty-gsb.stanford.edu/larcker/PDF/6%20Unobservable%20Variables.pdf
Omega, equation 2 (page 6): http://mailer.fsu.edu/~akamata/papers/MD_rel_paper.pdf

The notation is slightly different, but the meaning is the same. The key for my simple brain to remember was the denominator for Omega is the sample variance-covariance matrix, which is the same as the squared loadings + residual variance, which is the denominator in F&L composite reliability.

You can get composite reliability to equal omega but you first have to make sure you are getting the freely estimated factor loadings, so you have to fix the factor variance to 1 or whatever. Then you simply follow F&L's formula, and it will give you Omega.

library(MASS)
library
(lavaan)
library
(psych)

#Generate 3 multivariate normal variables, with differing degrees of correlation.
Sigma<-matrix(c(1,.6,.3,
               
.6,1,.4,
               
.3,.4,1),3,3)
data
<-as.data.frame(mvrnorm(n=1000,rep(0,3),Sigma))
colnames
(data)<-c("x1","x2","x3")

#Specify and estimate a simple CFA model, fixing the variance to 1.
model
<-"y=~NA*x1+x2+x3
        y~~1*y"

model_fit
<-cfa(model,data=data)
summary
(model_fit)

#extract from the results the residual variances and the loadings.
v
<-diag(inspect(model_fit,"coef")$theta)
l
<-inspect(model_fit,"coef")$lambda

#calculate composite reliability
cr
<-sum(l)^2/(sum(l)^2+sum(v))

#Calculate omega
omega
(data)$omega_h
#results differ slightly, probably due to differences in computation.

Now if you have a multidimensional construct, or you are doing more complicated paths between indicator variables, then it may not come out exactly as Omega in the psych package will get you. That's because Omega isn't really estimating your CFA - it's doing an EFA, some rotation, etc. (you can read the help file), so it won't give you the same thing. But I would argue that since you are really estimating a confirmatory model and you have delineated the relationships, that's the reliability you want to test, not the reliability that Omega as estimated by psych gives you.




On Tuesday, August 6, 2013 9:14:27 AM UTC-4, Sunthud Pornprasertmanit wrote:
I agree that users need to be careful. On one hand, the current practice of calculating alpha (or KR-20) is equivalent to the tau-equivalent model treating dichotomous variables as continuous variables. The normality assumption is violated. On the other hand, reliability is defined as the proportion of true score variance over observed score variance. Using polychoric correlation does not calculate totally based on "observed" score variance but "latent" score variance. Thus, I think both ways are not ideal solutions. If I have no choice, I may go conservative by using KR-20 (the current practice). Some guidelines are available for interpretation of the magnitude. I think the following article may provide the solution.

Raykov, T., Dimitrov, D. M., & Asparouhov, T. (2010). Evaluation of scale reliability with binary measures using latent variable modeling. Structural Equation Modeling, 17(2), 265-279.

Sunthud
On Tue, Aug 6, 2013 at 7:50 AM, kma...@aol.com <kma...@aol.com> wrote:
I could be missing something, but, for what it is worth, it seems to me that reliability estimates based on polychoric correlations require extreme care in interpretation.  Normally, one wishes to estimate the reliability of the observed scores.  It seems to me that reliability estimates based on polychorics estimate the reliability one would obtain were one able to observe the underlying continuous scores.  In many contexts, that may not be a counterfactual that has much practical import.

Keith
------------------------
Keith A. Markus
John Jay College of Criminal Justice, CUNY
http://jjcweb.jjay.cuny.edu/kmarkus
Frontiers of Test Validity Theory: http://www.routledge.com/books/details/9781841692203/



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

Hellen Geremias dos Santos

unread,
Mar 22, 2014, 8:38:09 AM3/22/14
to lav...@googlegroups.com
Dear Terrence Jorgensen and lavaan users

Is it possible to calculate confidence interval for Composite Reliability and Average Variance Extracted?

Any help would be much appreciated.
Hellen

robin ghertner

unread,
Mar 23, 2014, 11:52:27 PM3/23/14
to lav...@googlegroups.com
Hellen

Raykov has published extensively on this. This article by Padilla and Divers seems to be the latest, which summarizes Raykov's work as well as adding some new insight. The consensus method seems to use bootstrapping. I haven't done this yet so I can't tell you how to implement in lavaan or R. If you or anyone else works out the method, I'd be interested in seeing it. Likewise, if I find the time to do so, I'll share my results here.

Here's the P&D article:
http://digitalcommons.wayne.edu/cgi/viewcontent.cgi?article=1012&context=jmasm

Hellen Geremias dos Santos

unread,
Mar 25, 2014, 5:01:18 PM3/25/14
to lav...@googlegroups.com
Dear Robin

Many thanks for your advice on this. I will try to work out the method, if I do so, I'll share my results here too.

Kind regards, 
Hellen

Matthias Vorberg

unread,
Nov 4, 2020, 4:32:05 PM11/4/20
to lavaan

Hello everyone,

 

I would like to calculate Cronbach's Alpha and McDonald's Omega

 

Alpha I have always calculated with this function:

alpha(data5[c("Fd01","Fd03","Fd04")])

 

Since I additionally calculate McDonald's omega, I thought that I could skip the separate alpha calculation, because the alpha is also automatically given in the omega function.

 

By omega function I mean this one:

 

omegasubsetHa <- subset(data5, select = c("Fd01","Fd03","Fd04"))

omega(omegasubsetHa, nfactors=1)

 

So far, both calculation methods always resulted in the same value for alpha. But this was not the case in my current calculation, and here I got different values. Using the alpha function I got alpha = 0.42. Using the omega function I got alpha = 0.39.

 

My question is now: How can it be that for Cronbach's alpha two different values appear? And, if I should not have made a mistake, which value should I now trust and indicate in my study?

 

Below you find the output of the two analyses.

 

Thank you very much for looking at it!

 

 

---------------------------------------------

 

Output of the alpha function:

alpha(data5[c("Fd01","Fd03","Fd04")])

 

Reliability analysis   

Call: alpha(x = data5[c("Fd01", "Fd03", "Fd04")])

 

  raw_alpha std.alpha G6(smc) average_r  S/N   ase mean   sd median_r

      0.42      0.39     0.4      0.18 0.65 0.059  3.7 0.97    0.042

 

 lower alpha upper     95% confidence boundaries

0.3 0.42 0.53 

 

 Reliability if an item is dropped:

     raw_alpha std.alpha G6(smc) average_r    S/N alpha se var.r  med.r

Fd01     0.079     0.080   0.042     0.042  0.087    0.110    NA  0.042

Fd03    -0.085    -0.086  -0.041    -0.041 -0.079    0.130    NA -0.041

Fd04     0.694     0.694   0.531     0.531  2.268    0.037    NA  0.531

 

 Item statistics 

       n raw.r std.r  r.cor  r.drop mean  sd

Fd01 273  0.76  0.74 0.6057 0.36859  2.8 1.5

Fd03 272  0.80  0.78 0.6705 0.43515  3.5 1.5

Fd04 272  0.45  0.50 0.0062 0.00098  4.6 1.3

 

Non missing response frequency for each item

        1    2    3    4    5    6 miss

Fd01 0.18 0.32 0.22 0.12 0.10 0.07    0

Fd03 0.11 0.19 0.20 0.22 0.18 0.10    0

Fd04 0.02 0.07 0.10 0.18 0.33 0.31    0

 

 

---------------------------------------------------

 

Output of the omega function:

omegasubsetHa <- subset(data5, select = c("Fd01","Fd03","Fd04"))

omega(omegasubsetHa, nfactors=1)

 

 

Omega_h for 1 factor is not meaningful, just omega_t

Omega 

Call: omegah(m = m, nfactors = nfactors, fm = fm, key = key, flip = flip, 

    digits = digits, title = title, sl = sl, labels = labels, 

    plot = plot, n.obs = n.obs, rotate = rotate, Phi = Phi, option = option, 

    covar = covar)

Alpha:                 0.39 

G.6:                   0.41 

Omega Hierarchical:    0.59 

Omega H asymptotic:    1.02 

Omega Total            0.58 

 

Schmid Leiman Factor loadings greater than  0.2 

        g  F1*   h2   u2 p2

Fd01 0.53      0.28 0.72  1

Fd03 1.00      1.00 0.00  1

Fd04           0.00 1.00  1

 

With eigenvalues of:

  g F1* 

1.3 0.0 

 

general/max  Inf   max/min =   NaN

mean percent general =  1    with sd =  0 and cv of  0 

Explained Common Variance of the general factor =  1 

 

The degrees of freedom are 0  and the fit is  0.01 

The number of observations was  273  with Chi Square =  1.7  with prob <  NA

The root mean square of the residuals is  0.03 

The df corrected root mean square of the residuals is  NA

 

Compare this with the adequacy of just a general factor and no group factors

The degrees of freedom for just the general factor are 0  and the fit is  0.01 

The number of observations was  273  with Chi Square =  1.7  with prob <  NA

The root mean square of the residuals is  0.03 

The df corrected root mean square of the residuals is  NA 

 

Measures of factor score adequacy             

                                                 g F1*

Correlation of scores with factors            1.00   0

Multiple R square of scores with factors      1.00   0

Minimum correlation of factor score estimates 0.99  -1

 

 Total, General and Subset omega for each subset

                                                 g  F1*

Omega total for total scores and subscales    0.58 0.59

Omega general for total scores and subscales  0.59 0.59

Omega group for total scores and subscales    0.00 0.00

Warning messages:

1: In schmid(m, nfactors, fm, digits, rotate = rotate, n.obs = n.obs,  :

  Omega_h and Omega_asymptotic are not meaningful with one factor

2: In cov2cor(t(w) %*% r %*% w) :

  diag(.) had 0 or NA entries; non-finite result is doubtful

Patrick (Malone Quantitative)

unread,
Nov 4, 2020, 5:03:18 PM11/4/20
to lav...@googlegroups.com
It looks like omega() (from which package?) is returning standardized
alpha, which is .39 in your alpha() (from which package?) results.
> --
> 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/003b6098-acf5-4fa5-a31f-d452f9c46eeco%40googlegroups.com.



--
Patrick S. Malone, Ph.D., Malone Quantitative
NEW Service Models: http://malonequantitative.com

He/Him/His

Matthias Vorberg

unread,
Nov 5, 2020, 9:32:39 AM11/5/20
to lavaan
Thank you, Patrick!

Both functions alpha() and omega() are from the psych package. Do you think that both functions use a different algorithm for calculating Cronbach's Alpha? What kind of standardization could be in use here?


Patrick (Malone Quantitative)

unread,
Nov 5, 2020, 9:39:27 AM11/5/20
to lav...@googlegroups.com
It just means the variables are standardized before computing alpha.

On Thu, Nov 5, 2020 at 9:32 AM Matthias Vorberg <inform...@gmx.de> wrote:
>
> Thank you, Patrick!
>
> Both functions alpha() and omega() are from the psych package. Do you think that both functions use a different algorithm for calculating Cronbach's Alpha? What kind of standardization could be in use here?
>
>
> --
> 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/667d521f-16ba-4d4f-8f85-e0fd5eeb9b95o%40googlegroups.com.

Patrick (Malone Quantitative)

unread,
Nov 5, 2020, 9:43:04 AM11/5/20
to lav...@googlegroups.com
Or, said differently, it's computed from correlations instead of covariances.

Terrence Jorgensen

unread,
Nov 5, 2020, 3:05:42 PM11/5/20
to lavaan
Bill's omega in the psych package is different from the originally proposed omega(s) implemented in semTools:


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

Matthias Vorberg

unread,
Nov 6, 2020, 5:32:27 AM11/6/20
to lavaan
Thank you, Terrence! Okay, this means that there are two different omegas. But what about alpha (indeed, that was my question), where I get different values when using either the function alpha() or the function omega ()? Patrick suggested that it is calculated once from the covariances and once from the correlations. Or is it perhaps the same with alpha as with omega, that Cronbach's alpha is calculated once and another alpha when using the other function?


Terrence Jorgensen

unread,
Nov 6, 2020, 9:37:29 AM11/6/20
to lavaan
what about alpha

semTools(reliability) and psych::alpha() should yield the same alpha for continuous data:

example(reliability) # alpha for the total (all 9 items) = 0.76
## all these also yield 0.76
alpha
(HolzingerSwineford1939[paste0("x",1:9)])
alpha
(cov(HolzingerSwineford1939[paste0("x",1:9)]))

If you are analyzing ordinal items, then semTools will calculate alpha from the polychoric correlations, whereas psych will treat them as continuous (as it says on the help page, as well as a message printed to the screen whenever running semTools::reliability() on ordinal data)

example(lavTables)
reliability
(fit, return.total = TRUE)
## matches standard calculate on polychorics:
alpha
(lavInspect(fit, "sampstat")$cov)
## as opposed to treating binary data as numeric:
alpha
(HSbinary)

Matthias Vorberg

unread,
Nov 6, 2020, 3:15:37 PM11/6/20
to lavaan
Thank you, Terrence.

I will then re-check my alpha calculation with the function semTools(reliability). I did not know the different treatment as  continuous and ordinal data by the functions.

Reply all
Reply to author
Forward
0 new messages