본문 바로가기
사실 암산으로 대부분 어렵습니다. - MCMC (마코프체인 몬테카를로)

방금 전까지는 Conjugate Prior를 이용해서 암산으로 Posterior를 구할 수 있었습니다만, 실제로는 베이지안의 가장 큰 단점이 Conjugate 관계가 없는 경우에는 계산이 어렵다는 점입니다. 이게 그냥 어려운 게 아니고, 대단히 (졸라) 어렵습니다. 생각만 해도 머리가 지끈. 그래서 사용하는 것이 사후 분포를 시뮬레이션한 표본을 만들어 낸 후에 그것을 이용해서 대략~ 분포를 계산합니다. 그것을 가능하게 하는 것이 바로 MCMC라는 시뮬레이션 방법인데요, 이게 또 먼 소리냐 싶죠. 대략이라니. 어쨌든 수학을 도구로 사용하고 수식으로 명확하게 Derive가 안되면 직성이 안 풀리도록 깐깐한 사람에게는 으~~ 음? 이래도 되나 싶을 정도의 내용입니다. 수학적으로 모델링한 수식을 직접 계산해서 확률분포를 구하는 것이 아니라, 무작위 실험을 수없이 반복해 결과를 얻어내는 것이니까요.

MCMC는 Markov Chain Monte Carlo - 마코프체인 몬테카를로 -를 줄여서 부르는 용어입니다. 여기에서 마코프체인이라고 함은 바로 이전의 상태에 의존하여 사건이 일어나는 것을 의미하고요, 몬테카를로는 도박의 도시이름입니다, 도박처럼 무작위 샘플링을 한다는 의미로 통용됩니다. 두개를 합치고 보면 새로운 샘플은 이전의 상태에 의존해서 샘플링을 하면서도, 무작위로 한다는 뜻인데요. 이건 또 뭔가 싶습니다. 무작위이긴 한데 이전 상태에 관계되는 뭐 그런 건데 말이죠. 

일단 마코프체인이라는 게 이전의 상태가 현재의 상태에 영향을 미친다는 거라니까, 대충의 정의들을 살펴보면, "마지막의 샘플이 다음 샘플에 영향을 미친다" → "마지막의 샘플이 다음 샘플을 추천해 준다" 라는 표현이 있습니다. 뭔가 멋진데? 싶긴 하면서도 어떻게 추천해 준다는 것인가?라는 문제에 봉착하게 됩니다. 

이걸 조금 더 이해하기 쉽도록 절차적인 표현으로 바꿔쓰면 다음과 같습니다. 

MCMC는 
⓵ 랜덤 샘플을 고름 
⓶ 고른 샘플을 기반으로 다음의 샘플을 추천  (요게 마코프)
⓷ 추천된 다음 샘플이 특정한 조건을 만족하면 채택/ 그렇지 않으면 버림
⓸ 헤헷 ⓶로 돌아가서 무쟈게 반복  (요게 몬테카를로)

가만히 보면 ⓷에서 뭔가 채택하거나 버리거나 하니까, 이 부분에서 여러 가지 특화된 알고리즘이 적용될 수 있겠습니다. 요 부분에서 유명하면서 간단하게 사용되는 Metropolis 알고리즘을 아주 간략하게 소개한 후 MCMC를 통해 어떤 식으로 Posterior를 간편하게 알아내는지를 설명할 텐데, 이때! 알고리즘에 매몰되어서 원래 뭘 하려 했는지 잊으면 안 됩니다. 알고리즘은 알고리즘 일 뿐 오히려 잘 이용하겠다는 생각만으로도 충분합니다. "MCMC는 계산하기 어려운 사후 분포를 시뮬레이션을 통해 근사적으로 알아낼 수 있도록 문제를 해결해 준다" 정도로만 알고 있어도 충분하다고 생각합니다. 

일단 Metropolis 알고리즘은 말이죠, 별건 아니고, 다음 샘플을 추천할 때, 현재 샘플을 중심으로 제안분포라는 것을 설정하고, 그 분포에서 랜덤으로 다음 샘플을 샘플링 합니다. 이 샘플링한 것이 특정한 조건을 만족하면 이 샘플을 다음 샘플로 채택, 그렇지 않으면 버림이 되는데요, 그러면 어떤 특정한 조건을 만족하면 채택하게 되는지가 관건이겠죠.

자, 그럼 이 시점에 궁금한 점은 제안분포라는 것은 무엇이고, 특정조건이 무엇인가? 라는 점이 당연히 궁금해야 하겠습니다. 

