You are here

Trivariate Cholesky with thresholds -> Mx status RED

6 posts / 0 new
Last post
Charlotte's picture
Offline
Joined: 07/02/2012 - 11:13
Trivariate Cholesky with thresholds -> Mx status RED
AttachmentSize
Binary Data Cholesky_thresholds.R14.62 KB
Binary Data myFunctions.R3.44 KB

Dear all,

I would like to fit a trivariate Cholesky model in OpenMx. I have data of twins aged 7, 10, and 12 years old with longitudinal data for a part of the sample. My data are skewed, so I made three categories and I would like to fit a threshold model on them. I fixed the thresholds to be equal across ages, so I estimate two thresholds in total. The means and the variances of the phenotype increase with age though, so I would like to fix the mean at time point 1 to zero and the variance to one and freely estimate the means and variances of the subsequent ages. There are also quantitative and qualitative sex differences. When I run the script that is attached to this posting, I keep getting the following warning message:

Running CholACE
Warning message:
In model 'CholACE' Optimizer returned a non-zero status code 6. The model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch (Mx status RED)

I have changed the starting values many times and each time I get a different -2LL and totally different path estimates. I thought that there might be too few observations in some of the cells (especially for the longitudinal data), however the model does not run on simulated data either (see script).

There is no problem when I run univariate threshold models for each age separately.

Does anyone know what is going wrong here? Could it be that the variances need to be one for all ages?

Many thanks in advance!

Best regards
Charlotte

neale's picture
Offline
Joined: 07/31/2009 - 15:14
Two things

Hi Charlotte

I'm not sure what's going on, but there are a couple of things to try. One would be to use CSOLNP instead, it's a bit better with ordinal data. The other is to use my preferred trick for ordinal data, which is to fix the first two thresholds to zero and one. Then you can estimate means and variances, just the same as if fitting to continuous data.

CSOLNP does look better:

> mxOption(NULL,"Default optimizer","CSOLNP")
> CholAceFitCSOLNP   <- mxRun(CholAceModel, intervals=F)
Running CholACE
> CholAceFitCSOLNP@output$gradient
       Tm      incm        Tf      incf    am_1_1    am_2_1    am_3_1    am_2_2    am_3_2 
-7.77e-06 -2.92e-06  1.45e-06 -5.32e-06 -2.09e-04  1.29e-06  1.61e-06 -2.06e-06 -1.75e-06 
   am_3_3    cm_1_1    cm_2_1    cm_3_1    cm_2_2    cm_3_2    cm_3_3    em_1_1    em_2_1 
 7.07e-06 -3.40e-04 -3.43e-06  9.77e-08 -1.19e-06  1.40e-07  3.84e-06  5.46e-04 -5.78e-06 
   em_3_1    em_2_2    em_3_2    em_3_3    af_1_1    af_2_1    af_3_1    af_2_2    af_3_2 
-1.79e-06 -1.17e-06 -2.11e-06  6.46e-06 -2.18e-04  1.67e-06  4.49e-07  5.89e-06 -1.86e-06 
   af_3_3    cf_1_1    cf_2_1    cf_3_1    cf_2_2    cf_3_2    cf_3_3    ef_1_1    ef_2_1 
-1.38e-07  6.58e-04 -6.35e-06 -1.66e-06  4.05e-06 -3.57e-06 -6.05e-07 -1.09e-03  3.34e-06 
   ef_3_1    ef_2_2    ef_3_2    ef_3_3   rcdosP1   rcdosP2   rcdosP3 
-5.23e-06  4.67e-06 -1.34e-06  3.53e-06 -1.26e-03 -1.19e-04 -3.36e-05 
Charlotte's picture
Offline
Joined: 07/02/2012 - 11:13
Still not working...

Dear Mike,

Many thanks for your quick reply! I have tried three things now:
1- Adding mxOption(NULL,"Default optimizer","CSOLNP") to the syntax.
2- Fixing the thresholds to -1 and 1 (increment of 2) and freely estimating the means and variances –> I chose -1 and 1 rather than 0 and 1 based on the simulation script.
3- Both 1- and 2- at the same time.

