Confirmatory factor analysis problems using sem package (works in Amos)
Dear Solomon,
-----Original Message----- From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
On
Behalf Of S. Messing Sent: May-22-09 1:27 AM To: r-help at r-project.org Subject: [R] Confirmatory factor analysis problems using sem package
(works
in Amos) Hello all, I'm trying to replicate a confirmatory factor analysis done in Amos.
Occasionally in an ill-conditioned problem, one program will produce a solution and another won't. As a general matter, I'd expect Amos to be more robust than sem() since Amos is written specifically for SEMs, while sem() uses nlm(), a general-purpose optimizer.
The idea is to compare a one-factor and a two-factor model. I get the
following
warning message when I run either model: "Could not compute QR decomposition of Hessian. Optimization probably did not converge." I have no idea what to do here.
A general strategy is to set debug=TRUE in the call to sem() and see what happens in the optimization. Then there are several things that you can do to try to get the optimization to converge; see ?sem. In this case, however, I wasn't able to get a solution. The one-factor model is equivalent to a one-factor exploratory FA, which can be fit by ML using factanal():
factanal(factors=1, covmat=correl, n.obs=1100)
Call:
factanal(factors = 1, covmat = correl, n.obs = 1100)
Uniquenesses:
pvote jmposaff jmnegaff boposaff bonegaff
obama.therm mccain.therm oddcand.D evencand.D
0.100 0.496 0.497 0.277 0.397
0.129 0.312 0.466 0.585
Loadings:
Factor1
pvote -0.949
jmposaff 0.710
jmnegaff -0.709
boposaff -0.850
bonegaff 0.777
obama.therm -0.934
mccain.therm 0.830
oddcand.D 0.731
evencand.D 0.645
Factor1
SS loadings 5.744
Proportion Var 0.638
Test of the hypothesis that 1 factor is sufficient.
The chi square statistic is 1710.03 on 27 degrees of freedom.
The p-value is 0
As you can see, the one-factor model fits the data very poorly (as does a
two-factor EFA); I suspect, but am not sure, that this is the source of the
problem in sem(). I couldn't get a solution from sem() even when I used the
factanal() solution as start values.
I believe posters reported the same problem.
In almost all cases, the models haven't been properly specified, which is not the case here. Here, the model just doesn't fit the data.
It seems that the ability to invert the correlation matrix may have something to do with this error, but solve(correl) yields a solution.
No, the input correlation matrix is positive-definite. sem() would have complained if it were not:
eigen(correl, only.values=TRUE)
$values [1] 6.12561630 0.82418329 0.71616585 0.51263750 0.24467315 0.18248909 0.17024374 [8] 0.13905585 0.08493524 I'll keep your problem as a test case to see whether I can produce a solution, possibly using a different optimizer -- as I mentioned, sem() uses nlm(). Regards, John
Here are my correlation matrix and model specifications: --------------------------- R CODE BEGIN ------------------------------------------------ library(sem) correl <- matrix( c(1.0000000,-0.6657822,0.6702089,0.7997673,-0.7225454,0.8992372, -0.8026491,-0.6715168,-0.5781565,- 0.6657822,1.0000000,-0.5107568, -0.5030886,0.6971188,- 0.6306937,0.7200848,0.5121227,0.4496810, 0.6702089,-0.5107568,1.0000000,0.7276558,- 0.3893792,0.6043672, -0.7176532,-0.5247434,-0.4670362,0.7997673,- 0.5030886,0.7276558, 1.0000000,-0.6251056,0.8164190,-0.6728515,- 0.6371453,-0.5531964, -0.7225454,0.6971188,-0.3893792,- 0.6251056,1.0000000,-0.7760765, 0.6175124,0.5567924,0.4914176,0.8992372,- 0.6306937,0.6043672, 0.8164190,-0.7760765,1.0000000,-0.7315507,- 0.6625136,-0.5814590, -0.8026491,0.7200848,-0.7176532,- 0.6728515,0.6175124,-0.7315507,
1.0000000,0.5910688,0.5096898,-0.6715168,0.5121227,-
0.5247434, -0.6371453,0.5567924,- 0.6625136,0.5910688,1.0000000,0.8106496, -0.5781565,0.4496810,-0.4670362,- 0.5531964,0.4914176,-0.5814590, 0.5096898,0.8106496,1.0000000),
,nrow=9,ncol=9)
rownames(correl) = c("pvote", "jmposaff", "jmnegaff",
"boposaff","bonegaff",
"obama.therm", "mccain.therm",
"oddcand.D", "evencand.D" )
colnames(correl) = c("pvote", "jmposaff", "jmnegaff",
"boposaff","bonegaff",
"obama.therm", "mccain.therm",
"oddcand.D", "evencand.D" )
#One Factor Model:
model.all <- specify.model()
allmeasures -> pvote, b1, NA
allmeasures -> obama.therm, b2, NA
allmeasures -> mccain.therm, b3, NA
allmeasures -> jmposaff, b4, NA
allmeasures -> jmnegaff, b5, NA
allmeasures -> boposaff, b6, NA
allmeasures -> bonegaff, b7, NA
allmeasures -> oddcand.D, b8, NA
allmeasures -> evencand.D, b9, NA
allmeasures <-> allmeasures, NA,1
pvote <-> pvote, v1, NA
obama.therm <-> obama.therm, v2, NA
mccain.therm <-> mccain.therm, v3, NA
jmposaff <-> jmposaff, v4, NA
jmnegaff <-> jmnegaff, v5, NA
boposaff <-> boposaff, v6, NA
bonegaff <-> bonegaff, v7, NA
oddcand.D <-> oddcand.D, v8, NA
evencand.D <-> evencand.D, v9, NA
sem.all <- sem(model.all, correl, 1100)
summary(sem.all)
#Two Factor Model:
model.vi <- specify.model()
verbal -> pvote, b1, NA
verbal -> obama.therm, b2, NA
verbal -> mccain.therm, b3, NA
verbal -> jmposaff, b4, NA
verbal -> jmnegaff, b5, NA
verbal -> boposaff, b6, NA
verbal -> bonegaff, b7, NA
imp -> oddcand.D, b8, NA
imp -> evencand.D, b9, NA
imp <-> imp, NA, 1
verbal <-> verbal, NA, 1
pvote <-> pvote, v1, NA
obama.therm <-> obama.therm, v2, NA
mccain.therm <-> mccain.therm, v3, NA
jmposaff <-> jmposaff, v4, NA
jmnegaff <-> jmnegaff, v5, NA
boposaff <-> boposaff, v6, NA
bonegaff <-> bonegaff, v7, NA
oddcand.D <-> oddcand.D, v8, NA
evencand.D <-> evencand.D, v9, NA
imp <-> verbal, civ, NA
sem.vi <- sem(model.vi, correl, 1100)
summary(sem.vi)
--------------------------- R CODE END
------------------------------------------------
Thanks very much.
-Solomon
--
View this message in context: http://www.nabble.com/Confirmatory-factor-
analysis-problems-using-sem-package-%28works-in-Amos%29-
tp23664618p23664618.html
Sent from the R help mailing list archive at Nabble.com.
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.