제안분포는 다음 샘플을 고를 때, 현재 샘플을 중간에 놓고 Uniform 또는 Gaussian분포를 만듭니다. 그 분포에서 다음 샘플을 채취하는 것인데요, 그림으로 보면 간단합니다.

 

제안분포를 가우시안으로 둔다면, 현재의 샘플을 평균으로 한 가우시안 분포를 만들고 이 분포에서 랜덤 샘플링을 하여 next 샘플을 채취합니다. 제안분포란 다음 샘플을 채취하기 위한 분포입니다. 

그런 다음에, 채취한 다음번 샘플이 괜찮은가?를 판단해야 하겠는데, 샘플이 괜찮은가를 보는 것이 Likelihood 곱하기 Prior의 곱의 비를 확인하는 것입니다. 한마디로 Posterior의 비를 비교해서 사후확률이 될만한 녀석들만 샘플로 인정하는 거죠. 간단하죠?

자, 그럼 어떤 식으로 비교하느냐? 하면 현재와 다음 샘플에 대해서 Posterior를 구한 후에 비율을 비교해서 다음 샘플의 Posterior가 더 큰 경우에, 그리고 Uniform Distribution에서 Random 1 Value보다 이 비율이 큰 경우에 추천된 다음 샘플을 채택합니다. 

$$ \alpha = \cfrac{P(\theta_{next}| Sample)}{P(\theta_{curr}|Sample)} = \cfrac{P(\theta_{next}| x)}{P(\theta_{curr}|x)} > 1 $$

이런 셈 인데요, 딱 느낌 오죠? 새로운 샘플의 Posterior가 더 큰 경우에 채택하고요, 그런데 이 Posterior를 모르는 게 문제잖아요? 그러니까, 이걸 베이지안으로 보게 되면 Prior와 Likelihood의 곱으로 나타낼 수 있겠죠? 게다가 Sample에 대한 Evidence (분모는 Sample에 대한 x marginal 확률)이 같으니까 약분해서 간단하게 판단할 수 있겠습니다. 후후. 그르니깐~

$$ \alpha = \cfrac{P(x | \theta_{next})P(\theta_{new})/P(x)}{P(x|\theta_{curr})P(\theta_{curr})/P(x)} = \cfrac{P(x | \theta_{next})P(\theta_{new})}{P(x|\theta_{curr})P(\theta_{curr})} > 1 $$

가 되겠죠. 이렇게 구한 새로운 Sample을 채택하게 되면 이게 다시 현재의 샘플이 되고, 이 샘플을 기준으로 제안분포를 만들어서 거기에서 새로운 샘플을 또 구해서 같은 짓을 하여 아주 마~~아니~~ 하게 되면 이게 결국에는 Posterior 샘플을 만들어내는 시뮬레이션 결과가 됩니다. 생각보다 간단합니다. 

 

자, 이렇게 하면 계속 더 Posterior의 가장 밀도가 큰 쪽으로 이동해서 그곳에서만 샘플을 만들어 내겠죠? 그렇긴 한데, 1보다 큰 것들만 채택하면 좋겠는데, α비율이 1보다 작은 경우도 꽤나 나오거든요. 이것을 아주 버리는 것이 아니라, 균등분포(Uniform)에서 1개 랜덤 샘플링을 해서 그 값보다 크면 그 샘플도 무작정 버리지 말고 다음 샘플로 인정하자~ 뭐 그런 것입니다. 그러니까  0~1사이의 랜덤 샘플보다는 크면 α비율이 1보다 작은 경우일지라도 뭐 적당히 새로운 샘플의 Posterior가 그런대로 쓸만하다~ 고 인정하자 하는 것입니다. 그렇게 되면 결국, 그리고 어차피, 밀도가 큰 곳에서의 샘플이 많이 채취될 것이고 밀도가 낮은 곳에서의 샘플이 덜 채취될 것이라는 뭐 그런 것입니다. 

Metroplis 알고리즘을 베이지안 Posterior 추정에 사용하는 예제로 곧바로 돌입해서 깨방정을 떨어보시죠. 이게 말로 하니까 굉장히 어려워 보이는 데, 막상 구현한 코드를 보면 굉장히 쉽다고 느껴질 것인데요. 샘플링 시뮬레이션이 잘 되었는지 알아야 하니까, 우리가 정답까지 잘 알고 있는 Beta Prior Conjugate관계를 이용해서 한번 보시죠.

