You are here

SLSQP Confidence Intervals?

5 posts / 0 new
Last post
rabil's picture
Offline
Joined: 01/14/2010 - 16:47
SLSQP Confidence Intervals?
AttachmentSize
Plain text icon comparison_of_npsol_slsqp_lavaan.txt3.98 KB

I ran a very simple one-factor with 3 indicators measurement error model. I used the default optimizer (SLSQP) and then the NPSOL optimizer that I've been using for a long time. Either NPSOL has been outputting bad intervals for a long time or SLSQP is not working very well for confidence intervals. Please see the attached text file that compares SLSQP, NPSOL, and lavaan results for the same model and same data. SLSQP does not compute some bounds, even though NPSOL has no problem computing bounds. (Note that the results include factor loadings, intercepts and variances. The OpenMx results have sds rather than variances while the lavaan results show the variances.)

I also needed to re-run a much more complex model using OpenMx 2.2.6 and using the default optimizer SLSQP. It ran much more quickly than under 1.3 with NPSOL - but many intervals were not computed and the bounds that were computed did not match very well (the estimates appeared to be very similar).

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
additional diagnostics

If you use summary(model, verbose=TRUE) then you can obtain additional diagnostics for the CIs. NPSOL seems to work better because we don't even bother to check whether NPSOL obtained the correct deviance target for the CI (if often does not).

mhunter's picture
Offline
Joined: 07/31/2009 - 15:26
I'm sorry you're experiencing

I'm sorry you're experiencing some trouble here. One of the things the OpenMx development team has been finding is that a feeling of certainty about ones results is often a consequence of limited options. That is, if you only have one optimizer then you think your optimizer works fine. Once you have options, idiosyncrasies can show.

One thing worth mentioning about the lavaan results is their confidence intervals are based purely on the standard errors. Many SEM programs do this. OpenMx does not. The confidence intervals provided by OpenMx are profile likelihood confidence intervals. They are not necessarily symmetric, and I believe they require less stringent assumptions than SE-based intervals. In theory they should be close to SE-based confidence intervals, but only when all the assumptions are met.

You can get SE-based confidence intervals for OpenMx models with something like

require(OpenMx)
data(demoOneFactor)
manifests <- names(demoOneFactor)
latents <- c("G")
factorModel <- mxModel("One Factor", 
                       type="RAM",
                       manifestVars=manifests, 
                       latentVars=latents,
                       mxPath(from=latents, to=manifests),
                       mxPath(from=manifests, arrows=2),
                       mxPath(from=latents, arrows=2, free=FALSE, values=1.0),
                       mxData(observed=cov(demoOneFactor), type="cov", numObs=500))
factorRun <- mxRun(factorModel)
 
sum.param <- summary(factorRun)$parameters
param <- sum.param[,5]
SE <- sum.param[,6]
res.ci <- param + t(outer(c(-1, 0, 1), qnorm(.975)*SE))
rownames(res.ci) <- sum.param[,1]
colnames(res.ci) <- c('lower', 'estimate', 'upper')
res.ci
#                       lower   estimate      upper
#One Factor.A[1,6] 0.36667556 0.39715240 0.42762924
#One Factor.A[2,6] 0.46792656 0.50366144 0.53939631
#One Factor.A[3,6] 0.53716407 0.57724201 0.61731995
#One Factor.A[4,6] 0.65571296 0.70277421 0.74983547
#One Factor.A[5,6] 0.74397964 0.79625052 0.84852139
#One Factor.S[1,1] 0.03530140 0.04081423 0.04632707
#One Factor.S[2,2] 0.03252078 0.03802004 0.04351931
#One Factor.S[3,3] 0.03464878 0.04082720 0.04700562
#One Factor.S[4,4] 0.03270578 0.03938705 0.04606833
#One Factor.S[5,5] 0.02907724 0.03628709 0.04349694

This could be put into a function.

confidenceIntervals <- function(model){
    sum.param <- summary(model)$parameters
    param <- sum.param[,5]
    SE <- sum.param[,6]
    res.ci <- param + t(outer(c(-1, 0, 1), qnorm(.975)*SE))
    rownames(res.ci) <- sum.param[,1]
    colnames(res.ci) <- c('lower', 'estimate', 'upper')
    res.ci
}
rabil's picture
Offline
Joined: 01/14/2010 - 16:47
Thanks!

Thanks!

rabil's picture
Offline
Joined: 01/14/2010 - 16:47
I realize that lavaan is just

I realize that lavaan is just using a simple approximation and that it can be difficult to compute profile intervals. I guess what concerns me most is that the 3-indicator one factor model is about as simple as it gets and the dataset is very well behaved (although without the constraint on the variances, the one variance would have a negative estimate). Given the likely undue weight put on p-values and even confidence intervals in published research, it seems there should be very prominent disclaimers about their accuracy. Most researchers are not going to take the time to investigate this especially when it would put them at a relative disadvantage.