I fit a nuclear twin family model. AE model is the best. However, the CI of e makes me confused. Here are the output in free parameters and CI details section. The CI details looks strange, ranging from 0 to 2.053589e+06.
name matrix row col Estimate Std.Error A lbound ubound
1 AddGenPath MZntf.a 1 1 0.708885615 1.004657e-02
2 EnvPath MZntf.e 1 1 0.668286218 1.224575e-02
3 AMCopath MZntf.mu 1 1 0.590740014 1.836063e-02
CI details:
parameter side value fit diagnostic statusCode method
AddGenPath lower 5.859355e-01 11404.18 success OK neale-miller-1997
AddGenPath upper 8.259182e-01 11404.18 success OK neale-miller-1997
EnvPath lower -5.363198e+00 11403.31 alpha level not reached infeasible non-linear constraint neale-miller-1997
EnvPath upper 2.053589e+06 11402.99 alpha level not reached infeasible non-linear constraint neale-miller-1997
Then I add mxBootstrap(AE.NTF.Fit, 1000) to calculate CIs. Here are the results. It looks confusing that the std.error of e is 0.000000 and the CI is so close to the estimate.
name matrix. row col Estimate Std.Error 2.5% 97.5%
1 AddGenPath MZntf.a 1 1 0.708885615 0.06318805 0.56058213 0.81236401
2 EnvPath MZntf.e 1 1 0.668286218 0.00000000 0.66828622 0.66828622
3 AMCopath MZntf.mu 1 1 0.590740014 0.02719992 0.54368092 0.64904098
I would like to know what caused this. Which one is right? Is it acceptable to use CIs in bootstrap method? If not, what can I do to get the accurate CIs?
Printer-friendly version
I agree that those results seem odd. However, I'd need to see your full script to be able to tell what's going on.
Thanks for your reply. I have been working on it for many days but still could not find the solutions. Here is my script.
Which optimizer are you using? If you're using CSOLNP (the on-load default), then you should try either (1) switching optimizers or (2) placing the MxConstraints in either the MZ or DZ MxModel, but not both. You can switch optimizers with
mxOption(), e.g.mxOption(NULL,"Default optimizer","SLSQP")switches optimizers to SLSQP.What do you get from
qr(ADSFE.NTF.Fit$output$constraintJacobian)?What do you get from
mxVersion()?Aside from the confidence intervals...do the point estimates you get look reasonable?
How many phenotypes are you analyzing, anyhow? Your posted script doesn't show what value is assigned to
nv.Thanks for your quick reply. My openMx version is 2.14.11. I am analyzing one phenotype and nv<-1 in my script. Aside from the CIs, the point estimate looks reasonable. The estimate of E is 0.662516155 and the std.error is 1.503494e-02.
The output I get from
qr(ADSFE.NTF.Fit$output$constraintJacobian)$qr AddGenPath DomPath AMCopath FamilialPath VarPhen SibPath MZntf.VpCon[1,1] 2.0110369 1.0050667 0.000000000 0.00000000 -9.241115e-01 5.415512e-01 MZntf.deltaCon[1,1] 0.5352611 -1.1643035 0.000000000 0.00000000 1.070522e+00 -6.273514e-01 MZntf.qCon[1,1] 0.0000000 0.0000000 0.305701843 0.04914830 3.202266e-04 0.000000e+00 MZntf.xCon[1,1] 0.0000000 0.0000000 0.003955549 0.95232837 6.771227e-03 0.000000e+00 MZntf.wCon[1,1] 0.0000000 0.0000000 0.037370059 0.66355074 -2.248940e-03 1.233176e-29 DZntf.VpCon[1,1] 0.4620558 -0.3661017 0.000000000 0.00000000 7.404974e-14 1.755417e-16 DZntf.deltaCon[1,1] 0.5352611 0.6580158 0.000000000 0.00000000 -4.936649e-14 -9.486833e-01 DZntf.qCon[1,1] 0.0000000 0.0000000 0.706107522 -0.02135926 -1.677573e-02 7.375195e-16 DZntf.xCon[1,1] 0.0000000 0.0000000 0.003955549 0.24386946 -5.335934e-01 2.345862e-14 DZntf.wCon[1,1] 0.0000000 0.0000000 0.037370059 0.66355074 5.979116e-01 -2.628627e-14 EnvPath CovPhenGen LatentVarAddGen VarF CovFA mean MZntf.VpCon[1,1] 1.224478 -1.070522e+00 6.080053e-01 9.241115e-01 1.842139e+00 0 MZntf.deltaCon[1,1] -1.418477 -9.241115e-01 1.992176e-01 -1.070522e+00 3.024421e-02 0 MZntf.qCon[1,1] 0.000000 4.670787e-01 -1.412215e+00 -7.911097e-03 -7.474012e-02 0 MZntf.xCon[1,1] 0.000000 1.940195e-02 7.288240e-02 -4.875699e-01 -1.325505e+00 0 MZntf.wCon[1,1] 0.000000 9.909658e-03 1.835284e-02 1.327484e+00 -4.872886e-01 0 DZntf.VpCon[1,1] 0.000000 -5.206094e-16 -1.337556e-15 -9.300964e-14 3.421490e-14 0 DZntf.deltaCon[1,1] 0.000000 3.659091e-16 1.059826e-15 7.249687e-14 -2.655185e-14 0 DZntf.qCon[1,1] 0.000000 -1.479618e-17 -1.024771e-17 -2.017033e-16 1.122347e-16 0 DZntf.xCon[1,1] 0.000000 3.517239e-01 -9.219073e-18 -4.955415e-16 2.676370e-16 0 DZntf.wCon[1,1] 0.000000 5.862064e-02 1.338068e-01 -1.855310e-16 -7.219559e-17 0 $rank [1] 5 $qraux [1] 1.462056e+00 1.658016e+00 1.706108e+00 1.243869e+00 1.597912e+00 1.316228e+00 0.000000e+00 [8] 1.934267e+00 1.991007e+00 7.249711e-14 2.655195e-14 0.000000e+00 $pivot [1] 1 2 5 6 7 3 4 8 9 10 11 12 attr(,"class") [1] "qr"The output of free parameters in AE model
The beginning of the mentioned script
Here is the output for mxVersion().
Finally, I tried to place the MxConstraints only in MZ model and in DZ model by the following correction. These two results are consistent and the CIs could be estimated. However the result changes a lot compared to placing the MxConstraints both MZ and DZ MxModel. The best model is ADSE model and the CIs of a and s are including zero.
The code I changed
the results of ADSE model
Why are these two results so different? I am looking forward to your reply. Thanks forward.
I think the model might not be identified. What do you get if you run
mxCheckIdentification()on your MxModel objects?Also, I notice that you're not adjusting for age of parents and twins, nor for sex of twins. Doing so is generally advised. But, worry about that after trying
mxCheckIdentification().For that matter, your script also models parents and twins as having the same model-expected mean:
Is that really realistic for the trait you're studying?
You are right. The ADSFE model could not be identified. However, the submodels seemed to be identified. The output for the submodels is just like this. The status is “TRUE” and $non_identified_parameters is "None". Could it be seen as identified? If it is identified, is it acceptable that I set some path to zero to make submodels on a non-identified model? If not, should I model ADFE, ADSE and ASFE separately?
Also, the question I mentioned earlier still make me confused. Why I move constraints to MZmodel make such a difference?
You are right that it seems not feasible to model parents and twins as having the same model-expected mean. However it is always to do so in nuclear twin family study. I am new in this area and also have some confusion. I will consider it seriously and control sex and age factors after dealing with the identification problem.
I am looking forward to your reply.
Yes, the model is at least locally identified, at the current values of its free parameters.
Yes, that sounds reasonable to me. It's analogous to how, in the classical twin design, the ACDE model is unidentified, and therefore one must fix either the C or D component to zero in any given model.
That was probably an artifact of model unidentification.
Thanks for your reply. I still have some problems with the result. The path diagram of nuclear twin model is in the attachment.
Firstly, the estimates of m and w are significantly, but NEGATIVE. It seems difficult to explain it reasonably. I don’t know if I should report these negative values as the output shows or I should set these factors between 0 and 1 ahead of running the script.
Secondly, I tried to make ADFE, ASFE and ADSE model separately. However, the ADSE model do not have an equal degree of freedom as ASFE and ADFE model because ADSE model has no estimates of m, w and some nonlinear constraints (i.e., w1=w2, x1=x2). Therefore, I am not sure if it is acceptable for me to do this.
Your reply is highly appreciated.
Why is it a problem if they're negative? You haven't mentioned what trait you're studying.
I don't see that there's anything wrong in doing that.