Problem using Particle Filter in simple SI model

[Kok Ben Toh]

Hi,

I was hoping to use PMCMC on a more complicated compartmental SIR model. After a number of setbacks, I decided to start from the simple and use particle filter to evaluate a simple SI model (Susceptible -> Infected).

Here's my code:

library(nimble)

library(nimbleSMC)

#### Model

SI <- nimbleCode({

  S[1] <- iniState[1]

  I[1] <- iniState[2]

  y[1] ~ dpois(pdet * I[1] + 1e-5)

  for (t in 2:nT) {

    lambda[t] <- beta * I[t-1] / N * S[t-1] + 1e-5

    negdS[t] ~ dpois(lambda[t])

    

    S[t] <- S[t-1] - negdS[t]

    I[t] <- I[t-1] + negdS[t]

    

    y[t] ~ dpois(pdet * I[t] + 1e-5)

    X[1:2,t] <- c(S[t], I[t])

  }

  

  beta ~ dunif(0, 0.3) #Transimissibility

  pdet ~ dunif(0, 1) #Probability of detection

})

states <- c(10000, 5)

constants <- list(

  N = sum(states),

  nT = 50,

  iniState = states

)

y <- c(0, 3, 2, 2, 8, 9, 3, 6, 7, 5, 6, 10, 14, 9, 8, 16, 

       15, 20, 27, 29, 37, 33, 32, 52, 55, 67, 56, 68, 85, 85, 111, 

       124, 145, 169, 206, 203, 260, 277, 320, 335, 414, 438, 539, 621, 

       681, 783, 932, 1015, 1099, 1222)

SIMod <- nimbleModel(code = SI,

                     constants = constants,

                     data = list(y = y),

                     inits = list(beta = 0.2, pdet = 0.5))

my_AuxF <- buildAuxiliaryFilter(SIMod, 'I',

                                control = list(saveAll = TRUE))

SIModC <- compileNimble(SIMod)

my_AuxFC <- compileNimble(my_AuxF, project = SIMod)

logLik <- my_AuxFC$run(m = 10000)

logLik # Always -Inf

The code runs and my nimble model should work since I could simulate data using SIMod$simulate() and the simulated dataset looks correct. However, the log-likelihood evaluated by the algorithm is always -Inf. I am not sure if it was my model's problem, or my code's problem. I have tried using Bootstrap filter instead, and I have swapped out the Poisson distribution with Normal and I still get -Inf log-likelihood. Might anyone know what could have been the problem here?

Thank you!

Cheers,

Ben

[Chris Paciorek]

hi Ben,

It looks like the -Inf is coming from numerical underflow in our calculations of weights in the particle filtering algorithm. This occurs for this specific model (and not others on which we've used nimble's particle filters) because of the large counts in the Poisson likelihood (producing very small likelihood values) as well as because the structure of the model induces many dependent data nodes for each latent node (resulting in a sum of many large magnitude negative log-likelihoods). That said, the underlying issue is that our numerical approach for working with the weights and log-likelihood values could be much more robust. Thanks for posting your question so that we are able to recognize this issue.

We'll try to make it a priority to work on this in the next couple weeks and at least get a branch on our GitHub repository that addresses the issue (making a release would take a bit longer).

Side note: in your call to buildAuxiliaryFilter, the 'nodes' argument should be 'negdS', since that is the variable containing the latent **stochastic** nodes. 

-chris

--
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 on the web visit https://groups.google.com/d/msgid/nimble-users/52a853c5-cdde-4ecf-804c-4a61b43101a3n%40googlegroups.com.

[Chris Paciorek]

hi Ben,

I've reworked the numerical approach to avoid the underflow problem in a branch on our GitHub repository. That should address the problem you were facing, so please feel free to try again after installing the updated version from the branch:

remotes::install_github("https://github.com/nimble-dev/nimbleSMC", ref = "fix_underflow", subdir = "packages/nimbleSMC")

We'll get this into a release of nimbleSMC in coming weeks.

-chris

[Kok Ben Toh]

Seems to work now! Thank you!