컴퓨터의 연산능력을 이용하여 수식만으로 계산하기 어려운 문제가 있을 때, 데이터의 무작위 샘플을 얻고 이를 이용하여 해답을 구하는 방법이다.
예를 들면, 원의 원주율 $\pi$를 구한다고 하자. 샘플링만으로 이 값을 구하기 위하여, 정사각형 내에 반지름이 1인 내접하는 원을 그린 후, 무작위 샘플을 발생시킨다. 그리고 이 샘플이 원 내부에 있는지, 외부인지 판단 만으로 값을 구할 수 있다. 즉, $2^2 - \pi \cdot 1^2 = $ samples between two areas 로 구할 수 있다.
이때 고성능 컴퓨터를 이용하여 샘플을 많이 발생 시킬수록 더 정확한 값을 구할 수 있다.
Markov Chain
연속으로 이어지는 일련의 값이 주어질 때, 다음에 나오는 값은 바로 직전(이전) 값에만 영향을 받는 연속 값의 체인이다. 이 때 값은 상태(state)를 말하고 체인은 상태 값의 연속 시퀀스이다.
값들은 상호 독립적이지 않고 서로 상관성이 존재하는데 연속된 값 사이에서만 그 상관성이 존재하고 멀리 떨어진 값들 사이의 상관성은 무시된다.
어떤 수식의 값을 샘플링 만으로 근사하여 구하며, 이때 샘플을 발생시키는 방법을 상호 독립적이지 않고 어떤 관계를 가지도록 발생 시킨다는 의미로 사용된다. 따라서 이름이 MCMC 근사로 붙여졌다.
고차공간에서 확률분포 추정에 사용하며 왜 MCMC인지는 아래 그림 하단부(MH 알고리즘)에서 다시 설명한다.
MCMC(Markov Chain Monte Carlo)는 고차원 공간에서 수행해야 하는 적분이나 최적화 문제를 풀기 위해 사용된다. 주요 적용 분야는 기계학습, 물리, 통계학, 경제학, 결정이론 등이다.
고차공간에서 확률분포 $p(z)$의 추정은 매우 어렵다. 확률분포(probability density function: pdf)는
$p(z)={\tilde{p}(z) \over Z_p}$, where $Z_p = \int_Z {\tilde{p}(z) dz}$
와 같이 계산한다. 여기서 $\tilde{p}(z)$는 정규화가 안된 함수 값이고, $Z_p$는 정규화 상수(normalization constant)이다. 보통 $\tilde{p}(z)$는 쉽게 구할 수 있지만, $Z_p$의 경우는 $Z$ 공간 내의 가능한 모든 영역의 $z$에 대해 적분을 필요로 하며, 적분할 수식의 형태도 주어지지 않는 경우가 많아 어렵다. 이 때 sampling기법을 이용하는 MCMC는 $p(z)$의 추정에 유용한 도구이다.
MCMC가 사용되는 몇 가지 분야를 살펴 보자.
(1) Bayesian inference and learning: 미지의 변수 $x \in X$와 데이터 $y \in Y$가 주어질 때, 다루기 힘든(intractable) 적분의 존재는 Bayesian 통계학에서 자주 나타난다.
$\bullet$ 정규화: Bayes' 이론에서 분포의 정규화 항은 고차원 공간 $X$에 대한 적분을 필요로 한다.
$p(x|y)={p(y|x) \cdot p(x) \over p(y)}={p(y|x) \cdot p(x) \over \int_{X}p(y|x')p(x')dx'}$
$p(y)$를 직접 추정하든지, 주변확률을 이용하여 추정해야 한다. 어느 경우나, 정규화 상수 $p(y)$의 추정은 쉽지 않다.
분모 $p(y)$ 값 계산을 주변 확률을 이용하려 하려면 적분 인자인 변수 $x'$가 무엇이고 어떤이 있는지 알아야 하는데 불가능하다.
$\bullet$ 주변화(Marginalization): 합이나 적분을 통해 어떤 변수를 제거할 수 있다. 이 때, 공간 $Z$가 크다면 적분이 어려워 진다.
$p(x|y)=\int_{Z} p(x,z|y)dz$
$\bullet$ 기대값(Expectation): 기대값(또는 가중 평균)이 필요할 때 적분을 사용할 수 있다. 공간이 크면 적분이 어렵다.
$E_{p(x|y)}(f(x))=\int_{X}f(x)p(x|y)dx$
(2) Statistical mechanics(주1): 상태 $s$와 Hamiltonian $E(s)$를 가진 어떤 시스템의 partition 함수 $Z$를 계산할 때 필요하다.
$Z=\sum_{S} exp[-{E(s) \over kT}]$,
여기서
$k$는 Boltzmann상수, $T$는 시스템의 Temperature를 나타낸다. 상태 $S$의 공간이 아주 크다면 합의 계산은 풀기 어려운 문제가 된다.
(3) 최적화(optimization): 최적화는 무수히 많은 해답 공간 내에서 목적함수를 최소화하는 후보를 하나 골라 내는 것이다. 그런데 이 공간은 연속이고 경계가 없다. 답을 찾기 위해 모든 공간을 뒤지는 것은 어렵다.
(4) 우도(likelihood) 모델의 선정: 두 단계로 구성된다. 먼저, 가능한 각각의 모델에 대해 ML(maximum likelihood) 추정치를 찾는다. 다음, 그 중 하나를 선정하기 위해 비용(penalization) 항(MDL, BIC, or AIC)을 사용한다. 이 방법의 문제점은 모델 집합이 아주 클 경우이다.
Monte Carlo의 예
MC의 아이디어는 어떤 고차 공간에서 정의된 확률 분포 $p(x)$로부터 상호 독립적인 샘플들의 집합 $\{ x^{(i)} \}_{i=1}^N$을 추출하는 것이다. 고차 공간은 시스템의 다양한 구성을 나타내는 공간이거나, 사후확률(posterior)이 정의된 공간, 가능한 해답의 결합 공간이 될 수 있다.
샘플들은 공간의 확률 분포를 근사적으로 표현하기위해 사용된다.
$p_{N}(x)={1 \over N} \sum_{i=1}^N {\delta_{x^{(i)}}(x)}$
여기서 $\delta_{x^{(i)}}(x)$는 Dirac-delta함수로 $x$가 $x^{(i)}$에 위치할 때 1이 되고, 그렇지 않으면 0이 된다. 따라서 $x$가 공간의 확률 분포에 따라 샘플링된다면 $p_{N}(x)$는 공간의 $pdf$에 대한 히스토그램(histogram)이 된다.
(예제) $x \in \{1,2,3 \}, N=100, pdf=Uniform$일 때 $p_{N}(2)$를 구하시오.
$p_{N}(2)={1 \over 100} \sum_{i=1}^N {\delta_{x^{(i)}}(2)}$
이고 이것은 $x^{(i)}=2$인 샘플의 비율이다. 즉, $x=2$인 확률이고, $p_{N}(2) \approx {33 \over 100}$이다.
만일 $pdf=N(2,3^2)$라면 $P_N(2)$는 ${1 \over 3}$보다 휠씬 클 것이다. $x^{(i)}$는 $pdf$에 따라 샘플링되므로 $p_{N}(x)$가 근사적으로 그 $pdf$이다.
Monte Carlo를 이용하면 다루기 힘든(intractable) 함수의 적분도 $I_{N}(f)$로 근사할 수 있다.
$x$의 $pdf$에 따른 $f(x)$의 근사 적분을 이산 합을 통해 구성해 보면
$I_{N}(f)={1 \over N} \sum_{i=1}^N { f(x^{(i)})}$
가 된다.
이 때 $x$의 공간을 이산에서 연속으로 확장하면 $N$은 $\infty$가 되고 이 식은
$I(f)=\int_{X}f(x)p(x)dx$
이다.
Metropolis Hasting 알고리즘
MH는 대표적인 MCMC 알고리즘이다. 샘플링을 임으로 하지 않고 Markov 성질을 가지는 연속 값의 체인을 통해 수행한다.
위 그림을 이해하기 위해 그림의 내부에 표시된 1->2->3 순으로 읽어본다. MH을 사용하면 최소한의 sampling으로 분포를 근사화 할 수 있다.
그래프와 코드를 보고 간략히 설명한다. 그림에서 2개의 봉우리를 가진 어떤 분포가 있다. 이 분포를 데이터 샘플링을 통해 근사하고자 한다.
분포 $q$는 proposal 확률이라 부르며, 현재 위치 $x$에서 다음 위치 $x^*$을 샘플링하는 것에 관여한다. 연속 샘플링인 마코프 체인의 이동은 acceptance 확률 $A$에 따라 $x$에서 $x^*$로 움직이는데, reject되면 $x$에 그냥 머문다.
$q$가 Gaussian proposal이라면 $q(x^*|x^{(i)})=N(x^{(i)},\sigma)$이므로 이동은 대칭형 random walk이고 $q(x^*|x^{(i)})=q(x^{(i)}|x^*)$이다.
알고리즘을 보면 먼저, $u$를 하나 샘플링하는데 이 값은 0~1사이의 uniform분포에서 얻는다. 다음 $x^*$를 $q$분포에서 하나 샘플링한다.
분포 $q$는 여러가지 형태가 될 수 있으나 여기서는 Gaussian과 같은 대칭형 분포로 $q(x^*|x^{(i)})=q(x^{(i)}|x^*)$ 관계가 성립하므로 수식의 분자, 분모에 있는 $q$는 소거되어 샘플링의 기준치는 $A=min( 1, {p({x^*}) \over {p(x^{(i)})}} )$가 된다.
$x^*$가 $q(x^*|x^{(i)})$에 따라 샘플링된다는 말은 이전 샘플링 값인 $x^{(i)}$근처에서 샘플링된다고 볼 수 있다. Gaussian 분포의 분산 크기에 따라 샘플링 점과 현재 점 사이의 거리가 결정된다. 타겟 분포 $p$는 모르지만 위치 $x$에서 확률 값 $p(x)$는 측정 가능하다고 가정한다.
어쨌던 현재의 $u$는 이렇게 계산된 $A$보다 작을 수도 있고 클 수도 있다.
만일 $x^*$에서의 확률 값 $p(x^*)$가 $p(x^{(i)})$보다 크다면 $A$는 1이되고 $u<A$이므로 $x^*$는 무조건 accept된다. 즉, 이전 샘플 $x^{(i)}$근처에서 샘플링된 $x^*$의 확률값이 $x^{(i)}$보다 크다면 이 샘플은 분포가 커지는 방향에 있고 accept되어야 좋다. 이는 확률값이 높은 곳에서 더 많은 샘플이 얻어질 수 있다는 것을 표현하게 된다.
만일 $x^*$에서의 확률 값 $p(x^*)$가 $p(x^{(i)})$보다 작다면 $A$는 1보다 작아진다. 이 때는 $u$의 값에 따라 $u<A$가 되거나 $u>A$가 될 수도 있다. $u$의 값에 따라 $x^*$가 받아질 수도 있고 아닐 수도 있다.
이러한 방식으로 아주 많은 샘플을 취하면 이 샘플들의 빈도는 확률값의 크기에 따라 얻어지게 되고 높은 확률 값에서 많은 샘플이 얻어지므로 확률분포를 잘 근사하게 된다.
여기서 Monte Carlo는 샘플링으로 어떤 값을 근사한다는 의미로 사용되었고, Markov Chain은 샘플링 시, 현재 샘플링 값은 이전 샘플링 값과 독립적이지 않고 어떤 관계를 가지면서 얻어진다는 의미로 사용된다. 따라서 MCMC 근사가 된다.
위 그림은 분포 $q$의 분산 크기가 미치는 영향을 보여 준다. 분산 값에 따라 MCMC근사의 정밀도가 달라진다.
(주1) Statistical mechanics is a branch of mathematical physics that studies, using probability theory, the average behavior of a mechanical system where the state of the system is uncertain.
References
[1] http://arongdari.tistory.com/m/post/view/id/62
[2] C. Andrieu, An Introduction to MCMC for Machine Learning, Machine Learning, 2003.
