1. 배경

• 베이즈 정리를 이용하여 사후 분포를 구하는 과정은 다음과 같다.

$$ P(\theta|D) = \frac{P(D|\theta)P(\theta)}{\int P(D|\theta)P(\theta)\,d\theta} $$

• 문제는 파라미터인 $\theta$가 다차원인 경우, 적분을 수식으로 구하기가 쉽지 않다.
• 이 경우 적분을 정확히 구할 수는 없고 추정해야 하는데, 이 때 몬테카를로 적분이 쓰인다.

2. 함수의 적분값 구하기

• 임의의 함수 g(x)를 적분한다고 하자. 즉, $\int g(x)\,dx$ 를 구한다고 하자.
• 약간의 변형을 거쳐서 다음처럼 계산할 수 있다. 여기서 f(x)는 확률분포이다.

$$ \begin{align} \int g(x)dx = \int \frac{g(x)}{f(x)}f(x)dx = E_{x \sim f(x)}\left(\frac{g(x)}{f(x)}\right) \approx \frac{1}{N}\sum_{n=1}^{N}\frac{g(x_{n})}{f(x_{n})}

\end{align} $$

• 만약 표본을 추출하는 분포 f를 적분 영역 $dom(f)$를 포함하는 균등 분포로 설정하면, 위의 수식은 다음과 같아진다.

$$ \begin{align} \frac{1}{N}\sum_{n=1}^{N}\frac{g(x_{n})}{f(x_{n})} = \frac{1}{N}\sum_{n=1}^{N}\frac{g(x_{n})}{\frac{1}{(b-a)}} = (b-a)\sum_{n=1}^{N}\frac{g(x_{n})}{N} \end{align} $$

• 코드로 구현하면 다음과 같다. 예시로 사용한 함수는 $\int_{0}^{2}3x^{2}dx$이다.

import numpy as np
import math

def calc_integral_of_func(uniform_samples, func, start, end):
		func_output = [func(i) for i in uniform_samples]
		integral = (end - start) * sum(func_output) / len(func_output)
		
    return integral
    
a, b, N = 0, 2, 10000
uniform_samples = np.random.uniform(a, b, N)
func = labda x : 3*math.pow(x, 2)
calc_integral_of_func(uniform_samples, func, a, b)