Confidence intervals of ACE model fit
Posted on
ywsky01
Joined: 08/06/2020
| Attachment | Size |
|---|---|
| The principal code used for bivariate ACE model | 36.64 KB |
| The code used to carry out batch process | 5.29 KB |
| table output | 38.33 KB |
Hi all,
I'm new with OpenMx. I have been trying to use a code obtained from others to carry out a bivariate ACE model. However, I encountered a situation that most of the confidence intervals of the model fit are NA. It seems this topic has been discussed thoroughly on this forum. But I still cannot confirm my problem by myself. I tried switching the optimizer to “SLSQP” but didn't solve the problem. In addition to mxCI, I found there are several alternatives, discussed on this forum, to estimate confidence intervals. Before trying the other options, I want to make sure if it is okay to proceed without confidence intervals. As discussed at (https://openmx.ssri.psu.edu/node/4269), there might be three reasons for the missing value of confidence intervals. I wonder if the fact that the confidence interval is too narrow means the estimate is okay? Or how can I set the limit to verify it? As suggested in the abovementioned discussion, I cannot carry out the analysis. So wish somebody can help me to examine the problem. Thanks a lot.
I'm new with OpenMx. I have been trying to use a code obtained from others to carry out a bivariate ACE model. However, I encountered a situation that most of the confidence intervals of the model fit are NA. It seems this topic has been discussed thoroughly on this forum. But I still cannot confirm my problem by myself. I tried switching the optimizer to “SLSQP” but didn't solve the problem. In addition to mxCI, I found there are several alternatives, discussed on this forum, to estimate confidence intervals. Before trying the other options, I want to make sure if it is okay to proceed without confidence intervals. As discussed at (https://openmx.ssri.psu.edu/node/4269), there might be three reasons for the missing value of confidence intervals. I wonder if the fact that the confidence interval is too narrow means the estimate is okay? Or how can I set the limit to verify it? As suggested in the abovementioned discussion, I cannot carry out the analysis. So wish somebody can help me to examine the problem. Thanks a lot.
ps: I attach the code I used to this post.
**Below is my mxVersion:**
mxVersion()
OpenMx version: 2.17.3 [GIT v2.17.3]
R version: R version 3.6.2 (2019-12-12)
Platform: x86_64-w64-mingw32
Default optimizer: CSOLNP
NPSOL-enabled?: No
OpenMP-enabled?: No
**The Results I got:**
summary(CholAeFit,verbose=T)
Summary of CholAE
see attachment "ace output.r"
CI details; generic advice
modelsumm <- summary(myModel)
, you can extract the CI details table from it. Save it as a text file or .csv, and post it as an attachment. Its first 7 columns are the most important.
In the meantime, I can offer some generic advice.
First, a lot of your confidence limits have diagnostic "active box constraint". A simple thing you can do about that is to get rid of the `lbound` argument when you construct `pathA`, `pathC`, and `pathE`.
Second, you may want to consider using an alternate parameterization of your models--the "direct-symmetric" parameterization. For example, see 'twoACE5c25raFac.R' on this page. There is a lot of prior discussion on these forums about the two parameterizations. You currently use the "Cholesky" parameterization, in which variance components are the squares of one-headed path coefficients; the alternative is the "direct-symmetric" parameterization, under which the variance components are themselves free parameters of the model. There's a good chance that switching parameterizations will eliminate most of your problematic confidence limits, but your results may be much less interpretable.
Third, one of your models got a status RED. Try replacing `mxRun()` with
mxTryHard()for at least that model. Note that `mxTryHard()` behaves stochastically, so you might want to use `set.seed()` to make your results reproducible.Finally, do you know about the 'umx' package? You might find it really useful!
Log in or register to post comments
In reply to CI details; generic advice by AdminRobK
Thanks & further questions
I got rid of the lbound argument of pathA, pathC, and pathE, it helped me solve most of the cases without a value for the confidence intervals. However, there are still some confidence limits left as NA.
I tried mxTryHard().
I have tried to save the summary results to a text file without confirming if I performed the procedure correctly.
I have been trying to understand direct-symmetric & switching parameterizations. However, it is difficult for me to understand it thoroughly in a short time.
The reference to umx also opened a new area for me to study. I am still learning it.
I attached the results to this post. If possible, could you help me to take a further examination? I want to know if it is okay to proceed with current results.
Thanks for the informative comments. I am grateful for that.
Log in or register to post comments
In reply to Thanks & further questions by ywsky01
write.csv()
modelsumm <- summary(myModel)
cid <- modelsumm$CIdetails
write.csv(cid,"CSOLNP.csv",quote=TRUE)
That will create CSOLNP.csv in R's working directory.
Log in or register to post comments
In reply to write.csv() by AdminRobK
without CIdetails
I have CholAceModel &
< CholAceFit <- mxRun(CholAceModel, intervals = TRUE)>,
Performing
summary(CholAceModel, verbose=T) , or
summary(CholAceFit), $$
turns out that the model$CIdetails is NULL,
while $CI shows the values of lbound, estimate, ubound and note for each parameter.
Log in or register to post comments
In reply to without CIdetails by ywsky01
CIdetail
Log in or register to post comments
In reply to CIdetail by AdminRobK
$CIdetail
I attached the $CIdetail to this post. The name of the file indicates which optimizer was used to produce the results. Wish you can help me to figure out the problems. Thanks!
Log in or register to post comments
In reply to $CIdetail by ywsky01
Prop*
I don't know what to make of the off-diagonal 'Prop*' elements. Both optimizers seem to really have difficulty with them (are you able to install an NPSOL-enabled build of OpenMx on your system?). If you really need confidence intervals for those elements, you might need to use their standard errors to do so. You can get standard errors for them with
mxSE(), but you should run your model throughimxRobustSE()first if you can.Alternately, you could bootstrap.
Log in or register to post comments
In reply to Prop* by AdminRobK
Things going on interesting
I cannot install the optimizer NPSOL. At first, it seemed my R cannot connect to (http://openmx.ssri.psu.edu/getOpenMx.R). Actually, it is true. However, I tried to copy the command on the internet page to R command window to run it, it turned out that the OpenMx is not available (perhaps for the R version 3.6.2 or 3.5.2).
I am now studying SEM with youtube videos. It is really a pleasure to learn more about it. I will try the methods you mentioned, mxSE and bootstrap, to estimate the Prop* parameter.
It is really a pleasant environment giving kind messages and detailed suggestions.
All of a sudden, a question was coming out. Given the bivariate model I used, I wonder if the model only estimates the path from latent variable A1 to the observed measurement 1 (such as structural connectivity), the path from latent variable A1 to the observed measurement 2 ( such as functional connectivity), and the path from latent variable 2 to the observed measurement 2, but the path from latent variable A2 to the observed measurement 1?--This question came out of my mind while I was watching a video to understand Cholesky bivariate model.
Log in or register to post comments
In reply to Things going on interesting by ywsky01
Cholesky
That path is not estimated under the Cholesky parameterization. The model would be unidentified if that path were a free parameter.
Log in or register to post comments
In reply to CIdetail by AdminRobK
Another question-about model fitting abnormally exited
Log in or register to post comments
In reply to Another question-about model fitting abnormally exited by ywsky01
NaN legal in data?
Possibly that should be caught by OpenMx as invalid data rather than running for hours...?
Log in or register to post comments
In reply to NaN legal in data? by tbates
Thanks you!
Since I have a large dataset to run the bivariate model, it might be the reason for running hours. I need to run 30,135 bivariate models for 30, 135 pairs of features.
Log in or register to post comments
In reply to Another question-about model fitting abnormally exited by ywsky01
abnormal termination
Are you calculating confidence intervals for the ACE, AE, and CE models as in the script in your OP? There is a known issue concerning CI calculations being unnecessarily slow. Also, you can use
omxRunCI()to calculate CIs for a model you've already run.Since you're trying to automate running all your models as a batch, you should make your code more error-tolerant. For example:
CholCeFit <- try(mxRun(CholCeModel, intervals = TRUE))
if( !("try-error" %in% class(CholCeFit)) ){
# Code to cxtract results from CE model here.
}
Also, consider replacing `mxRun()` with `mxTryHard()`.
I doubt the presence of `NaN` in your data explains why your CE model terminated abnormally. Instead, what I think is happening is that, once you drop _A_ from the ACE model, the model-expected DZ covariance matrix is computationally singular, and the fitfunction is evaluated as non-finite. The optimizer can cope with non-finite fit in the middle of optimization, but not right at the beginning. You will need to tweak the values of the free parameters for your CE model before running it. Or, again, use `mxTryHard()` in place of `mxRun()`.
Log in or register to post comments
In reply to abnormal termination by AdminRobK
how to tweak the values of free parameters
Log in or register to post comments
In reply to CIdetail by AdminRobK
Generic advice
Log in or register to post comments
In reply to Generic advice by ywsky01
SEM & umx resources
1. The umx paper. Bates, T. C., Maes, H., & Neale, M. C. (2019). umx: Twin and Path-Based Structural Equation Modeling in R. Twin Res Hum Genet, 22(1), 27-41. doi:[10.1017/thg.2019.2](https://www.cambridge.org/core/journals/twin-research-and-human-genetics/article/umx-twin-and-pathbased-structural-equation-modeling-in-r/B9658AC0CDA139E540BFAC0C9D989623)
2. Examples in `?umx` package.
3. Online [umx tutorial](https://tbates.github.io) at https://tbates.github.io If anything needs improving, let me know: feedback good bad or ugly greatly valued!)
Log in or register to post comments
In reply to SEM & umx resources by tbates
Thanks a lot. I quite like
Log in or register to post comments
In reply to Generic advice by ywsky01
See Boulder Workshop Materials
The Boulder workshops have a lot of material, including presentations and accompanying scripts and exercise materials. A reasonable place to start is https://www.colorado.edu/ibg/international-workshop/2020-international-statistical-genetics-workshop/workshop-2020-file
Log in or register to post comments
In reply to See Boulder Workshop Materials by AdminNeale
Appreciated
Thanks for the sharing. I am trying to google if there are related videos but haven't found any. I will make use of the materials after some video study, which might help me to better understand the materials.
Log in or register to post comments