The parameters of model can not be recovered

Modeling
1

I simulated a set of data

library(tidyverse)

df <- tibble(
  x    = seq(from = -13,  to = 13, by = 1),
  z    = rnorm(n = 27, mean =  4.8 + 0.35 * x, sd = 2),
  y    = rpois(n = 27, lambda = exp(0.18 - 0.13 * z))
)
df %>% head(5)

# A tibble: 5 × 3
#      x      z     y
#  <dbl>  <dbl> <int>
#1   -13  1.57      1
#2   -12 -0.705     1
#3   -11  4.47      0
#4   -10  0.123     0
#5    -9 -0.471     3

I want to recover the parameters \beta s,

\begin{aligned} y_{i} &\sim \text{Poisson}\left( \lambda_i\right) \\ \lambda_i &=\exp \left(\beta_1 + \beta_2 z_i\right) \\ z_i & \sim \text{Normal}\left( \beta_3 + \beta_4 x_i, \sigma\right) \\ \end{aligned} 
data {
  int<lower=1> N;                                      
  vector[N] x;  
  array[N] int y;                          
}

parameters {
  vector[4] beta;
  real<lower=0> sigma;
}
model {
  vector[N] z;

  beta ~ normal(0, 20);                                                            
  sigma ~ exponential(1);

  z ~ normal(beta[3] + beta[4] * x, sigma);
  y ~ poisson_log(beta[1] + beta[2] * z);
}

but the model failed, pop up a message without samping

Running MCMC with 4 sequential chains...

I don’t know how to solve it. Thank you very much for your help.

code.Rmd (1.4 KB)

2

You need to move z to the parameters block.

To see why, recall that z ~ normal(beta[3] + beta[4] * x, sigma) is just a shorthand for target += normal_lupdf(z | beta[3] + beta[4] * x, sigma). But z is just an empty vector at this point; it doesn’t contain anything. You’re expecting z ~ normal(beta[3] + beta[4] * x, sigma) to define z to be a sample from the distribution on the right-hand side of the ~, but that’s not how Stan works.

The lack of an informative error message here is definitely not ideal, though. @WardBrian @stevebronder @rok_cesnovar , on my system if I run the rmd provided in the OP I can compile the model and run the line

fit <- model$sample(data = stan_data, refresh = 1)

with no warning messages or errors. Clearly it doesn’t sample; the only output to the console is

Running MCMC with 4 sequential chains...

But I only see an error if I then run, e.g.

> fit
Error: Fitting failed. Unable to print.

I have range checks turned off; but I don’t think that should matter?

> cmdstan_make_local()
[1] "STAN_THREADS=false"                         "STAN_NO_RANGE_CHECKS=true"                  "CXXFLAGS_OPTIM=-march=native -mtune=native"

@Wang_MinJie , as an aside, even once you fix this I don’t think you can expect this model to do a good job of parameter recovery. The parameters are not identified except via the prior (which is very weak). In particular, I can multiply beta[3], beta[4], and sigma by a constant and divide beta[2] by the same constant and the model will be completely unchanged.

5 Likes
3

Thank you for your help. The problem that has been troubling for a long time has been solved.
There is a problem with the model. I need to modify it