stan-dev / example-models

Example models for Stan
http://mc-stan.org/
771 stars 479 forks source link

Update models in Ch.06 of BPA #107

Closed ito4303 closed 7 years ago

ito4303 commented 7 years ago

Non-centered parameterization is applied for hierarchical models.

bob-carpenter commented 7 years ago

I didn't understand the non-centering in BPA/Ch.06/Mh.stan. You have (equivalently---I translated to use compound declare/define):

parameters {
  real<lower=0,upper=1> omega;      // Inclusion probability
  real<lower=0,upper=1> mean_p;     // Mean detection probability
  real<lower=0,upper=5> sigma;
  vector[M] sigma_raw;
  ...
transformed parameters {
  vector[M] eps = logit(mean_p) + sigma * sigma_raw;
  ...
model {
  ...
  sigma_raw ~ normal(0, 1);

I think it's partly confusing naming and partly an unusual approach to centering. I think what you really want to do if you want to parameterize in terms of log odds is to use

parameters {
  real mu_eps;
  vector[M] eps_raw;
  real<lower=0> sigma_eps;
  ...
transformed parameters {
  vector[M] eps = mu_eps + eps_raw * sigma_eps;
  ...
model {
  eps_raw ~ normal(0, 1);

You could alternatively draw eps_raw from a logistic(0, 1), which gives you the same distribution as a uniform(0, 1) variate that is log-odds transformed.

If you're going to change the models from the way they were coded in the book, I'd suggest not using those interval constrained parameters as they will mess up both computation and estimation if there would otherwise be any probability mass near the boundaries.

bob-carpenter commented 7 years ago

I'd also suggest rewriting the likelihood (which doesn't include sigma_raw) as:

for (i in 1:M) {
  if (y[i] > 0)
    target += bernoulli_lpmf(1 | omega)
                + binomial_logit_lpmf(y[i] | T, eps[i]);
  else 
    target += log_sum_exp(bernoulli_lpmf(1 | omega)
                            + binomial_logit_lpmf(0 | T, eps[i]),
                          bernoulli_lpmf(0 | omega));
}
ito4303 commented 7 years ago

Thank you for your comments. I fixed the mistakes and updated models using vectorization and compound declaration.