Monte Carlo 速成教程:Sampling 方法详解
Max Slater
计算机图形学、编程和数学
Monte Carlo 速成教程
Sampling
在前一章节中,我们假设可以均匀随机地对我们的域进行采样。然而,如何实际做到这一点并不明显 —— 事实上,一台确定性的计算机如何生成随机数呢?1
- Pseudo-Random Numbers
- Uniform Rejection Sampling
- Non-Uniform Rejection Sampling
- Inversion Sampling
- Changes of Coordinates
Pseudo-Random Numbers
幸运的是,Monte Carlo 方法并不需要真正的随机数。2 相反,我们可以使用 pseudo-random number generator (PRNG)。PRNG 产生一个确定性的数字流,看起来是均匀随机的: E = 0State: -8855663050364721856♺ 通过“看起来是均匀随机的”,我们指的是该序列表现出某些统计特性:
- Uniformity: 样本均匀分布。
- Independence: 之前的样本不能用于预测未来的样本。3
- Aperiodicity: 样本序列不重复。
确定性生成器无法完全实现这些属性,但在精确意义上可以非常接近。4 在这里,我们将使用 PCG 系列生成器,它们性能良好、体积小且统计鲁棒。 PRNG 为我们提供均匀随机的标量,但我们最终想要对复杂的高维域进行采样。 幸运的是,我们可以使用一些简单的算法来构建有趣分布的采样器。
Uniform Rejection Sampling
Rejection sampling 将简单域 DDD 的采样器转换为复杂域 Ω\OmegaΩ 的采样器,其中 Ω⊆D\Omega \subseteq DΩ⊆D。我们所需要的只是一个函数 accept\text{accept}accept 来指示点 x∈D\mathbf{x} \in Dx∈D 是否也包含在 Ω\OmegaΩ 中。
让我们为二维单位圆盘构建一个 rejection sampler。首先,我们将选择 D=[−1,1]×[−1,1]D = [-1,1]\times[-1,1]D=[−1,1]×[−1,1],它显然包含 Ω\OmegaΩ。我们可以使用 PRNG 来生成 [−1,1][-1,1][−1,1] 的均匀样本序列,表示为 ξi\xi_iξi。然后取每一对 Di=(ξ2i,ξ2i+1)D_i = (\xi_{2i},\xi_{2i+1})Di=(ξ2i,ξ2i+1) 即可提供 DDD 的样本。
其次,我们将定义 accept(x)\text{accept}(\mathbf{x})accept(x) —— 对于单位圆盘,我们可以检查 ∣∣x∣∣≤1||\mathbf{x}|| \le 1∣∣x∣∣≤1。 现在, rejection sampler:
def Ω():
x = D()
if accept(x):
return x
return Ω()
换句话说,对 DDD 进行采样,如果结果不在 Ω\OmegaΩ 中,那就再试一次! D:000 samplesΩ:000 samples\begin{align*} D : \phantom{00}0\text{ samples} \\ \Omega : \phantom{00}0\text{ samples}\end{align*}D:000 samplesΩ:000 samples 直观地说,rejection sampling 会过滤掉不在 Ω\OmegaΩ 中的样本。因此,如果我们从 DDD 的均匀样本开始,我们应该剩下 Ω\OmegaΩ 的均匀样本。 为了形式化我们的推理,让我们推导出采样器的 PDF,表示为 frejf_\text{rej}frej。 为了生成样本 x\mathbf{x}x,我们必须首先从 fDf_DfD 中采样它,然后接受它。 因此 frej(x)f_\text{rej}(\mathbf{x})frej(x) 等价于 fD(x)f_D(\mathbf{x})fD(x) 以接受 x\mathbf{x}x 为条件。 frej(x)=fD∣ accept(x)=P{accept ∣x}fD(x)P{accept}=1⋅1Vol(D)Vol(Ω)Vol(D)=1Vol(Ω) \begin{align*} f_\text{rej}(\mathbf{x}) &= f_{D\ |\ \text{accept}}(\mathbf{x}) \\ &= \frac{\mathbb{P}\left\{\text{accept}\ |\ \mathbf{x}\right\}f_{D}(\mathbf{x})}{\mathbb{P\left\{\text{accept}\right\}}} \tag{Bayes' rule}\\ &= \frac{1\cdot \frac{1}{\text{Vol}(D)}}{\frac{\text{Vol}(\Omega)}{\text{Vol}(D)}} = \frac{1}{\text{Vol}(\Omega)} \tag{$f_D$ is uniform} \end{align*} frej(x)=fD∣accept(x)=P{accept}P{accept∣x}fD(x)=Vol(D)Vol(Ω)1⋅Vol(D)1=Vol(Ω)1(Bayes’ rule)(fD is uniform) Vol\text{Vol}Vol 表示域的体积(在二维中,面积)。 因此,frejf_\text{rej}frej 确实是 Ω\OmegaΩ 上的均匀分布。5
Non-Uniform Rejection Sampling
正如我们在 Monte Carlo integration 中所看到的,rejection sampling 可以直接扩展到与非均匀分布一起使用。 假设 DDD 上的分布的 PDF 是 fD(x)f_D(\mathbf{x})fD(x),我们想用它从 Ω\OmegaΩ 中采样,其 PDF 为 fΩ(x)f_\Omega(\mathbf{x})fΩ(x)。 我们已经知道 Ω⊆D\Omega \subseteq DΩ⊆D,但我们还需要检查一个稍微严格的条件 —— PDF 之间的比率具有有限的上界,表示为 ccc。6 c=supx∈ΩfΩ(x)fD(x) \begin{align*} c = \sup_{\mathbf{x}\in\Omega}\frac{f_\Omega(\mathbf{x})}{f_D(\mathbf{x})} \end{align*} c=x∈ΩsupfD(x)fΩ(x) 上面,我们要求 Ω⊆D\Omega \subseteq DΩ⊆D,因为否则不可能对 Ω\OmegaΩ 的所有部分进行采样。 在这里,我们需要一个有限的 ccc,原因本质上相同 —— 我们要检查 Ω\OmegaΩ 的任何部分都没有被无限不频繁地采样。 一旦我们有了 ccc,我们只需要更新 accept\text{accept}accept。 现在,我们将以概率 fΩ(x)cfD(x)\frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})}cfD(x)fΩ(x) 接受样本 x\mathbf{x}x,该概率始终最多为 111。
def accept(x):
return random(0, 1) < f_Ω(x) / (c * f_D(x))
直观地说,我们通过以与 fΩ(x)fD(x)\frac{f_\Omega(\mathbf{x})}{f_D(\mathbf{x})}fD(x)fΩ(x) 成比例的概率接受 x\mathbf{x}x,将 x\mathbf{x}x 处的概率密度从 fD(x)f_D(\mathbf{x})fD(x) 转换为 fΩ(x)f_\Omega(\mathbf{x})fΩ(x)。 请注意,如果 fDf_DfD 是均匀的,我们直接以与 fΩ(x)f_\Omega(\mathbf{x})fΩ(x) 成比例的概率接受 x\mathbf{x}x。 例如,给定均匀 fDf_DfD 和 fΩ(x)∝11+∣∣x∣∣2f_\Omega(\mathbf{x}) \propto \frac{1}{1+||\mathbf{x}||^2}fΩ(x)∝1+∣∣x∣∣21: D:000 samplesΩ:000 samples\begin{align*} D : \phantom{00}0\text{ samples} \\ \Omega : \phantom{00}0\text{ samples}\end{align*}D:000 samplesΩ:000 samples 正如您所期望的那样,我们看到更多的接受样本靠近中心,其中 fΩ(x)f_\Omega(\mathbf{x})fΩ(x) 最大。 最后,让我们检查一下采样器的 PDF 实际上是 fΩ(x)f_\Omega(\mathbf{x})fΩ(x)。 与上面类似,PDF 等价于 fD(x)f_D(\mathbf{x})fD(x) 以接受 x\mathbf{x}x 为条件。 frej(x)=fD∣ accept(x)=P{accept ∣x}fD(x)P{accept}=fΩ(x)cfD(x)fD(x)∫DfΩ(x)cfD(x)fD(x) dx=fΩ(x)∫DfΩ(x) dx=fΩ(x) \begin{align*} f_\text{rej}(\mathbf{x}) &= f_{D\ |\ \text{accept}}(\mathbf{x}) \\ &= \frac{\mathbb{P}\left\{\text{accept}\ |\ \mathbf{x}\right\}f_{D}(\mathbf{x})}{\mathbb{P}\left\{\text{accept}\right\}} \tag{Bayes' rule}\\ &= \frac{\frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})} f_{D}(\mathbf{x})}{\int_D \frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})} f_D(\mathbf{x}) , d\mathbf{x}} \tag{$\mathbb{P}\left\{\text{accept}\right\} = \mathbb{E}\left[\frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})}\right]$}\\ &= \frac{f_\Omega(\mathbf{x})}{\int_D f_\Omega(\mathbf{x}) , d\mathbf{x}} \tag{Algebra}\\ &= f_\Omega(\mathbf{x}) \tag{$\int_{D\supseteq\Omega} f_\Omega = 1$} \end{align*} frej(x)=fD∣accept(x)=P{accept}P{accept∣x}fD(x)=∫DcfD(x)fΩ(x)fD(x)dxcfD(x)fΩ(x)fD(x)=∫DfΩ(x)dxfΩ(x)=fΩ(x)(Bayes’ rule)(P{accept}=E[cfD(x)fΩ(x)])(Algebra)(∫D⊇ΩfΩ=1) 在第二步中,我们通过计算接受 x\mathbf{x}x 在所有 x∈D\mathbf{x} \in Dx∈D 上的预期概率,来获得接受任意样本的概率。在第四步中,请注意我们将 fΩf_\OmegafΩ 定义为 Ω\OmegaΩ 之外的 000。
Sample Efficiency
许多实际问题仅使用 rejection sampling 和均匀 Monte Carlo integration 即可解决。 选择 DDD 作为包含 Ω\OmegaΩ 的框在任何数量的维度上都有效 —— 框总是易于采样的,因为每个维度都是独立的。7
[WebGL 2 not available]
但是,只有当 fΩf_\OmegafΩ 可以利用 fDf_DfD 中很大比例的概率密度时,rejection sampling 才是有效的。 fΩf_\OmegafΩ 的每个样本都需要几何数量的 fDf_DfD 样本,这些样本根据 P{accept}\mathbb{P}\{\text{accept}\}P{accept} 分布:
P{accept}=E[fΩ(x)cfD(x)]=∫DfΩ(x)cfD(x)fD(x) dx=1c∫DfΩ(x) dx=1c\begin{align*} \mathbb{P}\left\{\text{accept}\right\} &= \mathbb{E}\left[\frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})}\right] \\ &= \int_D \frac{f_\Omega(\mathbf{x})}{cf_D(\mathbf{x})} f_D(\mathbf{x}) , d\mathbf{x} \\ &= \frac{1}{c}\int_D f_\Omega(\mathbf{x}), d\mathbf{x} = \frac{1}{c} \end{align*}P{accept}=E[cfD(x)fΩ(x)]=∫DcfD(x)fΩ(x)fD(x)dx=c1∫DfΩ(x)dx=c1
由于我们有 1c\frac{1}{c}c1 机会接受每个样本,因此我们应该期望为 fΩf_\OmegafΩ 的每个样本生成 ccc 个 fDf_DfD 样本。 直观地说,当 ccc 很大时,这意味着 fDf_DfD 很少对 fΩf_\OmegafΩ 频繁采样的区域进行采样。
例如,当 Ω\OmegaΩ 没有覆盖 DDD 的大部分时,您可能不想使用 rejection sampling:
因此,我们需要设计一种更有效的采样算法。
Inversion Sampling
Inversion sampling 是一种对任何具有可逆累积分布函数 (CDF) 的一维分布进行采样的方法。 随机变量 XXX 的 CDF,表示为 FX(x)F_X(x)FX(x),衡量样本小于 xxx 的概率。 FX(x)=P{X<x}=∫−∞xfX(x′) dx′\begin{align*} F_X(x) &= \mathbb{P}\{X < x\}\\ &= \int_{-\infty}^x f_X(x^\prime), dx^\prime \end{align*}FX(x)=P{X<x}=∫−∞xfX(x′)dx′ 直观地说,CDF 将 xxx 映射到低于 xxx 的概率质量百分比: 因此,inverse CDF FX−1(p)F^{-1}X(p)FX−1(p) 将概率质量的百分比映射到相应的 xxx。8 我们可以通过从 [0,1][0,1][0,1] 中均匀采样 ppp 并计算 FX−1(p)F_X^{-1}(p)FX−1(p) 来定义一个 inversion sampler Inv\text{Inv}Inv。 为了表征采样器的行为,我们可以找到它的 PDF。 Inv=FX−1(p)\text{Inv} = F_X^{-1}(p)Inv=FX−1(p) 落入 dxdxdx 范围内的概率等价于 ppp 落入相应 dpdpdp 范围内的概率: 我们均匀地对 ppp 进行采样,因此 ppp 落入 dpdpdp 的概率是 dpdpdp 的长度。 然后,dxdxdx 上的平均概率密度为 dpdx\frac{dp}{dx}dxdp。 在极限情况下,dpdpdp 的长度与 FXF_XFX 的斜率成正比 —— 该比率是它的导数! fInv(x)=dpdx=dFX(x)dx=fX(x)\begin{align*} f\text{Inv}(x) &= \frac{dp}{dx}\\ &= \frac{dF_X(x)}{dx} \tag{$p = F_X(x)$}\\ &= f_X(x) \end{align*}fInv(x)=dxdp=dxdFX(x)=fX(x)(p=FX(x)) 因此,我们有 fInv=fXf_\text{Inv} = f_XfInv=fX。9 这意味着 inversion sampling 有效,但更严格地说,我们可以检查采样器的 CDF 是否与 FXF_XFX 匹配。 FInv(x)=P{Inv<x}=P{FX−1(Uniform(0,1))<x}=P{FX(FX−1(Uniform(0,1)))<FX(x)}=P{Uniform(0,1)<FX(x)}=FX(x) \begin{align*} F_\text{Inv}(x) &= \mathbb{P}\{\text{Inv} < x\}\\ &= \mathbb{P}\{F_X^{-1}(\text{Uniform}(0,1)) < x\}\\ &= \mathbb{P}\{F_X(F_X^{-1}(\text{Uniform}(0,1))) < F_X(x)\} \tag{$F_X$ non-decreasing}\\ &= \mathbb{P}\{\text{Uniform}(0,1) < F_X(x)\} \tag{Def. $F_X^{-1}$}\\ &= F_X(x) \tag{$F_X(x) \in [0,1]$} \end{align*} FInv(x)=P{Inv<x}=P{FX−1(Uniform(0,1))<x}=P{FX(FX−1(Uniform(0,1)))<FX(x)}=P{Uniform(0,1)<FX(x)}=FX(x)(FX non-decreasing)(Def. FX−1)(FX(x)∈[0,1]) 由于它们的 CDF 是等价的,因此我们确实有 Inv∼X\text{Inv} \sim XInv∼X。
Marginal Inversion Sampling
如前所述,inversion sampling 仅适用于一维分布。 幸运的是,我们可以通过迭代采样每个维度的 marginal distribution 将 inversion 扩展到更高的维度。 让我们推导出二维分布 fXY(x,y)f_{XY}(x,y)fXY(x,y) 的 inversion sampler。 首先,我们将定义 marginal distribution fX(x)f_X(x)fX(x),它计算所有 yyy 选择中 xxx 处的总概率密度。 此分布是一维的,因此我们可以使用 inversion sampling 来选择样本 XXX。 其次,我们将计算以 XXX 为条件的 marginal distribution fY(y)f_Y(y)fY(y),该分布必须与 fXY(X,y)f_{XY}(X,y)fXY(X,y) 成比例。 fY∣X(y)=f(X,y)∫−∞∞f(X,y) dy=f(X,y)fX(X)\begin{align*} f_{Y\ |\ X}(y) &= \frac{f(X,y)}{\int_{-\infty}^\infty f(X,y), dy}\\ &= \frac{f(X,y)}{f_X(X)} \end{align*}fY∣X(y)=∫−∞∞f(X,y)dyf(X,y)=fX(X)f(X,y) 最后,我们可以再次应用 inversion 来采样 YYY。 直观地说,fY∣Xf_{Y\ |\ X}fY∣X 根据我们对 XXX 的选择为 yyy 选择正确的分布。 我们将在下一节中更严格地探讨 inversion sampling 的工作原理。
Changes of Coordinates
虽然 marginal inversion sampling 可以构建任意高维分布,但在实践中通常没有必要。 这是因为 inversion sampling 是变换随机变量的更通用技术的特殊情况。
为了说明这一点,让我们尝试为单位圆盘定义一个均匀采样器。 与 rejection sampling 不同,我们首先需要选择域的 parameterization。 一个自然的选择是极坐标,其中 θ\thetaθ 是相对于 x 轴的角度,rrr 是距原点的距离。
Φ(r,θ)=(rcosθ,rsinθ)Φ−1(x,y)=(x2+y2,atan2(y,x))\begin{align*} \Phi(r,\theta) &= (r\cos\theta, r\sin\theta)\\ \Phi^{-1}(x,y) &= (\sqrt{x^2+y^2}, \text{atan2}(y,x)) \end{align*}Φ(r,θ)Φ−1(x,y)=(rcosθ,rsinθ)=(x2+y2,atan2(y,x))
因此,单位圆盘由 S=Φ(R)\mathcal{S} = \Phi(\mathcal{R})S=Φ(R) 描述,其中 R=[0,1]×[0,2π]\mathcal{R} = [0,1]\times[0,2\pi]R=[0,1]×[0,2π]。 为了生成 S\mathcal{S}S 的采样器,我们可以尝试将 R\mathcal{R}R 的均匀样本映射到 S\mathcal{S}S 上:
但是 R\mathcal{R}R 的均匀样本不会变成 S\mathcal{S}S 的均匀样本。 这是因为从极坐标到矩形坐标的转换并没有保留面积 —— 较小的半径包含的面积较小,但我们对所有 rrr 的权重都相同。
为了确定哪里出了问题,让我们找到此采样器的 PDF。 关键的观察是,样本 s\mathbf{s}s 落入圆形块 dSd\mathcal{S}dS 的充要条件是 r=Φ−1(s)\mathbf{r} = \Phi^{-1}(\mathbf{s})r=Φ−1(s) 落入相应的矩形 dRd\mathcal{R}dR 中。10
因此,采样任一区域的概率必须是等价的。
∫dSfS(s) ds=∫dRfR(r) dr \int_{d\mathcal{S}} f_\mathcal{S}(\mathbf{s}), d\mathbf{s} = \int_{d\mathcal{R}} f_\mathcal{R}(\mathbf{r}), d\mathbf{r} ∫dSfS(s)ds=∫dRfR(r)dr
在极限情况下,这些积分简化为相应的 PDF 乘以块的面积。
fS(s)⋅∣dS∣=fR(r)⋅∣dR∣ ⟹ fS(s)=fR(r)⋅∣dR∣∣dS∣\begin{align*} && f_\mathcal{S}(\mathbf{s})\cdot|d\mathcal{S}| &= f_\mathcal{R}(\mathbf{r})\cdot|d\mathcal{R}|\\ &\implies& f_\mathcal{S}(\mathbf{s}) &= f_\mathcal{R}(\mathbf{r})\cdot\frac{|d\mathcal{R}|}{|d\mathcal{S}|} \end{align*}⟹fS(s)⋅∣dS∣fS(s)=fR(r)⋅∣dR∣=fR(r)⋅∣dS∣∣dR∣
直观地说,面积比 ∣dR∣∣dS∣\frac{|d\mathcal{R}|}{|d\mathcal{S}|}∣dS∣∣dR∣ 告诉我们 dSd\mathcal{S}dS 在映射到 dRd\mathcal{R}dR 上时被压缩或拉伸了多少。 例如,如果 dSd\mathcal{S}dS 按两倍的系数缩小,则 dRd\mathcal{R}dR 必须包含两倍的概率密度。
最后,由于 Φ−1\Phi^{-1}Φ−1 将 S\mathcal{S}S 映射到 R\mathcal{R}R,因此面积缩放因子由其导数给出:11
fS(s)=fR(r)⋅∣dR∣∣dS∣=fR(Φ−1(s))⋅∣dΦ−1(S)dS∣=fR(Φ−1(s))⋅∣DΦ−1∣\begin{align*} f_\mathcal{S}(\mathbf{s}) &= f_\mathcal{R}(\mathbf{r})\cdot\frac{|d\mathcal{R}|}{|d\mathcal{S}|} \\ &= f_\mathcal{R}(\Phi^{-1}(\mathbf{s}))\cdot\left|\frac{d\Phi^{-1}(\mathcal{S})}{d\mathcal{S}}\right| \tag{$\mathcal{R} = \Phi^{-1}(\mathcal{S})$}\\ &=f_\mathcal{R}(\Phi^{-1}(\mathbf{s}))\cdot |D\Phi^{-1}| \end{align*}fS(s)=fR(r)⋅∣dS∣∣dR∣=fR(Φ−1(s))⋅dSdΦ−1(S)=fR(Φ−1(s))⋅∣DΦ−1∣(R=Φ−1(S))
其中 ∣DΦ−1∣|D\Phi^{-1}|∣DΦ−1∣ 表示 Φ−1\Phi^{-1}Φ−1 的 Jacobian 的行列式。12
Sampling via Change of Coordinates
现在我们知道了 fSf_\mathcal{S}fS 和 fRf_\mathcal{R}fR 之间的关系,我们可以选择一个不同的、非均匀的 fRf_\mathcal{R}fR,它将产生一个均匀的 fSf_\mathcal{S}fS。 我们新的 PDF 将需要消除 ∣DΦ−1∣|D\Phi^{-1}|∣DΦ−1∣ 的因子: ∣DΦ−1∣=∣[δrδxδθδxδrδyδθδy]∣=δrδx⋅δθδy−δθδx⋅δrδy=xx2+y2⋅xx2+y2+yx2+y2⋅yx2+y2=1x2+y2=1r\begin{align*} |D\Phi^{-1}| &= \left|\begin{bmatrix}\frac{\delta r}{\delta x} & \frac{\delta \theta}{\delta x}\\ \frac{\delta r}{\delta y} & \frac{\delta \theta}{\delta y}\end{bmatrix}\right|\\ &= \frac{\delta r}{\delta x}\cdot\frac{\delta \theta}{\delta y} - \frac{\delta \theta}{\delta x}\cdot\frac{\delta r}{\delta y} \\ &= \frac{x}{\sqrt{x^2+y^2}} \cdot \frac{x}{x^2+y^2} + \frac{y}{x^2+y^2} \cdot \frac{y}{\sqrt{x^2+y^2}}\\ &= \frac{1}{\sqrt{x^2+y^2}} = \frac{1}{r} \end{align*}∣DΦ−1∣=[δxδrδyδrδxδθδyδθ]=δxδr⋅δyδθ−δxδθ⋅δyδr=x2+y2x⋅x2+y2x+x2+y2y⋅x2+y2y=x2+y21=r1 与 1r\frac{1}{r}r1 成比例是有道理的 —— 我们行为不端的采样器在原点附近产生了太多的样本。 如果我们改为根据 fR(r,θ)=r2πf_\mathcal{R}(r,\theta) = \frac{r}{2\pi}fR(r,θ)=2πr 对 R\mathcal{R}R 进行采样,我们将最终得到一个均匀的 fSf_\mathcal{S}fS。13 fS(x,y)=fR(Φ−1(x,y))⋅∣DΦ−1∣=r2π⋅1r=12π\begin{align*} f_\mathcal{S}(x,y) &= f_\mathcal{R}(\Phi^{-1}(x,y))\cdot|D\Phi^{-1}| \\ &= \frac{r}{2\pi}\cdot \frac{1}{r} = \frac{1}{2\pi} \end{align*}fS(x,y)=fR(Φ−1(x,y))⋅∣DΦ−1∣=2πr⋅r1=2π1 在前一节中,我们应用了一维的坐标变换。 也就是说,通过取 Φ=FX−1\Phi = F_X^{-1}Φ=FX−1,我们将均匀单位分布转换为具有我们期望的 PDF。 fInv(x)=fU(FX−1(x))⋅∣D(FX−1)−1∣=1⋅∣DFX(x)∣=fX(x)\begin{align*} f_{\text{Inv}}(x) &= f_U(F_X^{-1}(x)) \cdot \left|D\left(F_X^{-1}\right)^{-1}\right|\\ &= 1 \cdot |DF_X(x)|\\ &= f_X(x) \end{align*}fInv(x)=fU(FX−1(x))⋅D(FX−1)−1=1⋅∣DFX(x)∣=fX(x) 在实践中,许多有用的分布都可以通过适当的坐标变换来有效地采样。 然而,这样做需要对域进行参数化,这有时是无法实现的。 在这种情况下,我们可以转向诸如 Markov Chain Monte Carlo 之类的方法,这将在以后的章节中讨论。
Footnotes
- 您个人可以[生成随机数](https://thenumb.at/Sampling/https:/calmcode.io/blog/inverse-turing-test)吗? 很难说。↩︎
- 另一方面,Cryptography 可能需要真正的随机数。 有硬件设备可以根据物理传感器噪声甚至量子过程提供样本。↩︎
- 虽然我们在第二章中假设了独立样本,但这不是 Monte Carlo 工作所必需的! 我们将在以后的章节中探讨相关随机数的用法。↩︎
- 例如,如果 PRNG 具有有限数量的内部状态,则该序列必然在它们耗尽后重复。 同样,从当前状态计算下一个样本的算法的存在意味着先前的样本提供有关未来样本的_一些_信息。 良好的 PRNG 确保计算上难以利用该信息。 有关更多分析,请参阅 PCG paper。↩︎
- 请注意,贝叶斯规则的非标准用法用于计算以离散事件为条件的概率密度。 此步骤应理解如下:在 x\mathbf{x}x 处,以 x\mathbf{x}x 被接受为条件的密度等价于接受 x\mathbf{x}x 的概率(即 1),乘以在 x\mathbf{x}x 处的密度,再除以我们接受任意样本的概率(即体积比率)。↩︎
- 请注意,ω(x)d(x)\frac{\omega(\mathbf{x})}{d(\mathbf{x})}d(x)ω(x) 可能没有最大值; 我们实际上使用 Ω\OmegaΩ 上的上确界。↩︎
- 也就是说,在更高的维度上,球体——以及大多数您想要采样的其他域——所包含的体积比立方体小得多! 维度诅咒再次出现。↩︎
- 假设 CDF 是单射的。 某些分布(例如,只要 fX=0f_X = 0fX=0 在某个区间上)没有单射 CDF,这在通常意义上是不可逆的。 幸运的是,我们仍然可以将 inversion sampling 应用于它们的广义逆,定义为 FX−1(p)=infx{FX(x)≥p}F_X^{-1}(p) = \inf_x\{F_X(x) \ge p\}FX−1(p)=infx{FX(x)≥p}。 直观地说,广义逆会跳过 fX=0f_X = 0fX=0 的区域。
↩︎
- 假设 CDF 是可微的。 正如第一章中所述,涉及 Dirac deltas 的分布可能具有不连续的(因此不可微的)CDF。 幸运的是,inversion sampling 仍然有效:广义逆在对应于不连续性的区间上只是常数。 使用 [non-st