Hello all
I am a big fan of marginalising out latent variables as in the change point example in the manual. But this time, I have approx 100 discrete latent variables to marginalise out, which is generating computational issues. I am looking for either a computational trick that I missed in the manual or an alternative model approach. Any help would be great!
I am considering a simple augmented count model with augmented likelihood
n+m \sim Multinomial(\pi).
The observed counts n are from 100 locations, and we currently model missing counts m via a Negative Binomial
m \sim NB(n, p)
(in base R parameterisation) where p is a vector of observation probabilities. The expected value of the augmented counts, E(n+m), is n/p.
To marginalise out all 100 latent variables, I could put a cap m_{max} on m. But the dimension of the joint latent space is huge, \prod_i m_{i,max} \approx 1e135. I think this precludes marginalisation as in the change point example in the manual (15.2 p214).
Is there another trick that I could exploit for the augmented count model above?
Or maybe there is an approach for modelling missing counts under strong information on the observation probabilities p, that is more tractable in Stan?
Thank you very much for your input,
olli
Just for completeness, here is the full Stan model with no attempt of marginalising the latent variables (i.e. this won t compile, but may clarify the problem).
data {
// counts
int<lower=1> N; // number of observations, approx 100
int COUNTS[N]; // observed counts
// sampling
real<lower=0> P_ALPHA[N]; // sampling probabilities prior hyperparameter alpha
real<lower=0> P_BETA[N]; // sampling probabilities prior hyperparameter beta
}
transformed data {
// dirichlet hyperparameters depending on size of data
real<lower=0> prior_weight;
vector<lower=0>[N] dirprior;
prior_weight= N; // poor man casting
prior_weight= 0.8/prior_weight;
dirprior= rep_vector(prior_weight, N);
}
parameters {
simplex[N_PAIRS] props; // target parameter
vector<lower=0, upper=1>[N] prob_obs; // sampling parameter
int MISSING[N]; // 100 latent variables, to be integrated out
}
transformed parameters {
vector<lower=0, upper=1>[N] nb2_mu; // transformation to NegBin2 mu parameter
vector<lower=0, upper=1>[N] nb2_phi; // transformation to NegBin2 phi parameter
int AUX_COUNTS[N]; // augmented data
nb2_mu= COUNTS .* (1-prob_obs) ./ prob_obs;
nb2_phi= COUNTS;
AUX_COUNTS= COUNTS + MISSING;
}
model {
prob_obs ~ beta(P_ALPHA, P_BETA); // prior on sampling parameter
props ~ dirichlet_log(dirprior); // Dirichlet prior on target parameter
MISSING ~ neg_binomial_2_rng(nb2_mu, nb2_phi);
AUX_COUNTS ~ multinomial(props); // augmented likelihood
}