presence of mxCI breaks Standard Error reporting

neale's picture
Project:OpenMx
Component:Code
Category:bug report
Priority:normal
Assigned:Unassigned
Status:closed
Description

Adding mxCI() to a model makes all SE's reported as NaN. I don't think it should.

data(demoOneFactor)
 
Rff <- mxMatrix(type="Stand",nrow=1,ncol=1,free=F,name="Rff")
L <- mxMatrix(type="Full",nrow=5,ncol=1,free=T,values=0.2,labels=paste("l",1:5,sep=""),name="L")
I <- mxMatrix(type="Iden",nrow=5,ncol=5,name="I")
C <- mxAlgebra(L %*% Rff %*% t(L),name="C")
U <- mxAlgebra(I-(I*C),name="U")
SD <- mxMatrix(type="Full",nrow=5,ncol=1,free=T,values=0.6,
               labels=c("sd1","sd2","sd3","sd4","sd5"),lbound=0,name="SD")
Sigma <- mxAlgebra( (C+U) * (SD%*%t(SD)), name="Sigma", dimnames=list(colnames(demoOneFactor),colnames(demoOneFactor)))
 
factorModelStand <- mxModel(
  "1FactorStd", Rff, L, I, C, U, SD, Sigma, 
  mxMLObjective(covariance="Sigma"),
#  mxCI("L"),
  mxData(observed=cov(demoOneFactor), type="cov", numObs=500)
)
run_fms <- mxRun(factorModelStand,intervals=T)
summary(run_fms)
 
 
factorModelStand <- mxModel(
  "1FactorStd", Rff, L, I, C, U, SD, Sigma, 
  mxMLObjective(covariance="Sigma"),
  mxCI("L"),
  mxData(observed=cov(demoOneFactor), type="cov", numObs=500)
)
run_fms <- mxRun(factorModelStand,intervals=T)
summary(run_fms)

Comments

tbates's picture

#1

This does seem wrong, and is still true with current beta

packageVersion("OpenMx")
[1] ‘999.0.0.3156’
 
data(demoOneFactor)
m1 <- mxModel("1FactorStd",
	mxMatrix(type = "Stand", nrow = 1, ncol = 1, free = F, name = "Rff"),
	mxMatrix(type = "Full", nrow = 5, ncol = 1, free = T, values = 0.2, labels = paste0("l", 1:5), name = "L"),
	mxMatrix(type = "Iden", nrow = 5, ncol = 5, name = "I"),
	mxAlgebra(L %*% Rff %*% t(L), name = "C"),
	mxAlgebra(I-(I*C), name = "U"),
	mxMatrix(type = "Full", nrow = 5, ncol = 1, free = T, values = 0.6, labels = paste0("sd",1:5), lbound = 0, name = "SD"),
	mxAlgebra( (C + U) * (SD%*%t(SD)), name = "Sigma", dimnames = list(colnames(demoOneFactor), colnames(demoOneFactor))),
	mxMLObjective("Sigma"),
	mxData(cov(demoOneFactor), type="cov", numObs=500)
)
m1 <- mxRun(m1,intervals=T)
summary(m1)$parameters[,1:6]
#    name matrix row col  Estimate   Std.Error
# 1    l1      L   1   1 0.8913096 0.009723307
# 2    l2      L   2   1 0.9325548 0.006408531
# 3    l3      L   3   1 0.9438468 0.005514373
# 4    l4      L   4   1 0.9623626 0.004007500
# 5    l5      L   5   1 0.9725557 0.003275464
# 6   sd1     SD   1   1 0.4455831 0.014104647
# 7   sd2     SD   2   1 0.5400880 0.017096124
# 8   sd3     SD   3   1 0.6115843 0.019359291
# 9   sd4     SD   4   1 0.7302594 0.023115874
# 10  sd5     SD   5   1 0.8187199 0.025916025
 
m2 <- mxRun(mxModel(m1, mxCI("L"), name= "withCI"), intervals = T)
a = summary(m2)
round(a$CI,2)
#               lbound estimate ubound
# withCI.L[1,1]   0.87     0.89   0.91
# withCI.L[2,1]   0.92     0.93   0.94
# withCI.L[3,1]   0.93     0.94   0.95
# withCI.L[4,1]   0.95     0.96   0.97
# withCI.L[5,1]   0.97     0.97   0.98
 
