蒙特卡洛积分与重要性采样详解
最近在看有關蒙特卡洛積分的內(nèi)容,發(fā)現(xiàn)網(wǎng)上很多博主寫的證明過程跳步較為嚴重,而且過程晦澀,不太容易理解。我在自己閱讀國外相關教材附錄后發(fā)現(xiàn)證明蒙特卡洛積分方法并不難,利用的僅是概率論的基本知識,現(xiàn)整理下來與大家分享。
那么什么是蒙特卡洛積分?簡而言之就是,在求積分時,如果找不到被積函數(shù)的原函數(shù),那么利用經(jīng)典積分方法是得不到積分結(jié)果的,但是蒙特卡洛積分方法告訴我們,利用一個隨機變量對被積函數(shù)進行采樣,并將采樣值進行一定的處理,那么當采樣數(shù)量很高時,得到的結(jié)果可以很好的近似原積分的結(jié)果。這樣一來,我們就不用去求原函數(shù)的形式,就能求得積分的近似結(jié)果。
一、前提知識:
由概率論基本知識,假設一連續(xù)型隨機變量$X$的樣本空間為$D$,其概率密度分布函數(shù)為$p(x)$,則其數(shù)學期望為
\begin{equation}
E\left[X\right]=\int_D {xp(x){\rmze8trgl8bvbq}x}
\end{equation}
若另一連續(xù)隨機變量Y滿足$Y=f(X)$,則$Y$的數(shù)學期望$E[Y]$可由下式給出
\begin{equation}
E\left[Y\right]=\int_D {f(x)p(x){\rmze8trgl8bvbq}x}
\end{equation}
二、蒙特卡洛積分與重要性采樣
根據(jù)以上敘述,假設這里我們要計算一個一維積分式
\begin{equation}
A=\int_a^b {f(x){\rmze8trgl8bvbq}x}
\end{equation}
根據(jù)經(jīng)典的方法,我們需要求得$f(x)$的原函數(shù)$F(x)$,才能解出這個積分結(jié)果,但如果$f(x)$的原函數(shù)形式復雜,或者根本求不出來,總之在不知道$F(x)$的具體形式的情況下,如果我們還想計算這個積分,怎么辦?這時候我們就需要借助蒙特卡洛積分(Monte Carlo Integration)方法。蒙特卡洛積分方法告訴我們,為求得積分結(jié)果,可以構(gòu)造
\begin{equation}
\label{eq:FN}
F_N = \frac{{b - a}}{N}\sum\limits_{i = 1}^N {f(X_i )}
\end{equation}
其中的每一個$X_i(i=1,2,3,...,N)$為$[a,b]$之間的均勻連續(xù)隨機變量,所有的$X_i$組成一個隨機變量集合。下面證明,$F_N$的數(shù)學期望即為我們要求的結(jié)果$A$。
\begin{equation}
\begin{split}
\label{eq:proof1}
E\left[ {F_N } \right] &= E\left[ {\frac{{b - a}}{N}\sum\limits_{i = 1}^N {f(X_i )} } \right] \\
&=\frac{{b - a}}{N}\sum\limits_{i = 1}^N {E\left[ {f(X_i )} \right]} \\
&= \frac{{b - a}}{N}\sum\limits_{i = 1}^N {\int_a^b {f(x)p(x){\rmze8trgl8bvbq}x} } \\
&= \frac{{b - a}}{N}\sum\limits_{i = 1}^N {\int_a^b {f(x)\frac{1}{{b - a}}{\rmze8trgl8bvbq}x} } \\
&= \frac{{b - a}}{N}\frac{1}{{b - a}}\sum\limits_{i = 1}^N {\int_a^b {f(x){\rmze8trgl8bvbq}x} } \\
&= \frac{1}{N}\sum\limits_{i = 1}^N {\int_a^b {f(x){\rmze8trgl8bvbq}x} } \\
&= \int_a^b {f(x){\rmze8trgl8bvbq}x}
\end{split}
\end{equation}
以上證明過程表明,若我們根據(jù)式(\ref{eq:FN})來構(gòu)造一個新的隨機變量$F_N$,則$F_N$的期望就是積分$\int_a^b{f(x){\rmze8trgl8bvbq}x}$的結(jié)果,隨著$N$的增加,$F_N$就越逼近理論上$A$的值,即$F_N$是$A$的一個無偏估計。這樣我們實際上就避開了求$f(x)$的原函數(shù)$F(x)$的過程。整個求積分近似值的過程則可以用文字描述如下:首先從區(qū)間$[a,b]$上對均勻分布的隨機變量$X$連續(xù)取樣$N$次,得到$N$個取樣值$\{x_1,x_2,x_3,...,x_N\}$,對每個取樣值$x_i(i=1,2,3,...,N)$計算$f(x_i)$得到$\{f(x_1),f(x_2),f(x_3),...,f(x_N)\}$,再計算它們的和$\sum\limits_{i=1}^{N}{f(x_i)}$,最后乘系數(shù)$\frac{b-a}{N}$即可得到對理論值的一個估計。
進一步對上述過程分析,我們發(fā)現(xiàn)這里的$X$被規(guī)定為與原積分區(qū)間相同的均勻分布隨機變量。那么對于與原積分區(qū)間相同,但卻不是均勻分布的一般隨機變量,上述結(jié)論是否仍成立?結(jié)論是,如果這里的隨機變量$X$的概率密度分布函數(shù)已知,那么上述結(jié)論還是成立的。假設其概率密度分布函數(shù)為$p(x)$,證明如下
類似式(\ref{eq:FN}),這次我們構(gòu)造
\begin{equation}
\label{eq:FN2}
F_N = \frac{1}{N}\sum\limits_{i = 1}^N {\frac{{f(X_i )}}{{p(X_i )}}}
\end{equation}
再構(gòu)造隨機變量
$$Y=\frac{f(X)}{p(X)}$$
式(\ref{eq:FN2})的所有量都是已知的。則$F_N$的數(shù)學期望$E\left[F_N\right]$即為
\begin{equation}
\begin{split}
\label{eq:proof2}
E\left[ {F_N } \right] &= E\left[ {\frac{1}{N}\sum\limits_{i = 1}^N {\frac{{f(X_i )}}{{p(X_i )}}} } \right] \\
&= \frac{1}{N}\sum\limits_{i = 1}^N {E\left[ {\frac{{f(X_i )}}{{p(X_i )}}} \right]} \\
&= \frac{1}{N}\sum\limits_{i = 1}^N {E\left[ {Y_i } \right]} \\
&= \frac{1}{N}\sum\limits_{i = 1}^N {\int_a^b {\left[ {\frac{{f(x)}}{{p(x)}}} \right]p(x){\rmze8trgl8bvbq}x} } \\
&= \frac{1}{N}\sum\limits_{i = 1}^N {\int_a^b {f(x){\rmze8trgl8bvbq}x} } \\
&= \int_a^b {f(x){\rmze8trgl8bvbq}x}
\end{split}
\end{equation}
我們發(fā)現(xiàn)式(5)其實是式(7)的特殊形式。與之前的要求不同,這里僅要求隨機變量$X$的概率密度分布函數(shù)$p(x)$已知且在$X$的樣本空間內(nèi)$p(x)\ne0$。綜合上述敘述,我們可以得到蒙特卡洛積分方法如下
\begin{equation}
\int_D {f(x){\rmze8trgl8bvbq}x} = \mathop {\lim }\limits_{N \to \infty } \frac{1}{N}\sum\limits_{i = 1}^N {\frac{{f(X_i )}}{{p(X_i )}}}
\end{equation}
我們可以進一步地將求解一維積分的方法擴展到求解高維積分中去??紤]求解如下積分
\begin{equation}
\int_{z_1 }^{z_2 } {\int_{y_1 }^{y_2 } {\int_{x_1 }^{x_2 } {f(x,y,z){\rmze8trgl8bvbq}x{\rmze8trgl8bvbq}y{\rmze8trgl8bvbq}z} } }
\end{equation}
只需要在積分域$\int_{z_1 }^{z_2 } {\int_{y_1 }^{y_2 } {\int_{x_1 }^{x_2 } } } $定義的盒型區(qū)域內(nèi)選取概率密度分布為
$$
p(x) = \frac{1}{{\left( {x_2 - x_1 } \right)\left( {y_2 - y_1 } \right)\left( {z_2 - z_1 } \right)}}
$$
的均勻隨機變量$X$,則積分結(jié)果為
$$
\frac{{\left( {x_2 - x_1 } \right)\left( {y_2 - y_1 } \right)\left( {z_2 - z_1 } \right)}}{N}\sum\limits_{i = 1}^N {f\left( {X_i } \right)}
$$
現(xiàn)在我們分析一下隨著樣本數(shù)量$N$的增加,估計值$F_N$的方差$\sigma ^2 \left[ {F_N } \right]$的變化情況,以便得出蒙特卡洛積分方法的收斂速度特性。在式(7)的基礎上,我們繼續(xù)計算$\sigma ^2 \left[ {F_N } \right]$如下
\begin{equation}
\begin{split}
\label{eq:proof3}
\sigma ^2 \left[ {F_N } \right] &= \sigma ^2 \left[ {\frac{1}{N}\sum\limits_{i = 1}^N {\frac{{f(X_i )}}{{p(X_i )}}} } \right] \\
&= \frac{1}{{N^2 }}\sum\limits_{i = 1}^N {\sigma ^2 \left[ {\frac{{f(X_i )}}{{p(X_i )}}} \right]} \\
&= \frac{1}{{N^2 }}\sum\limits_{i = 1}^N {\sigma ^2 \left[ {Y_i } \right]} \\
&= \frac{1}{{N^2 }}\left( {N\sigma ^2 \left[ Y \right]} \right) \\
&= \frac{1}{N}\sigma ^2 \left[ Y \right]
\end{split}
\end{equation}
所以
\begin{equation}
\sigma \left[ {F_N } \right] = \frac{1}{{\sqrt N }}\sigma \left[ Y \right]
\end{equation}
上式告訴我們,估計值的不穩(wěn)定來源于隨機變量$Y$的取值不穩(wěn)定。換句話說,如果$\frac{{f(X_i )}}{{p(X_i )}}$因不同$X_i$的取值變化地越劇烈,就會造成$Y$的方差較大,則會造成估計值的收斂速度越慢。這啟示我們,若$p(x)$的形狀越接近$f(x)$,則有益于最終結(jié)果的收斂。上述思想即為“重要性采樣”方法,即對積分值有重要貢獻($f(x)$較大)的被積函數(shù)區(qū)間,我們以較大概率生成處于這個區(qū)間附近的隨機變量,用于快速逼近理論值。這也就是為什么我們要引入任意分布隨機變量的蒙特卡洛積分方法,而不滿足于利用均勻分布隨機變量來求蒙特卡洛積分的原因。
利用蒙特卡洛方法,我們可以得到任意一個積分的結(jié)果,但是代價就是我們得不到其理論值,我們得到的只是一個對理論值的估計,估計值與理論值之間的誤差可以通過增加樣本數(shù)來減小,但收斂速率僅為$O\left( {\sqrt N } \right)$,也就是說,若想將誤差降為現(xiàn)在的一半,我們需要再多計算$4$倍的計算量才可以達到。即便如此,原始的蒙特卡洛積分方法也不失為是一種經(jīng)典有效的方法。
?
轉(zhuǎn)載于:https://www.cnblogs.com/time-flow1024/p/10094293.html
總結(jié)
以上是生活随笔為你收集整理的蒙特卡洛积分与重要性采样详解的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 福大软工1816:项目测评
- 下一篇: 河南省第十一届ACM程序设计竞赛 修路