curious lack of warning/error
Posted on
carey
Joined: 10/19/2009
| Attachment | Size |
|---|---|
| curiousProblem.R | 1.99 KB |
| dz.csv | 1.58 KB |
Forums
In fitting a bivariate threshold model using a definition variable, OpenMx converged but the final result gave the following predicted covariance matrix
InitT1 InitT2
InitT1 1.00000 1.48088
InitT2 1.48088 1.00000
Obviously not positive definite. Yet I did not get any warnings or errors.
Script and data attached.
Greg
PS Could the problem be in
PS Could the problem be in omxMnor? Add the following line to the code in my previous post:
omxMnor(preCov, c(0, 0), c(-Inf, -Inf), c(t, t))
Log in or register to post comments
In reply to PS Could the problem be in by carey
It's a bug!
I believe this is a bug. omxMnor() seems to know what to do with correlation matrices only. I don't know what it is doing otherwise. Nor am I sure what it is doing with the non-pd input.
> omxMnor(matrix(c(1,0,0,1),2,2), c(0, 0), c(-Inf, -Inf), c(0,0))
[,1]
[1,] 0.25
> omxMnor(matrix(c(1,0,0,1),2,2), c(0, 0), c(-Inf, -Inf), c(1.282,1.282))
[,1]
[1,] 0.8101416
> omxMnor(matrix(c(1,90,90,1),2,2), c(0, 0), c(-Inf, -Inf), c(1.282,1.282))
[,1]
[1,] 0.5056822
> omxMnor(matrix(c(1,90,1,90),2,2), c(0, 0), c(-Inf, -Inf), c(1.282,1.282))
[,1]
[1,] 0.5370363
Log in or register to post comments
In reply to It's a bug! by neale
Well maybe it's a bug
The following do check out ok against classic Mx. I suspect that a) OpenMx only pays attention to the triangle below the diagonal of the covariance matrix (second example); and b) that it does something to curtail out-of-range values for the correlation so that optimization doesn't come to a grinding halt. I've logged a bug report: http://openmx.psyc.virginia.edu/issue/2013/02/omxmnor-problem
and also fixed a couple of minor errors in the example script (reading & mxFactoring the data)
> omxMnor(matrix(c(1,.5,.5,1),2,2), c(0, 0), c(-Inf, -Inf), c(0,0))
[,1]
[1,] 0.3333333
> omxMnor(matrix(c(1,.3,.5,1),2,2), c(0, 0), c(-Inf, -Inf), c(1, 1))
[,1]
[1,] 0.7281473
> omxMnor(matrix(c(1,.3,.3,1),2,2), c(0, 0), c(-Inf, -Inf), c(1, 1))
[,1]
[1,] 0.7281473
I don't see anything wrong with the standardization code (in omxAlgebraFunctions.c) called before the call to SADMVN to get the integral:
void omxStandardizeCovMatrix(omxMatrix* cov, double* corList, double* weights) {
// Maybe coerce this into an algebra or sequence of algebras?
if(OMX_DEBUG) { Rprintf("Standardizing matrix."); }
int rows = cov->rows;
for(int i = 0; i < rows; i++) {
weights[i] = sqrt(omxMatrixElement(cov, i, i));
}
for(int i = 0; i < rows; i++) {
for(int j = 0; j < i; j++) {
corList[((i*(i-1))/2) + j] = omxMatrixElement(cov, i, j) / (weights[i] * weights[j]);
}
}
}
Log in or register to post comments
In reply to It's a bug! by neale
Any fix soon? Bet many users
Any fix soon? Bet many users struggling with IFAIL=6 messages for threshold models will appreciate it while (shudder) some with models that converged may have to redo their (shudder even more, published) results.
In the interim, suggest using constraints to make certain the matrix is pd.
Log in or register to post comments