Hi everyone,

I am trying to do a joint analysis between 2 binary variables and 1 continuous variable I am having some problems to fix the variance for the 2 binary variables equal 1.

I think I should change something here

matUnv <- mxMatrix( type="Unit", nrow=nvo, ncol=1, name="Unv1" )

matOc <- mxMatrix( type="Full", nrow=1, ncol=nv, free=FALSE, values=oc, name="Oc" )

var1 <- mxConstraint( expression=diag2vec(Oc%&%V)==Unv1, name="Var1" )

NVO=2

NV=3

However, I have done several changes and it doesn’t work since in some tries the matrix are not conformable and in others it doesn’t fix the variance to 1. How could I fix the variance for the two binary variables to 1?

Thank you so much in advance.

Juan J Madrid-Valero

Without the context of the full script, I can't really answer your question.

Thanks Rob for your prompt response here you can find the full script:

At least in the script you posted,

`Oc%&%V`

is going to evaluate to a 1x1 matrix, but`matUnv`

is 2x1, and the two sides of the comparator in an MxConstraint have to be of the same dimensions. A simple way to fix the variance of the two binary variables to 1 is to replace those three lines with:You could have also just used two MxConstraint objects (one for each variance).

Since your analysis involves an MxConstraint, I recommend using SLSQP as the optimizer.

You could also avoid using an MxConstraint entirely if you change the parameterization of the model to direct variance components, and identify the scale of the ordinal variables by fixing their nonshared-environmental variances to 1.

Thank you so much Rob,

I didn’t know why the matrix was ( expression=diag2vec(Oc%&%V)) 1x1.

