PMCMC Implementation in NIMBLE to Multi-Parameter, Multi-State SSMs

[Joe Harris]

Dear Nimble Developer,

I have successfully implemented PMCMC estimation in NIMBLE for models with multiple parameters but a single latent state.

I am now planning to extend this approach to state-space models with multiple latent states, which may also include multiple observation variables. However, I have encountered several issues. For example, when compiling the model and building the MCMC, I receive the following error:

buildMCMC Error in findLatentNodes(model, nodes, timeIndex) : findLatentNodes: all latent nodes must come from the same variable

Interestingly, if I run a standard full Bayesian MCMC with multiple parameters and multiple latent states, it works without problems. I suspect the issue arises when I use a sampler like:

type = "RW_PF", control = list(filter = pf, Np = Np, latents = "x")

Could you help me identify the cause of this problem? Ideally, it would be extremely helpful if you could provide an example of implementing PMCMC for a state-space model with multiple parameters and multiple latent states in NIMBLE.

Finally, do I need to define custom NIMBLE functions to implement PMCMC for such models, or can this be accomplished using the existing built-in functions and samplers?

Thank you very much for your time and guidance.

Kind regards

Joe

[Perry de Valpine]

Hi Joe,

Thanks for the question. I'm glad to hear of this application.

The requirement is that non-scalar latent states must be arranged along a dimension of a larger object, for example as rows of a matrix. It is written that way to simplify the management of sets of particles. By default, the longest dimension of the object will be assumed to index time, but this can be controlled by the control list element "timeIndex" (e.g control = list(timeIndex=2))

One could imagine a more general implementation to handle more arbitrary cases, but that's not the way it is at the moment.

You can check and explore how cases will be handled by calling `findLatentNodes` directly. From the code in buildBootstrapFilter (or one of the others), it is called early (after extracting control list elements and inserting defaults), as:

    nodes <- findLatentNodes(model, nodes, timeIndex)  

I see findLatentNodes is not exported, so you can check with

`nimbleSMC:::findLatentNodes(model, nodes, NULL)`,

where `NULL` would be the value of timeIndex if it wasn't provided, to get default handling of the time index.

Here is an example of a toy model designed just to understand how latent nodes will be set up. The model is meaningless and I gave it just one real parameter (rho) and then set up state vectors of length 2 (x[1:2, t]) and a single observation for each one (but there could be arbitrary observations).

m <- nimbleModel(

 nimbleCode({

  rho ~ dnorm(0, 1)

  x[1, 1] ~ dnorm(0, sd=1)

  x[2, 1] ~ dnorm(0, sd = 1)

  for(t in 2:T) {

    x[1, t] ~ dnorm(rho * x[1, t-1], sd=1)

    x[2, t] ~ dnorm(x[1, t], sd=1)

    y[t] ~ dnorm(x[2, t], sd=1)

   }

 }),

 constants = list(T = 10)

)

nimbleSMC:::findLatentNodes(m, "x")

#  [1] "x[1:2, 1]"  "x[1:2, 2]"  "x[1:2, 3]"  "x[1:2, 4]"  "x[1:2, 5]"  "x[1:2, 6]"  "x[1:2, 7]"  "x[1:2, 8]"

#  [9] "x[1:2, 9]"  "x[1:2, 10]"

I see that between the JSS paper (https://www.jstatsoft.org/article/view/v100i03), the nimble User Manual, and the help pages in nimbleSMC, only the latter tries to explain the case of non-scalar latent states and that is quite brief. I'll file an issue suggesting to expand that and/or provide a vector example in the future. We could also consider exporting `findLatentNodes` so it is more visible for debugging and understanding the system.

If you're still stuck after absorbing this, please post an example.

HTH

Perry

--
You received this message because you are subscribed to the Google Groups "nimble-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nimble-users...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/nimble-users/87eaec4e-da37-4bd6-ae5d-9bf08e20fb19n%40googlegroups.com.