RW_llFunction_block

[Hamze Dokoohaki]

Hi Nimble team

I do have a process-based model with 5 parameters that is computationally expensive to run. So over a gird of these 5 parameters,  I developed a GP model to estimate the log-likelihood given a set of observations.  Now given a GP model that simulates the log-likelihood I'd like to optimize these 5 parameters. 

My nimble model looks like this :

```

code <- nimbleCode({

  p1 ~ T(dnorm(0.047, 0.5), 0.035, 0.055) 

  p2 ~ T(dnorm(-1.17, 1), -1.25, -1.05) 

  p3 ~ dunif(200, 500)

  p4 ~ T(dnorm(23, 15), 15, 30)

  p5 ~ T(dnorm(50, 30), 30, 80)

  sigma ~ dunif(0, 5)

  sigma2 ~ dunif(0, 5)

  

  ssq <- predic_NIMBLE_func(1, p1, p2, p3, p4 , p5)

  ssq2 <- predic_NIMBLE_func(2, p1, p2, p3, p4 , p5)

})

```

In the above code `predic_NIMBLE_func` is `nimbleRcall` which returns the log-likelihood. I have tested the  `predic_NIMBLE_func` and it works properly. 

Here is also my `nimbleFunction` for the  llFunction in `RW_llFunction_block` sampler. 

```

llFun <- nimbleFunction(

  setup = function(model) { },

  run = function() {

    

    p1 <- model$p1

    p2 <- model$p2

    p3 <- model$p3

    p4 <- model$p4

    p5 <- model$p5

    

    ssq <- model$ssq

    ssq2 <- model$ssq2

    

    sigma <- model$sigma

    sigma2 <- model$sigma2

 

    ll <- # estimate the likelihood 

    

    returnType(double())

    return(ll)

  }

)

RllFun <- llFun(Rmodel)

RllFun$run()

mcmcConf <- configureMCMC(Rmodel, nodes = NULL)

mcmcConf$addSampler(target = c("p1", "p2", "p3", "p4", "p5",

                               "sigma","sigma2","cc"),

                    type = "RW_llFunction_block",

                    control = list(llFunction = RllFun,

                                   adaptInterval=20, 

                                   adaptive=TRUE,

                                   includesTarget = FALSE))

```

So my problem is that, it looks like the MCMC doesn't use the defined `llFun` function mainly because the posteriors are exactly the same as priors.  I get the same posterior even when I set the `ll` to a constant in `llFun` function. 

It appears that I have an issue in the way I have defined my model. I was wondering if you could let me know if you see a major problem in my model. 

[Chris Paciorek]

hi Hamze,

I'm not seeing anything obviously wrong in what you are doing. Here's a suggestion for how to debug this:

Run the uncompiled version of the MCMC and step through the calculations of the sampler. You can do:

mcmc <- buildMCMC(mcmcConf)    # uncompiled MCMC

debug(mcmc$samplerFunctions[[1]]$run)

mcmc$run(5)      

And use R's debugger to step through what the sampler is doing. 

Let us know either way what you find out or if you have further questions.

-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/81c0f9a0-4407-4d0d-b125-1957fa2bcc48n%40googlegroups.com.

[Hamze Dokoohaki]

Hi Chris, 

Thanks for double checking the code. Sure, I'll give it a try and will share what I find for future reference. 

[Daniel Turek]

One additional comment: perhaps if your "llFun" nimbleFunction which you wrote, actually returns the *likelihood", as the comment in your code suggests:

  ll <- # estimate the likelihood 

In contrast to the *log-likelihood*, as is required.  If it mistakenly returns the (non-log) likelihood value, then realistically the likelihoods are always so small, e.g. 10^-10 or whatever, that regardless of the parameter values being used to calculate the likelihood, a likelihood of 10^-10 versus 10^-6 is essentially indistinguishable when it's used for (additive) log-likelihood calculations (being added onto the log-likelhood values of the parameter priors).  This theory is consistent with both (1) the comment from your code, and (2) the fact the sampling only follows the parameter priors.

Please let us know if this solves the problem, or if more detail would be helpful.  Or, if you're able to share a fully reproducible example, then one of us could dig further into what's happening.  Anyway, great that you're writing your own nimbleFuctions, and using the RW_llFunction_block sampler.

Daniel

To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/824dd687-7277-4d90-a90c-65893615c9c5n%40googlegroups.com.

[Hamze Dokoohaki]

Thank you Daniel; Really good point. Here is how I'm trying to combine my two likelihoods. Assuming residuals are gaussian and a fixed number of observations (5):

```   

    l.lai  <- ((5) * log(sigma)) - (ssq/(sigma*2) )

    l.ndvi <- ((5) * log(sigma2)) - (ssq2* (sigma2/2) )

    ll <-   l.lai  + l.ndvi 

```
Given the initial values the ll is around -5.

[Daniel Turek]

Looks like you're using the log-likelihood, which is correct.

Another thought. I see you're adding the RW_llFunction_block sampler to a bunch of valid model parameters ("p1", "p2", "p3", "p4", "p5","sigma","sigma2"), and also to the dubious model parameter "cc", which I don't see anywhere else.  It's plausible that the RW_llFunction_block sampler is trying to calculate the prior log-likelihood for these nodes (including "cc"), and the log-likelihood for "cc" is NaN (or similar) since "cc" doesn't exist in the model, rendering the sampling ineffective.  I haven't explored what would happen in this case, but just noticed that in your code.  Maybe try removing "cc" in the call to addSampler() and see how it goes?

Daniel

To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/29aa01b9-15ca-4d08-a8e1-80e8ab0701dcn%40googlegroups.com.