helske / KFAS

KFAS: R Package for Exponential Family State Space Models
64 stars 17 forks source link

Overdispersed distributions #40

Closed sschloss1 closed 3 years ago

sschloss1 commented 5 years ago

I'm interested in running an overdispersed binomial model. I was looking at the example in the vignette, and it looks like the overdispersed Poisson was handled by adding an additional random term to the state. I'm not sure that this is really an overdispersed Poisson, as my understanding is that the overdispersion should be at the observation level, not in the state. When I tried a formulation similar to the Poisson model in the vignette for my binomial model, I also got similar results for the standard (not overdispersed) and "overdispersed" formulations. If anything, the model with the additional state-level error term was actually more "wiggly" than the overdispersed model, which is the opposite of what should be happening with an overdispersion parameter.

For non-Gaussian models, would it be possible to add a random coefficient at the observation level so that it does not affect the state? I think would be equivalent to adding an H matrix to the non-Gaussian models.

helske commented 5 years ago

Hmm, I'm not sure I follow. The trick in the vignette is just like "adding extra H matrix", it is just augmented to state equation (additional state of form alpha_t ~ N(0, sigma^2)) and it then goes through the link function automatically (which you have to do in order to satisfy the distributional constraints like positivity in Poisson). You can do similar thing in Gaussian case as well, where you can either define H in observation equation or move the error term to state equation (see ?transformSSM function).

But, the thing is that you have to be careful in analyzing the results in this "augmented error" case. Which I just realized I hadn't written out in the vignette (and nobody has commented the issue before!): If you just grab the mu from the output of KFS, it gives you exp(mu+eps) in the 6.5 example, not exp(mu) what you want. So what you have to do is pick up the correct states from the output, for example like exp(coef(out_poisson, states = "level")) in the example (which is actually slightly biased due to exp(E(mu|y)) not being equal to E(exp(mu|y)). You could use importance sampling for asymptotically exact result if needed).