Also, I have run the script and the results are pretty feasible but I got this message “In mxTryHard(model = model, greenOK = greenOK, checkHess = checkHess, : 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 tried tons of starting values but I have not been able to remove the “problem”

By the way I have used the optimizer that you told me and the function MxtryhardOrdinal with extra tries.

Juan J Madrid-Valero

Here the full message:

Retry limit reached

Solution found

Start values from best fit:

0.0873764505129859,-0.00302648728148331,0.00387006292591793,0.446376724057658,0.765174696854374,1.24795735809318,2.27769791884809,-0.336483766998107,-0.275906439822281,0.416526067571617,0.471557647293468,0.00245202168140257,0.202289271880198,0.352135547711807,0.0718736069163078,0.00290227641677533,0.00372693059844336,0.00118080806141676,3.1768820509539,-0.106534386756017,-0.147029899250398,0.76021560181588,0.40208133001671,0.716267708855051

Warning message:

In mxTryHard(model = model, greenOK = greenOK, checkHess = checkHess, :

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)

`Oc`

is 1x3. Operator`%&%`

is the quadratic-product operator (see the documentation for`mxAlgebra()`

). Thus,`Oc%&%V`

is equivalent to`Oc%*%V%*%t(Oc)`

, which is (1x3) times (3x3) times (3x1), which evaluates to a 1x1 product.Which version of OpenMx are you using? As of v2.8.3, status Red should occur less frequently in ordinal analyses than it used to. Nonetheless, status Red is sometimes unavoidable in ordinal analyses, in that it can still occur even when the optimizer has found a local minimum. You could try using CSOLNP instead of SLSQP. If that still doesn't help, then it might be worth reparameterizing to eliminate the MxConstraint, which I alluded to in another post.

Thanks Rob,

I am using R version 3.4.3 and I have tried with the CSOLNP optimizer but it doesn't work either. I am going to try to change the parameterization of the model to direct variance components. I am not sure how I should do it. Could you please give me some guidance or where I could find information about that?

Juan J Madrid-Valero

What about your OpenMx version?

Let's see...one way, involving minimal changes to your script, would be as follows. First, remove the line that creates

`pathE`

. Then, define`covE`

as:Then, remove the lines that define

`matUnv`

,`matOc`

, and`var1`

. Finally, remove`pathE`

,`matUnv`

,`matOc`

, and`var1`

from the definition of`pars`

. What these changes do is identify the latent scale of the ordinal variables by fixing theirnonshared-environmental variance, rather than theirtotal variance, to 1.0. As the MxConstraint is no longer required, use the CSOLNP optimizer.While I'm thinking about it--consider removing the lower bounds for

`pathA`

and`pathE`

. You can put them back in if the optimizer keeps reaching really crazy solutions.Sorry! I am using OpenMx version: 2.8.3 [GIT v2.8.3] ; R version: R version 3.4.3 (2017-11-30).

I have made all the changes and now It works and I don't have a red mx status but I have this message:

Retry limit reached

Solution found

Running final fit, for Hessian and/or standard errors and/or confidence intervals

Warning messages generated from final fit for final fit for Hessian/SEs/CIs

Start values from best fit:

0.086752181736907,-0.00393385954682701,0.00457737378618429,0.481260848285712,0.526771126020844,1.19685352917708,2.29257177590332,-0.392479607148631,-0.315671611477633,0.782950339036923,0.501766203589714,-0.266923949247171,10.0657799256565,-0.474366335537433,-0.565466421429395,0.490967007947968

the resuts are the same that with the mx status red.

Thank you so much!!!

Do you think that these results are reliable?

by the way is it possible to know the influences shared by the three phenotypes (all together), not only two by two phenotypes?

and also I would like to know what is the easiest way to compare two correlations in the saturated model

I have tried creating these objects:

CICORMZF <-mxAlgebra(expression =cov2cor(MZf.covMZf), name="CORMZF")

CICORDZF <-mxAlgebra(expression =cov2cor(DZf.covDZf), name="CORDZF")

CICORMZM <-mxAlgebra(expression =cov2cor(MZm.covMZm), name="CORMZM")

CICORDZM <-mxAlgebra(expression =cov2cor(DZm.covDZm), name="CORDZM")

CICORDZO <-mxAlgebra(expression =cov2cor(DZo.covDZo), name="CORDZO")

rAcompare <- mxAlgebra(cbind(CORDZF[1,2] - CORDZM[1,2], name="rAcompare")

and adding them to the model but it doesn't work.

It is a 5group saturated model and I would like to compare the correlations between males and females.

I have thought to equal the correlations in the model using omxsetparameter but in the model I don't have the correlations for the five groups only in the submodels so I am not sure how to proceed.

Thanks!!!

How do you mean?

Yes, that can be done. Look into the Common-Pathway and Independent-Pathway models. Both of them do what you're describing, but in different ways.

One thing wrong I notice immediately is that there's an unbalanced parenthesis in

. I don't think you need the column-bind, so try changing it to

. Otherwise, that block on syntax looks OK to me, but it's hard to say without seeing it in the context of the full script, and without any specifics beyond "it doesn't work." Since you're referencing the covariance matrices by "modelname-dot-objectname", I assume you're placing these algebras in the "container" MxModel, right?

Sorry! When I said if these results are reliable I wanted to say if with this message:

"Retry limit reached

Solution found

Running final fit, for Hessian and/or standard errors and/or confidence intervals

Warning messages generated from final fit for final fit for Hessian/SEs/CIs

Start values from best fit:

0.086752181736907,-0.00393385954682701,0.00457737378618429,0.481260848285712,0...."

I could be sure that the soltuion is "good" I mean since retry limit have been reached I don't know if I should try other starting values or something else (although there is a solution)

Here the 5 group script I am working with:

I would like to add to the model the correlations for each group (MZM,MZF,DZM,DZF,DZO) I have this correlations in the submodels but not in the model. I mean when I run summary(fit) this values are only in the Confident invervals.

Therefore I can't use omxsetparemateres to equal the correlations and test the model (for example corMZF equal to corMZM)

I think I should create five matrix for the correlations but I don't know how to do this. Also I tried adding CICORMZF, CICORMZM... to the model and submodels but I only have the correlations in the submodels.

Aren't you working with two binary phenotypes? If so, the saturated model still needs to have thresholds, and constraints imposed to identify the scale of the latent liability.

Yes sorry, I am working with two binary variables and one continous. This script is for the continous variable. I am testing if there are sex differences (in the continous variable). It is a 5 group saturated univariate model but I don't know how to add the correlations for the 5 groups to the model (actually I have them in the submodels)

What was the warning generated during the final fit? Did you get status Red or not? If you're really concerned, you could increase the number of extra tries. You could also try my Nelder-Mead optimizer, but I doubt it would work any better than CSOLNP.

Ah, I see. The reason for why that doesn't work is unrelated to whether the correlations are in the submodel or in the container (so, you don't need to create new objects for the correlation matrices; the ones you have in the submodels already will do). The reason for why that doesn't work is that the correlations are not statistical parameters of the MxModel object. Parameters are elements of MxMatrices (or values of MxPaths), but the correlations are elements of MxAlgebras. That is, they are not parameters themselves, but are instead

functionsof parameters (informally called "algebras").`omxSetParameters()`

and the 'free parameters' table in`summary()`

output only deal with parameters, not algebras. In contrast, confidence intervals can be calculated for algebras. MxConstraints can also be imposed on algebras.If you want to fit models in which the correlations are constrained across groups in some way, you basically have two options. One option is to reparameterize your MxModels so that the correlations are in fact parameters rather than algebras. For advice on doing that, see this other thread. The other option, which would involve less work, is to use MxConstraints. For instance, if you want to constrain the MZF and MZM correlations to be equal to one another, make an object like this,

and put it in the container MxModel. I recommend using SLSQP when MxConstraints are involved.

You could instead request confidence intervals on the differences in the correlations. For instance,

, and put them into the container MxModel.

Hi Rob! thank you so much again!

1-I think the warnings are for the CI since I have NA in some of them. I don't have mx status red now.

The message is "Retry limit reached

Solution found

Running final fit, for Hessian and/or standard errors and/or confidence intervals

Warning messages generated from final fit for final fit for Hessian/SEs/CIs

Start values from best fit:

0.086752181736907,-0.00393385954682701,0.00457737378618429,0.48126084828...."

I have tried changing the optimizer and 100 extra tries but I have the same resuls and always retry limit is reached.

In this model we have fixed the nonshared-environmental variances to 1 but why not for A and C?

2-Where can I find information about Common-Pathway and Independent-Pathway models? I need to use another script right?

3- With regards to the correations for the five group univariate model it works perfect now! thank you!!!

If you run

`summary()`

with argument`verbose=TRUE`

, you can see diagnostics for the confidence intervals. Also,`mxTryHardOrdinal()`

runs with default argument`exhaustive=TRUE`

, so by default it will always use up all of its extra tries.You could always look around in the materials from the 2016 Boulder Workshop, e.g. here. See also:

Posthuma, D. (2009). Multivariate genetic analysis. In Kim, Y-K. (Ed.),

Handbook of Behavior Genetics(pp. 47-59). New York: Springer Science+Business Media, LLC.Neale, M.C., & Cardon, L.R. (1992).

Methodology for genetic studies of twins and families. Dordrecht, The Netherlands: Kluwer Academic Publishers