蒙地卡羅積分

在數學中，蒙地卡羅積分（Monte Carlo integration）是一種使用亂數進行數值積分的技術。它是一種特殊的蒙地卡羅方法，可對定積分進行數值計算。其他演算法通常在規則網格上評估被積函數，而^[1]蒙地卡羅隨機選擇被積函數評估的點。 ^[2]該方法對於高維積分特別有用。 ^[3]

蒙地卡羅積分示意圖。在該範例中，D是和正方形內切的圓，域E則是正方形。因為正方形的面積 (值是4) 很容易計算，所以圓的面積 (π*1.0²) 可以通過圓內的點 (40個紅點) 與總點數 (50個點) 的比率(0.8) 來估計，得出圓面積的近似值 4*0.8 = 3.2 ≈ π。

進行蒙地卡羅積分的方法有多種，例如均勻採樣、分層採樣、重要性採樣、順序蒙地卡羅（也稱為粒子濾波器）和平均場粒子方法。

概述

在數值積分中，梯形法則等方法使用確定性方法。另一方面，蒙地卡羅積分採用非確定性方法：每種實現都提供不同的結果。在蒙地卡羅中，最終結果是帶有各自誤差線的正確值的近似值，並且正確值很可能在這些誤差線內。

考慮函數$${\displaystyle I=\int _{\Omega }f({\overline {\mathbf {x} }})\,d{\overline {\mathbf {x} }}}$$其中 Ω 是${\displaystyle \mathbb {R} ^{m}}$的一個子集，它的體積是：$${\displaystyle V=\int _{\Omega }d{\overline {\mathbf {x} }}}$$

樸素的蒙地卡羅方法是在Ω上均勻採樣點： ^[4] 給定N個均勻樣本，$${\displaystyle {\overline {\mathbf {x} }}_{1},\cdots ,{\overline {\mathbf {x} }}_{N}\in \Omega ,}$$I可以被近似為：$${\displaystyle I\approx Q_{N}\equiv V{\frac {1}{N}}\sum _{i=1}^{N}f({\overline {\mathbf {x} }}_{i})=V\langle f\rangle .}$$ 這是因為大數法則確保以下結論成立：$${\displaystyle \lim _{N\to \infty }Q_{N}=I.}$$ 我們使用I給出Q_N的I估計， Q_N的誤差線（error bars）可以使用變異數的不偏估計通過樣本變異數來估計。$${\displaystyle \mathrm {Var} (f)\equiv \sigma _{N}^{2}={\frac {1}{N-1}}\sum _{i=1}^{N}\left(f({\overline {\mathbf {x} }}_{i})-\langle f\rangle \right)^{2}.}$$由此可得：$${\displaystyle \mathrm {Var} (Q_{N})={\frac {V^{2}}{N^{2}}}\sum _{i=1}^{N}\mathrm {Var} (f)=V^{2}{\frac {\mathrm {Var} (f)}{N}}=V^{2}{\frac {\sigma _{N}^{2}}{N}}.}$$只要以下序列$${\displaystyle \left\{\sigma _{1}^{2},\sigma _{2}^{2},\sigma _{3}^{2},\ldots \right\}}$$是有界的，該變異數就會逐漸減小到零，即 1/N 。 那麼，Q_N的誤差估計就是：$${\displaystyle \delta Q_{N}\approx {\sqrt {\mathrm {Var} (Q_{N})}}=V{\frac {\sigma _{N}}{\sqrt {N}}},}$$

我們看到，Q_N的誤差是自身的${\displaystyle {\tfrac {1}{\sqrt {N}}}}$，這是平均值的標準誤差乘以${\displaystyle V}$。該結果不依賴於積分的維數，這是蒙地卡羅積分相對於大多數指數依賴維數的確定性方法所具有的優勢。 值得注意的是，與確定性方法不同，誤差的估計不是嚴格的誤差界限；隨機抽樣可能無法揭示被積函數的所有重要特徵，從而導致誤差的低估。

雖然樸素蒙地卡羅適用於簡單的範例，但對確定性演算法的改進只能通過使用特定於問題的採樣分布的演算法來實現。通過適當的樣本分布，可以利用以下事實：幾乎所有高維被積函數都是非常局部化的，並且只有小子空間對積分有顯著貢獻。 ^[5]蒙地卡羅文獻的很大一部分致力於開發改進誤差估計的策略。特別是，分層採樣（將區域劃分為子體）和重要性採樣（從非均勻分布中採樣）是此類技術的兩個範例。

例子

使用蒙地卡羅方法最典型的例子是估計π。

我們考慮以下函數：$${\displaystyle H\left(x,y\right)={\begin{cases}1&{\text{if }}x^{2}+y^{2}\leq 1\\0&{\text{else}}\end{cases}}}$$即正方形集合Ω=[−1,1]×[−1,1]的體積是V = 4。 因此，使用蒙地卡羅積分計算π值的粗略方法是在Ω上選取N個亂數並計算：$${\displaystyle Q_{N}=4{\frac {1}{N}}\sum _{i=1}^{N}H(x_{i},y_{i})}$$右圖中，相對誤差${\displaystyle {\tfrac {Q_{N}-\pi }{\pi }}}$被認為是一個關於N的函數，由此確認${\displaystyle {\tfrac {1}{\sqrt {N}}}}$ 。

C 範例

請記住，以下的程式需要真正的亂數生成器。

int i, throws = 99999, insideCircle = 0;
double randX, randY, pi;

srand(time(NULL));

for (i = 0; i < throws; ++i) {
  randX = rand() / (double) RAND_MAX;
  randY = rand() / (double) RAND_MAX;
  if (randX * randX + randY * randY < 1) ++insideCircle;
}

pi = 4.0 * insideCircle / throws;

Python範例

from numpy import random
import numpy as np

throws = 2000
inside_circle = 0
i = 0
radius = 1
while i < throws:
    # Choose random X and Y centered around 0,0
    x = random.uniform(-radius, radius)
    y = random.uniform(-radius, radius)
    # If the point is inside circle, increase variable
    if x**2 + y**2 <= radius**2:
        inside_circle += 1
    i += 1

# Calculate area and print; should be closer to Pi with increasing number of throws
area = (((2 * radius) ** 2) * inside_circle) / throws
print(area)

參見

• 擬蒙地卡羅方法
• 蒙地卡羅輔助場
• 統計物理中的蒙地卡羅方法
• 蒙地卡羅方法
• 減小變異數（英語：Variance reduction）