Binomial_rng fails due to probability parameter being nan

Modeling
1

I want to perform a prior predictive check but have troubles with what I think is fairly easy model. I am confused why it does not work in stan, since I have successfully implemented the same prior predictive sampling procedure in R without any problems.

My stan code is the following:

data {
  vector[2] mu;
  vector<lower = 0>[2] prev_beta;
  cov_matrix[2] Sigma;
  int<lower=0> n;
  }
  
  model{
    // empty model block
  }

generated quantities {
  vector[2] logit_Se_Sp;
  real<lower=0, upper=1> prev;
  int<lower=0,upper=n> y;
  real<lower=0,upper=1> Se = inv_logit(logit_Se_Sp[1]);
  real<lower=0,upper=1> Sp = inv_logit(logit_Se_Sp[2]);
  real<lower=0,upper=1> p = prev*Se + (1-prev)*(1-Sp);
  logit_Se_Sp = multi_normal_rng(mu, Sigma);
  prev = beta_rng(prev_beta[1],prev_beta[2]);
  y = binomial_rng(n, p);
}

I use rstan to pass the priors:

sampling(prior_predictive_model, data=list(n=1000, mu=c(3,3), Sigma=matrix(c(1,0,0,1), nrow=2), prev_beta=c(1,1)), chains=4, iter=2000, refresh=500, cores=4, algorithm="Fixed_param")

The error message that I receive:

Exception: binomial_rng: Probability parameter is nan, but must be in the interval [0, 1]

I cannot make sense of this error message because as I read the code, all the parameters are correctly constrained. Any advice is highly appreciated.

2

This is because you don’t have any parameters or any priors. If you want to do a prior predictive check you need to declare the parameters of interest in the parameters block and then assign the desired prior(s) in the model block. Then, you can use these parameters in the generated quantities block for your prior predictive check

3

Thanks a lot! Indeed, the problem was that I misunderstood how the prior predictive check needs to be structured. I must admit that I find the section on prior predictive checks in the User Guide rather confusing.

For completeness, the following model does what I want:

data {
  vector[2] mu;
  cov_matrix[2] Sigma;
  vector<lower = 0>[2] prev_beta;
  int<lower=0> n;
}

parameters{
  vector[2] logit_Se_Sp;
  real<lower=0, upper=1> prev;
}

transformed parameters{
  real<lower=0,upper=1> Se = inv_logit(logit_Se_Sp[1]);
  real<lower=0,upper=1> Sp = inv_logit(logit_Se_Sp[2]);
  real<lower = 0, upper = 1> p = prev * Se + (1 - prev) * (1 - Sp);
}

model{
  logit_Se_Sp ~ multi_normal(mu, Sigma);
  prev ~ beta(prev_beta[1],prev_beta[2]);
}

generated quantities {
  int y = binomial_rng(n, p);
}