I got these results:
1- The gradient indeed looks much better. However, the estimates change dramatically when I change the starting values. Actually, the estimates equal the starting values apart from the ones that are related to the variance that I have fixed to one (a11, c11, e11). Standard errors are not included in the output (due to the statement fixing the variance of the first measurement to one, I guess).
2- It took very long for this model to converge. Many standard errors are not estimated (“NA”). The gradient doesn’t look good either (some >10). The estimates are not equal to the starting values this time. However, the estimates and the fit change even when the starting values stay the same.
3- For option 3, the gradient gets really, really bad (e.g., -14213.93885). The estimates are not on the starting values. Again, the estimates and the fit change even when the starting values stay the same.

For all three approached, I get a Mx status Red-warning.

I don’t know what it means, but I also get this:
[0] here in findMax

I have attached the three scripts. Any advice is very much appreciated!

Best regards
Charlotte

Charlotte's picture
Offline
Joined: 07/02/2012 - 11:13
Any ideas?

Dear all,

I am still stuck with this problem... Does anyone have a suggestion on what else I could try?

Best
Charlotte

neale's picture
Offline
Joined: 07/31/2009 - 15:14
Integration precision

Hi Charlotte

A couple of things:

  1. mxTryHard() - check out the arguments to this function and tune it to your application (the default for scale may be a bit large). I got:
> CholAceHard@output$gradient
           Tm          incm            Tf          incf        am_1_1        am_2_1        am_3_1 
-2.620126e-07 -7.971401e-07  4.074519e-07 -9.769963e-08  2.605360e-05 -1.265654e-07 -2.753353e-07 
       am_2_2        am_3_2        am_3_3        cm_1_1        cm_2_1        cm_3_1        cm_2_2 
 1.754152e-07  7.771561e-09  3.697043e-07 -2.787326e-05  2.109424e-08 -8.881784e-08  1.032507e-07 
       cm_3_2        cm_3_3        em_1_1        em_2_1        em_3_1        em_2_2        em_3_2 
 4.551914e-08  3.874678e-07  3.296252e-06  7.771561e-09  7.882583e-08 -5.551115e-08  2.220446e-08 
       em_3_3        af_1_1        af_2_1        af_3_1        af_2_2        af_3_2        af_3_3 
 2.509104e-07  2.262301e-05  1.532108e-07  3.630429e-07 -3.397282e-07 -4.440892e-08  9.880985e-08 
       cf_1_1        cf_2_1        cf_3_1        cf_2_2        cf_3_2        cf_3_3        ef_1_1 
-8.655299e-06  2.364775e-07  1.254552e-07 -2.664535e-07 -1.287859e-07  1.798561e-07 -1.315725e-05 
       ef_2_1        ef_3_1        ef_2_2        ef_3_2        ef_3_3       rcdosP1       rcdosP2 
 2.375877e-07  7.105427e-08 -2.164935e-07  7.771561e-09  2.109424e-07 -3.202327e-05 -2.648992e-06 
      rcdosP3 
-1.071032e-05 
  1. Increase the integration precision (see http://openmx.psyc.virginia.edu/thread/3868 or options()$mxOption )

I don't think that fixing the first two thresholds to -1 and 1 or 0 and 1 would make that much difference, except that where you start the mean and what the starting values generate as variance.

In my experience, starting parameters where the covariance matrix is diagonal (i.e. all E) can be a good strategy. It keeps things nicely positive definite and it is easy to ensure that there are no wild outliers (the likelihood is the product of the univariate likelihoods). It is always a good idea to choose reasonable starting values. Optimization isn't an exact science, so very poor starting values are likely to cause problems.

If you want SE's, it is necessary to eliminate the mxConstraint() statements. Or use bootstrapping. Eliminating mxConstraints can usually be done by fixing a diagonal E parameter or two to (e.g.) 1. Be careful to avoid bounds on parameters as well (this can interfere with things if you are working in unstandardized space).

Finally, I don't know what the Find warning message means, but will ask the optimizer coders.

RobK's picture
Offline
Joined: 04/19/2011 - 21:00
The message "here in findMax"

The message "here in findMax" is mxLog output. I think it happens when CSOLNP encounters a non-finite fitfunction value.