Semi-conjugate multinomial logistic regression with polya-gamma

[Ethan]

Hi all,

I have been working on implementing the polya-gamma logistic regression sampler for (possibly) multinomial observations. I am using the stick-breaking representation the multinomial distribution as discussed in this paper:

The polya-gamma part of the sampler seems to work fine, but I am getting the following (rather uninformative) error message about the regression coefficient sampler when I try to compile the MCMC:

> ## Build and compile MCMC 

> mcmc <- buildMCMC(mcmcConf) 

> cmcmc <- compileNimble(mcmc, project = model, showCompilerOutput = TRUE) Compiling 

    [Note] This may take a minute. 

    [Note] On some systems there may be some compiler warnings that can be safely ignored. 

Error in exists(as.character(varOrNodeExpr), x, inherits = FALSE) : 

      first argument has length > 1

The sampler appears to work with uncompiled using the uncompiled MCMC.

I have attached my code here, which contains four different files:

1. conjugate_logit_model.R - the main R file to compile the model and MCMC
2. polyaGammaRNG.R - the R file to sample from the polya-gamma distribution
   
3. pg_logit_conjugate_samplers.R - the R file that contains the samplers for the multinomial logistic regression model
   
4. pg_logit_helper_funs.R - contains helper functions (e.g., computation of probabilities based on stick-breaking probabilities)
   

The first file contains a reproducible example using a real data set. The only thing that needs to be changed is the directory at the beginning of the script.

Thanks in advance!

Best,

Ethan

[Chris Paciorek]

Hi Ethan,

You've run into an obscure limitation of how we handle model node access when you try to use "model[[omegaNames[i]]]" in one of your user-defined samplers. The model[[x]] syntax needs to have x be a single variable or node, and one can't use model[[x[i]]. This is documented in Section 15.4.2.1 of the manual, but it's easy to overlook. And the error message is completely unhelpful, which is our fault. I'm making a GitHub issue for us to improve the error message.

Instead, if you use:

values(model,omegaNames[i])[k]

it should work.

-Chris

P.S. Also, as a side note, you should take a look at this warning:

> source(file.path(dir, 'pg_logit_conjugate_samplers.R'))
Warning message:
In .checkFieldsInMethod(def, fieldNames, allMethods) :
  local assignment to field name will not change the field:
    Km1 <- gammaDims[2]
 Did you mean to use "<<-"? ( in method "run" for class "coefs_stick_sampler")

--
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/03f187f4-211b-4e52-9a34-b9299afb9defn%40googlegroups.com.

[Ethan]

Thanks a lot!