Calculating log_lik for occupancy (mixture) models

[mbjoseph]

Looks like you have everything you need here to compute the log likelihood:

@mattwilliamson13

// Likelihood
  for (i in 1:T) {
    if (sum_y[i]) { // Occupied and observed
      target += bernoulli_logit_lpmf(1 |  logit_psi[i])
        + bernoulli_logit_lpmf(y[i, 1:last[i]] | logit_p[i, 1:last[i]]);
    } else {        // Never observed
                            // Occupied and not observed
      target += log_sum_exp(bernoulli_logit_lpmf(1 | logit_psi[i])
                            + bernoulli_logit_lpmf(0 | logit_p[i, 1:last[i]]),
                            // Not occupied
                            bernoulli_logit_lpmf(0 | logit_psi[i]));
    }
  }

If you want a vector log_lik with T elements, you can store the log likelihood values in the transformed parameters block, and increment the log density in the model block as follows:

...
transformed parameters{
  ...
  vector[T] log_lik;

  for (i in 1:T) {
    if (sum_y[i]) {
      // observed, present
      log_lik[i] = bernoulli_logit_lpmf(1 |  logit_psi[i])
                   + bernoulli_logit_lpmf(y[i, 1:last[i]] | logit_p[i, 1:last[i]]);
    } else {
      // not observed
      log_lik[i] = log_sum_exp( // either:
                     // present and not detected
                     bernoulli_logit_lpmf(1 | logit_psi[i])
                       + bernoulli_logit_lpmf(0 | logit_p[i, 1:last[i]]),
                     // or absent
                     bernoulli_logit_lpmf(0 | logit_psi[i]));
    }
  }
  ...
}
model {
  ...
  // increment the log density
  target += sum(log_lik);
  ...
}
...