I'm trying to compute the expected covariance matrix using mxGetExpected() with a RAM model that includes several mxMatrix objects whose purpose is facilitate algebraic equating of some parameters. When I try to use mxGetExpected() to get the expected covariance, it returns the error:
Error: The following error occurred while evaluating the subexpression '1/(1 + circumplex.v)' during the evaluation of 'z' in model 'circumplex' : non-conformable arrays
By contrast umx::umxExpCov() returns the matrix as expected.
Would it be possible to adjust mxGetExpected() to accommodate such models?
The code below reproduces the error.
mat <- diag(6) vars = c("R","I","A","S","E","C") rownames(mat) <- vars colnames(mat) <- vars N <- 1000 latentvars <- paste("C",vars,sep="") angleHolders <- paste("P",vars[-1],sep="") angles <- paste("p",vars[-1],sep="") startAngles <- 2 * pi * 2:length(vars) / length(vars) k = length(vars) k1 = length(vars) - 1 test.model <- umxRAM("circumplex", data=mxData(mat,type="cov",numObs=N), umxPath(var = vars, fixedAt=0), # To constrain uniquenesses for each variable to equal, # replace these two lines with commented out following line umxPath(var = latentvars, labels=paste0("v",vars)), # labels="v" mxMatrix("Full", nrow=k, ncol=1, free=TRUE, labels=paste("v",vars,sep=""), lbound=0, name="v"), #umxPath(var = latentvars, labels="v"), # If want to estimate communalities and scaling parameters separately, # replace following three lines with the commented out line below; # always do this if the data are supplied as a covariance matrix with unstandardized data mxAlgebra(sqrt(1/(1+v)), name="z"), #umxPath(latentvars, to=vars, free=FALSE, labels="z[1,1]"), # When variances constrained to equal umxPath(latentvars, to=vars, free=FALSE, labels=paste0("z[",1:k,",1]")), # When variances estimated independently # Estimate factor loadings/angles mxMatrix("Full", nrow=k1, ncol=1, free=TRUE, values=startAngles, labels=angles, name="Angles"), mxAlgebra(cos(Angles), name="COS"), mxAlgebra(sin(Angles), name="SIN"), umxPath("Fc", to=latentvars[1], fixedAt=1, labels=paste0("Fc",latentvars[1])), umxPath("Fc", to=latentvars[-1], free=FALSE, # values=startloadingsC, labels=paste0("COS[",1:k1,",1]")), umxPath("Fs", to=latentvars[-1], free=FALSE, # values=startloadingsS, labels=paste0("SIN[",1:k1,",1]")), umxPath(var = c("Fc","Fs"), freeAt=.5, labels="beta1"), # Estimate variance attributable general factor and circumplex mxAlgebra(1-beta1, name="beta0"), # Constrain beta0 to 1 - beta1 umxPath("gamma", to=latentvars, fixedAt=1, labels=paste("gamma",vars,sep="")), umxPath(var = "gamma", labels="beta0[1,1]", free=FALSE), # Use constraint to define beta0 # Compute additional results mxAlgebra(Angles * 180/pi, name="Degrees"), mxAlgebra(1-2*beta1, name="mc180"), mxAlgebra(1/(1+v), name="communalities"), mxAlgebra(sqrt(1/(1+v)), name="communalityIndices") ) mxGetExpected(test.model, "covariance") umxExpCov(test.model)
I can't get your example script to run, but I'm pretty sure I know what is causing the problem. I'm getting
from your
umxRAM
statement.The issue is with the algebra
mxAlgebra(sqrt(1/(1+v)), name="z")
. There's a difference between how the frontend evaluates algebras and how the backend evaluates them. Officially 'v' is a k x 1 matrix, but '1' is 1 x 1 matrix. In the frontend (i.e. when you trymxGetExpected
) we're saying you can't add a k x 1 matrix to a 1 x 1 matrix. This is correct but kind of rude. The backend is being friendlier. If you replace '1' in that and similar algebras by appropriately sized matrices (e.g.mxMatrix('Unit', nrow=k, ncol=1, name='oneVec')
), then it will work.To show that the issue isn't really with
mxGetExpected
try usingmxEval(z, test.model, compute=TRUE)
. You should get the same error as thrown bymxGetExpected
.Okay, that explains the behavior, and I can try to swap in unit matrices instead of scalar 1s. Thanks!
But it is really counterintuitive to me that scalar operations on matrices work fine in some parts of the program but not others. Based on the general behavior in R and that scalar operations on matrices work for model estimation, I would just expect these to work without issue. Is there any possibility that the behavior of the front-end functions could be updated to respect scalar operations?
I glad that helps, and I certainly agree that this shouldn't be a problem users ever have to face. We are looking into possibilities for eliminating this problem. I just wanted to give you a workaround in the meantime.
Thanks, really appreciate it!