Is enkf() working well? There is no example on the documentation about enkf()

[marius190]

Hi Kingaa,

I am using the enkf() to estimate some models with the aim to estimate more complex model for meteoroligical data assimilation problems. Do you think it works when having big models lets say 100 variables? And is any way to define this variables in one time, without being necesarry to define it manually givin name to each variable.

I have a problem that I can't figure out what happens. The enkf() gives me back just a simulation of the hiden state variable xt (when I look for example at pred.mean(kf)). I doesn't matter the observation vector I give, and the vriance of it, always the model prediction return back the simulation of xt...

I wonder if you can make some example about how to use enkf() because I found no one on the doc or internet.

[kingaa]

Why is it that you expect something other than what you observe? For your simple problem, what is the exact theoretical expectation?

enkf is, by all present indications, working well. If you find evidence otherwise, I would be happy to see it so that I can fix the problem.

The reason there are no examples is that no one has yet contributed any examples. I would be happy to include contributed examples in the documentation.

I am closing this issue. In future, @marius190, please do not open up multiple issues on very closely related topics. Also, please spend a bit of time making sure that your questions are well phrased, that the codes you furnish produce the results you claim they do, and that your issue text is well formatted. Future readers will thank you.

[marius190]

I would expect the prediction (pred.mean()) would be something like that attached image. enkf() is no taking into account the observation vector. You can give a "dumb" y vector, for example y=1:100 that the result would be the same. I suspect that the kalman gain is computed as 0: x^a=x^f+K(y-Hx^f)=x^f, f for forecast, a for analysis.

PD: I think is a different issue, nevertheless I will take in account your advice. I am not native english speaker, so to "well phrase" is quite a challenge ;).

[kingaa]

Yes, the Kalman gain in this case is 0. This is because your state process is deterministic. There is no variance in the forecast. Míralo:

### Non Harmonic oscillator example
## Stochastic latent state process,
## Normal measurement noise.

R <- g^2 # note g is s.d., R is g^2

OscilatorData$y <- 1:100 

# *Stochastic* X process
oscilator.proc.sim <- function (x1, x2, omega, lambda, delta.t, ...) {
  Mk <- matrix(
    c(2+omega*omega-lambda*lambda*x1*x1,1,-1,0),
    nrow=2,ncol=2
  )
  x <- rnorm(n=2,mean=Mk %*% c(x1,x2),sd=1)
  setNames(x,c("x1","x2"))
}

kf <- enkf(OscilatorData,times="Time",t0=1,
  Np=10, R=R, h=h, params=theta,
  rprocess=discrete_time(step.fun = oscilator.proc.sim, delta.t=1))

## Note the value of Np.

plot(kf)
logLik(kf)

## is the filter mean identical to the simulation?
all(filter.mean(kf)==t(x))

And thanks for working on posing your questions well: the extra effort is noted!