Impute a distribution

Modeling
1

I am doing an applied project. I am wondering if there is any way to impute a distribution in Stan.

Separating a larger model into several smaller models can be either a computational approximation (e.g., replace a multilevel model with a two stage regression), or even a theoretical necessity (e.g., the calibration model and regression model can be assumed as independent in a measurement-error model, or the ''do" operator in causal inference.).

To start with, consider a model

data{
real y;
}
parameter{
real<lower=0>  sigma;
}
model{ 
y~normal (beta, sigma);
}

Instead of making inference on the exact posterior distribution of beta given y, suppose in the ideal situation I have already fitted another calibration model and know beta is exactly 1.
Then I can impute beta into the data block

data {
...
real beta;
}

More realistically, I may only know the distribution of beta from the separately-fitted calibration model, say I know p(beta | calibration model )= N(0,1). And I assume beta is independent of y, such that p(beta | y)=N(0,1). Conceptually, I am imagining something like

data {
...
real beta;
beta~normal(0,1);
}

In other words, it imputes the value of beta by a pre-fixed distribution. If beta is discrete, I know how to do that, as I can marginalize out beta by finite mixture. When beta is continuous, I have to calculate the likelihood as an integration. Setting a prior beta~N(0,1) in the model block does not achieve the same effect.

Arguably, this is not the most Bayesian approach, for I should always update beta after I see y. Nevertheless at least as an approximation inference, is there anyway to impute the distribution (in the posterior distribution)?

1 Like
2

Can’t you still marginalise it out by Monte Carlo (with a lot of log-sum-exps)? Maybe with the new fancy MPI stuff if the number of imputed draws is big?

3

Integrating the likelihood Monte-Carloly is one solution, but probably not scalable with the dimension of distribution to be imputed especially when it is done manually. The bias-variance trade-off again: a slow Monte-Carlo integration vs a biased point estimation.
And yes, MPI should give a more elegant approach.

More thoughts: I can always mix HMC with Gibbs sampling (via a loop outside Stan), but that seems to miss the main point. When I don’t impute distributions and make exact HMC sampling, the computation cost is d^0.25, which is in fact less than Gibbs. Shouldn’t I do even better if I know a priori the posterior distributions of some parameters?

4

I’m not really sure what you’re looking for here. You need to integrate out the parameter and if all you’ve got is a bag of samples, you’ll need to do Monte Carlo or some variant. If you’ve got a continuous function (like beta ~ normal(0,1)) there either is or will soon be an integrate_1d function that you could use.

5

Looks like you would like to have “cut” functionality (see. e.g.
Better together? Statistical learning in models made of modules). WinBUGS has this (maybe JAGS, too), but it’s not possible in Stan directly.

You can reduce the required number of Monte Carlo draws by various quadrature rules up to few dimensions, CCD up to 10-15 dimensions, and quasi Monte Carlo for more dimensions.

1 Like