functions {
  /*
  * Alternative to poisson_log_rng() that  avoids potential numerical 
  problems during warmup
  */
  int poisson_log_safe_rng(real eta) {
    real pois_rate = exp(eta);
    if (pois_rate >= exp(20.79))
      return -9;
    return poisson_rng(pois_rate);
  }
}

data {
  // The input data are: 
  int<lower=1> N;
  int<lower=0> complaints[N]; 
  vector<lower=0>[N] traps;
}

parameters {
 // Our model  accepts two parameters 'alpha' and 'beta'
 real alpha;
 real beta1; 
}

model {
   // weakly informative priors:
   alpha ~ normal(log(4), 1);
   beta1 ~ normal(-0.25, 1);
   // likelihood: Consider poisson_log(eta), which is  more efficient and stable alternative to poisson(exp(eta))
   complaints ~ poisson_log(alpha + beta1 * traps); 
  // Alternatively
  //complaints ~ poisson(exp(alpha + beta1 * traps));
} 

// Equivalently
// model {
//   target += normal_lpdf(alpha|log(4),1);
//   target += normal_lpdf(beta1|-0.25,1);
//   target += poisson_log_lpmf(complaints|alpha + beta1 * traps);
// } 

//Or
// model {
// beta1 ~ normal(-0.25, 1);
// alpha ~ normal(log(4), 1);
// target += poisson_log_lpmf(complaints|alpha + beta1 * traps);
// }

generated quantities {
  // sample predicted values from the model for posterior predictive checks
  int<lower=0> y_rep[N];
   
  vector[N] log_lik; 
 
  // for each obs we simulate a new value of y using the simulated parameters
  // the result will be a matrix with nrow = iter and ncol = n
  for (n in 1:N) {
    real eta_n = alpha + beta1 * traps[n];
    y_rep[n] = poisson_log_safe_rng(eta_n);
    log_lik[n] = poisson_log_lpmf(complaints[n]| eta_n);  
  }
}