Monte Carlo Crash Course: Sampling

Summary

Monte Carlo Crash CourseSamplingIn the previous chapter, we assumed that we can uniformly randomly sample our domain. However, it’s not obvious how to actually do so—in fact, how can a deterministic computer even generate random numbers?Pseudo-Random NumbersFortunately, Monte Carlo methods don’t need truly random numbers. Instead, we can use a pseudo-random number generator (PRNG). A PRNG produces a deterministic stream of numbers that look uniformly random:By “look uniformly random,” we mean the sequence exhibits certain statistical properties:Uniformity: samples are evenly distributed.Independence: previous samples cannot be used to predict future samples.Aperiodicity: the sequence of samples does not repeat.Deterministic generators cannot fully achieve these properties, but can get pretty close, in a precise sense. Here, we will use the PCG family of generators, which are performant, small, and statistically robust.PRNGs give us uniformly random scalars, but we ultimately want to sample complex, high-dimensional domains. Fortunately, we can build up samplers for interesting distributions using a few simple algorithms.Uniform Rejection SamplingRejection sampling transforms a sampler for a simple domain $$D$$ into a sampler for a complex domain $$\Omega$$, where $$\Omega \subseteq D$$. All we need is a function $$\text{accept}$$ that indicates whether a point $$\mathbf{x} \in D$$ is also contained in $$\Omega$$.Let’s build a rejection sampler for the two-dimensional unit disk. First, we’ll choose $$D = [-1,1]\times[-1,1]$$, which clearly encloses $$\Omega$$. We may use a PRNG to produce a sequence of uniform samples of $$[-1,1]$$, denoted as $$\xi_i$$. Taking each pair $$D_i = (\xi_{2i},\xi_{2i+1})$$ then provides samples of $$D$$.Second, we’ll define $$\text{accept}(\mathbf{x})$$—for the unit disk, we may check that $$||\mathbf{x}|| \le 1$$. Now, the rejection sampler:def Ω(): x = D() if accept(x): return x return Ω() In other words, sample $$D$$, and if the resul...