자. 예시 시작. 어떤 Posterior 분포를 구하려고 하는데요, Prior가 Beta(1,1)이고, Likelihood는 추가로 100회 시도에 70성공, 30 실패인 경우의 Posterior를 MCMC로 시뮬레이션해 보려고 합니다. 이게 우리 이미 암산으로 풀어봤잖아요? 사후 분포는 Beta(71,31)이 되어야 하겠군요. Metropolis 알고리즘을 이용한 MCMC 시뮬레이션결과가 이 분포를 따르는지 한번 확인해 보겠습니다. 

import numpy as np
from scipy.stats import beta
from scipy.stats import binom

def prior(a, b, p): # prior를 베타로 선언하고
    return beta.pdf(p, a, b)

def likelihood(n, x, p): # Likelihood를 Binomial로 선언하고 
    return binom.pmf(x, n, p)

 

이런 식으로 Prior(베타), Likelihood(Binomial)를 선언하고요. 여기에 각각의 분포에 대한 Parameter 설정을 한 후에 MCMC를 해 보도록 하겠습니다. 여기에서 Priror와 Likelihood가 argument로 받는 p가 θ(theta)가 되겠습니다. 길지 않으니까 너무 부담 갖지 말아 주세요.

# Prior Parameter ← Beta(1,1)
a = 1
b = 1

# Likelihood Parameter ← Binomial(100,0)
n = 100
x = 70  # 70 성공/ 30 실패 → 71/31이 됨 , B(71,31)

# 시뮬레이셔 샘플 총 겟수
n_samples = 2000
theta_samples = []  # 샘플을 모아둘 리스트

theta_curr = np.random.uniform(0, 1, 1)[0]  # 현재 샘플의 초기값
for _ in range(n_samples):
    
    # 제안 분포 Gaussian(theta,0.01)에서 다음 샘플 1개 추출
    theta_next = np.random.normal(theta_curr, 0.01, 1)[0] 
    
    # 목표 분포 판단을 위한 Posterior의 비 
    alpha = min(1, 
                ((likelihood(n, x, theta_next) * prior(a, b, theta_next)) / 
                (likelihood(n, x, theta_curr) * prior(a, b, theta_curr))))

    # Uniform(1,1)에서 threshold u_th 추출
    u_th = np.random.uniform(0, 1, 1)[0]
    
    # alpha와 u_th를 비교해 alpha가 u보다 크면 후보 값을 채택하여 샘플 list에 추가. 
    if u_th < alpha:
        theta_samples.append(theta_next)
        # 추출된 샘플을 현재 샘플로 갱신
        theta_curr = theta_next

 

이렇게 시뮬레이션하면 어떻게 될지 대~충 상상이 가능할 텐데요. 후후. 시뮬레이션 결과를 그려보면, 점선은 우리가 암산으로 구한 베타사후분포, 막대그래프는  베타사후분포에 대한 시뮬레이션 결과입니다. 주의할 점은 시뮬레이션하여 얻은 샘플들에 대해서는  hist(density=True, bins=30)정도로 그려주면 예쁘게 나옵니다. 

 

어때요 결과를 보니까 거~의 비슷하지요? 그러니까 너무 어렵게 생각 말고 계산이 어려우면 MCMC를 이용해서 사후분포를 시뮬레이션해서 대~략 구할 수 있다는 점만 알고 있으면 딱 좋은 이해의 정도입니다.

사실 마지막에 Uniform으로 샘플의 적절성을 판단하긴 했는데, 원래 제안분포가 가우시안이었기 때문에 가우시안으로 판단해도 됩니다.

    # Uniform(1,1)에서 threshold u 추출
    u = np.random.uniform(0, 1, 1)[0]

이 부분을 

    # Gaussian(0, 0.01)에서 threshold u 추출
    u = np.random.normal(0, 0.01, 1)[0]

이런식으로 해도 적절성 평가를 할 수 있습니다.

그런데, 가만히 보면 컴퓨터를 이용해서 뭔가 해답을 찾을 때에는 일단 한번 해보고 그것을 평가한 후에 다음번에는 어떻게 해야 할지에 대한 알고리즘에 따라 행동을 하게 하는 것이 기본 컨셉입니다. 이 개념을 알고 있으면 기계학습에 대해서는 반은 끝낸 것이나 다름없습니다. 

그러니까 이 경우에는 평가를 베이즈룰의 분자, 즉, likelihood × prior의 크기로 평가를 하여 더 괜찮은 샘플을 구한다고 보면 틀림없습니다. 

친절한 데이터 사이언스 강좌 글 전체 목차 (정주행 링크) -



댓글





친절한 데이터 사이언스 강좌 글 전체 목차 (정주행 링크) -