a$parameters[,1:6]
#    name matrix row col  Estimate Std.Error
# 1    l1      L   1   1 0.8913096        NA
# 2    l2      L   2   1 0.9325548        NA
# 3    l3      L   3   1 0.9438468        NA
# 4    l4      L   4   1 0.9623626        NA
# 5    l5      L   5   1 0.9725557        NA
# 6   sd1     SD   1   1 0.4455831        NA
# 7   sd2     SD   2   1 0.5400880        NA
# 8   sd3     SD   3   1 0.6115843        NA
# 9   sd4     SD   4   1 0.7302594        NA
# 10  sd5     SD   5   1 0.8187199        NA

jpritikin's picture

#2

Verify that SVN 3209 fixes it and close the bug.

mhunter's picture

#3

This works.

require(OpenMx)
data(demoOneFactor)
m1 <- mxModel("1FactorStd",
 mxMatrix(type = "Stand", nrow = 1, ncol = 1, free = F, name = "Rff"),
 mxMatrix(type = "Full", nrow = 5, ncol = 1, free = T, values = 0.2, labels = paste("l", 1:5, sep=''), name = "L"),
 mxMatrix(type = "Iden", nrow = 5, ncol = 5, name = "I"),
 mxAlgebra(L %*% Rff %*% t(L), name = "C"),
 mxAlgebra(I-(I*C), name = "U"),
 mxMatrix(type = "Full", nrow = 5, ncol = 1, free = T, values = 0.6, labels = paste("sd",1:5, sep=''), lbound = 0, name = "SD"),
 mxAlgebra( (C + U) * (SD%*%t(SD)), name = "Sigma", dimnames = list(colnames(demoOneFactor), colnames(demoOneFactor))),
 mxMLObjective("Sigma"),
 mxData(cov(demoOneFactor), type="cov", numObs=500)
 )
#Warning message:
#In mxMLObjective("Sigma") :
#  Objective functions like mxMLObjective() have been deprecated in favor of expectation and fit functions.
# Please use mxExpectationNormal(covariance= , ...) instead, and add a call to mxFitFunctionML()
# See examples at help(mxExpectationNormal)
m1 <- mxRun(m1,intervals=T)
#Running 1FactorStd 
summary(m1)$parameters[,1:6]
#   name matrix row col  Estimate   Std.Error
#1    l1      L   1   1 0.8913093 0.009723289
#2    l2      L   2   1 0.9325546 0.006408525
#3    l3      L   3   1 0.9438466 0.005514378
#4    l4      L   4   1 0.9623625 0.004007501
#5    l5      L   5   1 0.9725556 0.003275467
#6   sd1     SD   1   1 0.4455827 0.014104545
#7   sd2     SD   2   1 0.5400874 0.017095994
#8   sd3     SD   3   1 0.6115837 0.019359145
#9   sd4     SD   4   1 0.7302586 0.023115682
#10  sd5     SD   5   1 0.8187190 0.025915812
m2 <- mxRun(mxModel(m1, mxCI("L"), name= "withCI"), intervals = T)
#Running withCI 
a = summary(m2)
round(a$CI,2)
#              lbound estimate ubound
#withCI.L[1,1]   0.87     0.89   0.91
#withCI.L[2,1]   0.92     0.93   0.94
#withCI.L[3,1]   0.93     0.94   0.95
#withCI.L[4,1]   0.95     0.96   0.97
#withCI.L[5,1]   0.97     0.97   0.98
a$parameters[,1:6]
#   name matrix row col  Estimate   Std.Error
#1    l1      L   1   1 0.8913093 0.009723313
#2    l2      L   2   1 0.9325546 0.006408544
#3    l3      L   3   1 0.9438466 0.005514383
#4    l4      L   4   1 0.9623625 0.004007503
#5    l5      L   5   1 0.9725556 0.003275467
#6   sd1     SD   1   1 0.4455827 0.014104607
#7   sd2     SD   2   1 0.5400874 0.017096070
#8   sd3     SD   3   1 0.6115837 0.019359224
#9   sd4     SD   4   1 0.7302586 0.023115792
#10  sd5     SD   5   1 0.8187190 0.025915938

AdminHunter's picture

#4

Status:active» closed

Thanks, Joshua!