Using MCMC to evaluate the expected value of a high-dimensional function

I am working on a research project that is related to optimization and recently had an idea to use MCMC in this setting. Unfortunately, I am fairly new to MCMC methods so I had several questions. I'll start by describing the problem and then asking my questions.

Our problem boils down to estimating the expected value of a cost function $c(\omega)$ where $\omega = (\omega_1,\omega_2,...\omega_h)$ is an $h$-dimentional random variable with a density $f(\omega)$.

In our case, a closed form version of $c(\omega)$ does not exist. This means that we have to use Monte Carlo methods to approximate the expected value. Unfortunately, it turns out that estimates of $E[c(\omega)]$ that are generated using MC or QMC methods have too much variance to be useful within in a practical setting.

One idea that we had to use an importance sampling distribution to generate sample points that will produce a low variance estimate of $E[c(\omega)]$. In our case, the ideal importance sampling distribution, $g(\omega)$, has to be roughly proportional to $c(\omega)f(\omega)$. Seeing how $g(\omega)$ is known up to constant, I am wondering whether I can use MCMC along with the proposal distribution $c(\omega)f(\omega)$ to eventually generate samples from $g(\omega)$ .

My questions here are:

• Can MCMC be used within this setting? If so, what MCMC method would be appropriate? I am working in MATLAB, so I have a preference to anything that already has a MATLAB implementation.

• Are there any techniques that I can use to speed up the burn-in period for MCMC. And how can I tell that the stationary distribution has been reached? In this case, it actually takes a fair bit of time to calculate $c(\omega)$ for a given $\omega$.

• typo in the last paragraph? $C(w)$ instead of $c(w)$? – mpiktas Mar 17 '11 at 7:30
• You don't have to use Monte Carlo methods for high-dimensional integrals just because a closed form doesn't exist. You could use methods based on sparse grids. I don't know the pros and cons myself though. – onestop Mar 17 '11 at 8:26
• @onestop - you are right sparse grids work a little better in my case, though I'd also like to test this idea out. – Berk U. Mar 17 '11 at 8:53
• It sounds a bit like finance. – Wok Mar 17 '11 at 9:13
• It could be applied to finance. We are using it for the optimal management of power systems though :) – Berk U. Mar 17 '11 at 10:41

I would always remember, that MCMC is just a numerical integration tool (and a rather inefficient one at that). It is not some magic/mystical thing. It is very useful because it is reasonably easy to apply. It does not require much thinking compared to some other numerical integration techniques. For instance, you do not have to do any derivatives. You only have to generate "random numbers".

However, like any numerical integration method, it is not a universal catch all tool. There are conditions when it is useful, and conditions when it isn't.

It may be wiser to set up another technique. Depending on how big $h$ is, and how fast your computer is, and how much time you are prepared to wait for results. A uniform grid may do the job (although this requires small $h$ or a long amount of waiting). The "job" is to evaluate the integral - the equation does not care what meaning you or I attach to the result (and hence it does not care whether we obtained the result randomly or not).

Additionally, if your estimates of $\omega$ are quite accurate, the $f(\omega)$ will be sharply peaked and closely resemble a delta function, so the integral is effectively substituting $\omega\rightarrow\omega_{max}$.

Another numerical integration technique is using a taylor series under the integral. $f(\omega)\approx f(\omega_{max})+(\omega-\omega_{max})f'(\omega_{max})+\frac{1}{2}(\omega-\omega_{max})^{2}f''(\omega_{max})+\dots$

This is a useful strategy when the moments of $\omega$ are easily obtained.

Edwin Jaynes has a nice quote on this:

whenever there is a randomised way of doing something, there is a non-randomised way which yields better results, but requires more thinking

One "more thinking" way is to use "stratified MCMC" to do the integral. So rather than "randomly" pick a spot on the whole parameter space: divide it up into "strata". These "strata" should be picked so that you get a good range of the high part of the integral. Then randomly sample within each strata. But this will require you to write your own code I would imagine (i.e. more thinking).

• Thank you so much for your post. I very much like your stratified MC idea. I was wondering if you knew about any resources on it so that I could read up on it? Anything with an example in action would be much appreciated. – Berk U. Mar 18 '11 at 9:31
• This is an idea that I have had, so I don't know of any resources. It could be put into the MCMC framework, by adding another layer which defines the strata. And one of my main weaknesses is doing stuff - I'm a thinker at heart! I know this doesn't help you very much and I apologise for that. But surely I'm not the only person who has thought of this. Anyone who has done "design based" statistics and "Bayesian statistics" would surely have thought of this. – probabilityislogic Mar 18 '11 at 11:00

There isn't any indication that your variables here are correlated so I dont know why you would use MCMC as opposed to regular Monte Carlo. There are many different sampling methods including the mentioned stratified sampling (Latin hypercube) and QMC. Sparse quadrature methods are very good if the dimension of the problem is not too high (not more than 10) since sparse quadrature grids grow geometrically (curse of dimensionality).

But it sounds like you are on the right track with respect to importance sampling. The key here is to choose a biased distribution that has large probability concentrated near your region of interest and that it has thicker tails than the nominal distribution.

I'd like to add that this is an open research problem so if you can come up with something good it would be of great interest to the community!

• Thanks for this! Could you elaborate a little more about the relationship between the effectiveness of MCMC and the correlation between the variables? I had always thought that MCMC would provide an advantage over exisiting numerical integration techniques when the variables are correlated - though I cannot provide a good reason for why it would. – Berk U. Mar 22 '11 at 3:40

Since no one seemed to actually answer the question directly: yes you can use MCMC to sample from $g(\omega)$. MCMC can be used to sample from any distribution where the distribution is known only up to a constant of proportionality.

In addition, you may want to look up variance reduction techniques in the MC integration field. A great self contained set of resources are the free book chapters available from Art Owen at Stanford. Specifically chapters 8, 9, and 10.

There you will find in-depth treatments of adaptive sampling, recursion, and other techniques.