Problem with definition variables with state space expectation?

Am I doing something wrong? I was surprised when my model seemed to be ignoring my definition variables, but then I also notice that mxEval interprets them all as NA... this is using a build from yesterday, 1st Feb 2015, on 64 bit windows.
data(demoOneFactor)
nvar <- ncol(demoOneFactor)
varnames <- colnames(demoOneFactor)
demoOneFactorInputs <- cbind(demoOneFactor, V1=rnorm(nrow(demoOneFactor)))
ssModel <- mxModel(model="State Space Inputs Manual Example",
mxMatrix("Full", 1, 1, TRUE, .3, name="A"),
mxMatrix("Full", 1, 1, TRUE, values=1, name="B"),
mxMatrix("Full", nvar, 1, TRUE, .6, name="C", dimnames=list(varnames, "F1")),
mxMatrix("Zero", nvar, 1, name="D"),
mxMatrix("Diag", 1, 1, FALSE, 1, name="Q"),
mxMatrix("Diag", nvar, nvar, TRUE, .2, name="R"),
mxMatrix("Zero", 1, 1, name="x0"),
mxMatrix("Diag", 1, 1, FALSE, 1, name="P0"),
mxMatrix("Full", 1, 1, FALSE, labels="data.V1", name="u"),
mxData(observed=demoOneFactorInputs, type="raw"),
mxExpectationStateSpace("A", "B", "C", "D", "Q", "R", "x0", "P0", u="u"),
mxFitFunctionML()
)
ssRun <- mxRun(ssModel)
mxEval(data.V1, ssRun, compute=T, defvar.row=1)
Weird
What you're writing should work, but I'm replicating the same behavior. When you use the pre-run model, these evaluations work.
mxEval(data.V1, ssModel, compute=T, defvar.row=1)
mxEval(u, ssModel, compute=T, defvar.row=1)
And when you generate data based on a model with non-zero weights in the "B" and/or "D" matrices, the expectation correctly estimates the parameters. A model like this is in our test suite. http://openmx.psyc.virginia.edu/svn/trunk/models/passing/StateSpaceInputs.R
So, I'm very confident that internally the definition variables/inputs are being used correctly. There does seems to be a problem in evaluating them with mxEval after a model has been run.
I've created a Issue for this bug and I'll be fixing it soon. http://openmx.psyc.virginia.edu/issue/2015/02/mxeval-run-state-space-model
Thanks for finding it!
Log in or register to post comments
In reply to Weird by mhunter
Fixed
This ended up being a data sorting issue. It's now fixed in the latest svn (r4208) and will be part of the next binary release. Thanks again for finding this!
Log in or register to post comments
In reply to Fixed by mhunter
Ok great :) But since this
Ok great :) But since this was not the problem I was having, I worked out what was. According to my understanding, although you say that this:
x[t+1] = A x[t] + B u[t] + q[t]
is the state equation, it 'seems' to be functioning as this:
x[t+1] = A x[t] + B u[t+1] + q[t]
The following code hopefully shows either a) what I mean, or b) where my thinking is wrong!
A<-diag(.5,1)
B<-matrix(c(0,2.7),nrow=1)
Q<-diag(.1,1)
Tpoints<-200
u<-matrix(c(rep(0,Tpoints),sample(c(0,0,0,0,0,0,10),size=Tpoints,replace=T)),byrow=T,nrow=2)
q<-rnorm(200,0,sqrt(2))
x<-rep(NA,200)
x[1]<-1.5
for(i in 2:Tpoints){
x[i]<-A*x[i-1] + B%*%u[,i-1,drop=F] + q[i-1]
}
plot(x,type='l')
data<-cbind(x,u[2,])
colnames(data)<-c('y','input')
model<-mxModel(
mxData(data,type='raw'),
mxMatrix(name='A',free=T,nrow=1,ncol=1),
mxMatrix(name='B',free=T,nrow=1,ncol=2),
mxMatrix(name='C',values=1,free=F,nrow=1,ncol=1,dimnames=list('y','x')),
mxMatrix(name='D',values=0,free=F,nrow=1,ncol=2),
mxMatrix(name='Q',free=T,nrow=1,ncol=1),
mxMatrix(name='x0',free=T,nrow=1,ncol=1),
mxMatrix(name='P0',values=1,free=F,nrow=1,ncol=1),
mxMatrix(name='Q',free=T,nrow=1,ncol=1),
mxMatrix(name='R',values=0,free=F,nrow=1,ncol=1),
mxMatrix(name='u',values=c(1,0),labels=c(NA,'data.input'),free=F,nrow=2,ncol=1),
mxExpectationStateSpace('A', 'B', 'C', 'D', 'Q', 'R', 'x0', 'P0', 'u', dimnames='y'),
mxFitFunctionML()
)
fit<-mxRun(model)
summary(fit)
#now shift the u matrix in time
data[2:(Tpoints),'input']<-data[1:(Tpoints-1),'input']
model<-mxModel(
mxData(data,type='raw'),
mxMatrix(name='A',free=T,nrow=1,ncol=1),
mxMatrix(name='B',free=T,nrow=1,ncol=2),
mxMatrix(name='C',values=1,free=F,nrow=1,ncol=1,dimnames=list('y','x')),
mxMatrix(name='D',values=0,free=F,nrow=1,ncol=2),
mxMatrix(name='Q',free=T,nrow=1,ncol=1),
mxMatrix(name='x0',free=T,nrow=1,ncol=1),
mxMatrix(name='P0',values=1,free=F,nrow=1,ncol=1),
mxMatrix(name='Q',free=T,nrow=1,ncol=1),
mxMatrix(name='R',values=0,free=F,nrow=1,ncol=1),
mxMatrix(name='u',values=c(1,0),labels=c(NA,'data.input'),free=F,nrow=2,ncol=1),
mxExpectationStateSpace('A', 'B', 'C', 'D', 'Q', 'R', 'x0', 'P0', 'u', dimnames='y'),
mxFitFunctionML()
)
fit<-mxRun(model)
summary(fit)
Log in or register to post comments
In reply to Ok great :) But since this by CharlesD
Documentation fix
You're right. The state equation that's used has the "instantaneous" input effects for the dynamics (and the measurement).
x[t] = A x[t-1] + B u[t] + q[t]
not the lagged input effects.
x[t+1] = A x[t] + B u[t] + q[t]
If you want the lagged input effects you'll have to create a lagged variable yourself. I'll update the documentation to reflect the true state of things.
Log in or register to post comments
In reply to Documentation fix by mhunter
Ok, good to know then, all
Ok, good to know then, all working fine now.
Log in or register to post comments