This post is about some results from “Bias Properties of Budget Constrained Simulations“, by Glynn & Heidelberger and published in Operations Research in 1990. I have found these results extremely useful, and our latest manuscript on unbiased MCMC recalls them in detail. Below I go through some of the results and describe the simulations that lead to the above figure.
Consider the following setting: you have a time budget (say, one day), you have access to a number of machines, and you have Monte Carlo simulations to run. Each simulation produces a random variable with expectation , in a random time denoted by . The variable could, for instance, be the output of a rejection sampler, or of some Markov chain-based algorithm such as coupling from the past, some unbiased MCMC procedure, or some sequential Monte Carlo sampler for normalizing constant estimation. In any case you’re interested in .
You could first ask each machine to produce precisely estimators. But you would not be sure how long it would take for the calculations to finish. Also, it could be wasteful: suppose that machine finishes its calculations faster than machine : you would prefer machine to keep on producing estimators while waiting for to finish.
So instead, you prefer to ask each machine to produce as many estimators as possible within the given time budget, starting a new estimator whenever the previous one completes. But what exactly do you do when the budget expires? Do you interrupt the machines? Or do you wait for them to complete their on-going calculations? What else? These are the types of question that Glynn & Heidelberger help to answer.
Let’s introduce some notation (exactly that of Glynn & Heidelberger). Suppose that a machine produces in sequence, in times denoted by . Let us denote by the number of completed estimators at time . Define and for . Thus denotes the estimator obtained if we interrupt the machine at time and average over the available estimators; if no estimator is available yet, then .
Equation (1) in Glynn & Heidelberger states (assuming ):
Thus the estimator is biased for the object of interest , and the bias diminishes with . The proof of that result starts by noting the conditional exchangeability of the pairs given . That conditional exchangeability is used to compute, for any , . Finally, a multiplication by , a sum over and the application of Fubini’s theorem conclude the proof.
An alternative estimator could be , obtained if we wait for on-going calculations to complete. That estimator is biased too, and Theorem 2 in Glynn & Heidelberger gives explicitly:
The bias of that estimator goes to zero as , as indeed stated in Corollary 10 of the paper.
It turns out that one can obtain an unbiased estimator of with a mild modification of the above estimators. Define if and if . That is, if at least one estimator is already available, interrupt the on-going calculation at time ; otherwise wait for the first estimator to complete. Then . This is stated in Corollary 7 of Glynn & Heidelberger, and is a consequence of the first equation mentioned above.
Now, what about the figure above? It represents the bias and the expected completion time of the three estimators: , , and . The plots show bias (left) and expected completion time (right) as a function of the budget. In this example , with a mean , and , thus creating dependencies between estimators and compute times. The three estimators are represented by different colours; I’ll let you guess which one is which!