Vectorizing Mixture Poisson model

Modeling
1

which one is a right way of vectorizing the mixture of Poisson model (truncated at (y > 1)).

    target += log_mix(pi[i],
                poisson_lpmf(y[i]|mu[1,i])-log1m_exp(poisson_lccdf(1|mu[1,i])),
                poisson_lpmf(y[i]|mu[2,i])-log1m_exp(poisson_lccdf(1|mu[2,i])));

equivalently

   target += log_mix(pi[i],
              poisson_lpmf(y[i]|mu[1,i])-log1m_exp(-mu[1,i] + log(1+mu[1,i])),
              poisson_lpmf(y[i]|mu[2,i])-log1m_exp(-mu[2,i] + log(1+mu[2,i])));

OR

target += log_mix(pi[i],
                  poisson_lpmf(y[i]|mu[1,i])-log1m_exp(poisson_lpmf(0|mu[1,i]))
                                            -log1m_exp(poisson_lpmf(1|mu[1,i])),
                  poisson_lpmf(y[i]|mu[2,i])-log1m_exp(poisson_lpmf(0|mu[2,i]))
                                            -log1m_exp(poisson_lpmf(1|mu[2,i])));

and after we vectorize it in the right way, is it possible to define the generated quantities in this way

generated quantities{
   real y_hat_1[N];
   real y_hat_2[N];
   for (i in 1:N){
     y_hat_1[i] = poisson_rng(mu[1,i]);   // should T[2,] included here or not?
     y_hat_1[i] = poisson_rng(mu[2,i]);
   }
}
2

If you’re trying to write:

p(y_i | \mu_1, \mu_2) = \pi p_\text{truncated}(y_i | \mu_1) + (1 - \pi) p_\text{truncated}(y | \mu_2)

And:

p_\text{truncated}(y_i | \mu) = \frac{p(y_i | \mu)}{P(y_i > 1 | \mu)}

I think it should be:

for(i in 1:size(y)) {
  target += log_mix(pi[i],
                  poisson_lpmf(y[i]|mu[1,i])-poisson_lccdf(1|mu[1,i]),
                  poisson_lpmf(y[i]|mu[2,i])-poisson_lccdf(1|mu[2,i]));
}

I checked log_mix in the docs and it doesn’t list a vectorized version yet. Do you intend to have a different mixing proportion for every y[i]? Or should that just be a single mu in the log_mix?

4

Yes, I intend to have a different mixing proportion for every y[i]. If we state the truncation like this, what will be the generated quantities? Is it make sense…?

generated quantities{
   real y_hat_1[N];
   real y_hat_2[N];
   for (i in 1:N){
     y_hat_1[i] = poisson_rng(mu[1,i]);   // should T[2,] included here or not?
     y_hat_1[i] = poisson_rng(mu[2,i]);
   }
}

OR

generated quantities{
   real y_hat_1[N];
   for (i in 1:N){
          y_hat_1[i] = log_mix(pi[i],
                                 poisson_rng(mu[1,i]),   
                                 poisson_rng(mu[2,i]));
}
}
5

It depends on what you’re looking for. There are a bunch of different things you could generate and look at.

As I understand it, your generative process is to first randomly pick one of two distributions and then sample from the distribution you picked. This looks like:

int y_hat[N];
for(i in 1:N) {
  int which = bernoulli_rng(pi[i]);
  y[i] = poisson_rng(mu[which + 1, i]);
}
1 Like