for textmining

Gibbs Sampling

|

이번 글에서는 깁스 샘플링(Gibbs Sampling)에 대해 간단히 살펴보도록 하겠습니다. 이번 글 역시 고려대 강필성 교수님 강의와 위키피디아, ‘밑바닥부터 시작하는 데이터과학(조엘 그루스 지음, 인사이트 펴냄)’, 그리고 이곳을 정리했음을 먼저 밝힙니다.

개요

깁스 샘플링은 두 개 이상의 확률변수의 결합확률분포로부터 일련의 표본을 생성하는 확률적 알고리즘입니다. 결합확률분포나 그에 관련된 확률 계산을 근사하기 위해 널리 사용되고 있습니다. 깁스 샘플링은 마코프 연쇄 몬테카를로 방법(Markov Chain Monte Carlo;MCMC) 알고리즘의 한 예라고 합니다. 깁스 샘플링을 이해해보기 위해 몬테카를로 방법과 MCMC를 먼저 알아보겠습니다.

몬테카를로 방법

몬테카를로 방법(Monte Carlo Method)이란 랜덤 표본을 뽑아 함수의 값을 확률적으로 계산하는 알고리즘을 가리킵니다. 수학이나 물리학 등에 자주 사용되며 계산하려는 값이 닫힌 형식으로 표현되지 않거나 복잡한 경우에 그 값을 근사적으로 계산하려고 할 때 쓰입니다. 몬테카를로 방법의 대표적인 사례가 바로 원주율 계산입니다. 아래 그림(출처)을 보겠습니다.

원주율 계산을 위한 몬테카를로 방법은 다음과 같습니다.

(1) [0,1] x [0,1]에서 점 $(x,y)$를 무작위로 뽑는다.

(2) 이 점이 중심이 (0,0)이고 반지름이 1인 원에 속하는지 계산한다. (즉 $x^2+y^2≤1$을 만족하면 빨간색 점, 그렇지 않으면 파란색 점으로 분류)

(3) 위의 두 과정을 충분히 반복한다.

위 그림에서 4분 원의 넓이는 $π/4$입니다. 전체 점 개수를 빨간색 점 개수로 나눈 비율이 이 값에 근사합니다. 이를 통해 우리는 $π$값이 얼마인지 대략적으로 추정할 수 있습니다.

마코프 연쇄

마코프 연쇄(Markov Chain)마코프 가정(Markov assumption)을 따르는 이산 시간 확률 과정을 가리킵니다. 마코프 가정은 러시아 수학자 마코프가 1913년경에 러시아어 문헌에 나오는 글자들의 순서에 관한 모델을 구축하기 위해 제안된 개념입니다. 특정 시점의 상태 확률은 단지 그 이전 상태에만 의존한다는 것이 핵심입니다. 즉 한 상태에서 다른 상태로의 전이(transition)는 그동안 상태 전이에 대한 긴 이력(history)을 필요로 하지 않고 바로 직전 상태에서의 전이로 추정할 수 있다는 이야기입니다. 마코프가정은 아래와 같이 도식화됩니다.

그런데 특정 조건을 만족한 상태에서 마코프 연쇄를 반복하다 보면 현재 상태의 확률이 직전 상태의 확률과 같아지게, 즉 수렴하게 됩니다. 이렇게 평형 상태에 도달한 확률 분포를 정적분포(Stationary Distribution)라고 합니다.

마코프 연쇄 몬테카를로 방법

MCMC란 마코프 연쇄에 기반한 확률 분포로부터 표본을 추출하는 몬테카를로 방법입니다. 다음과 같습니다.

MCMC 알고리즘은 우리가 샘플을 얻으려고 하는 목표분포를 Stationary Distribution으로 가지는 마코프 체인을 만든다. 이 체인의 시뮬레이션을 가동하고 초기값에 영향을 받는 burn-in period를 지나고 나면 목표분포를 따르는 샘플이 만들어진다.

깁스 샘플링

깁스 샘플링은 MCMC의 일종인데요. 몬테카를로와 MCMC와의 차이점은 이렇습니다. 몬테카를로 방법은 모든 샘플이 독립(independent)이고 생성될(뽑힐) 확률 또한 랜덤입니다. 반면 마코프 연쇄에 기반한 MCMC는 다음번 생성될(뽑힐) 샘플은 현재 샘플의 영향을 받습니다. 깁스 샘플링은 다음번 생성될 표본은 현재 샘플에 영향을 받는다는 점에서는 MCMC와 같지만, 나머지 변수는 그대로 두고 한 변수에만 변화를 준다는 점이 다릅니다.

예를 들어보겠습니다. 3개의 확률변수의 결합확률분포 $p(x_1,x_2,x_3)$로부터 1개의 표본을 얻으려고 할 때 깁스 샘플링의 절차는 다음과 같습니다.

(1) 임의의 표본 $X^0=(x_1^0,x_2^0,x_3^0)$을 선택한다.

(2) 모든 변수에 대해 변수 하나만을 변경하여 새로운 표본 $X^1$을 뽑는다.

실제 사용시에는 처음 수집한 표본 $X^0$은 버리고 $X^1$만 쓰게 됩니다. (2)를 자세하게 쓰면 다음과 같습니다.

(ㄱ) 현재 주어진 표본 $X^0$의 두번째, 세번째 변수 $x_2^0$, $x_3^0$를 고정시킨다.

(ㄴ) 첫번째 기존 변수 $x_1^0$를 대체할 새로운 값 $x_1^1$을 다음과 같은 확률로 뽑는다. p($x_1^1$|$x_2^0$,$x_3^0$)

(ㄷ) 첫번째 변수 $x_1^1$, 세번째 변수 $x_3^0$ 를 고정시킨다.

(ㄹ) 두번째 기존 변수 $x_2^0$를 대체할 새로운 값 $x_2^1$을 다음과 같은 확률로 새로 뽑는다. p($x_2^1$|$x_1^1$,$x_3^0$)

(ㅁ) 첫번째 변수 $x_1^1$, 두번째 변수 $x_2^1$를 고정시킨다.

(ㅅ) 세번째 기존 변수 $x_3^0$를 대체할 새로운 값 $x_3^1$을 다음과 같은 확률로 새로 뽑는다. p($x_3^1$|$x_1^1$,$x_2^1$)

(ㅇ) 최종적으로 구한 $X^1$은 다음과 같다. $X^1=(x_1^1,x_2^1,x_3^1)$

(2)에서 새로운 값을 뽑는 데 쓰인 조건부 확률은 결합확률분포 $p(x_1,x_2,x_3)$에 비례한다고 합니다. 표본의 앞부분은 초기 상태 $X^0$에 크게 의존하지만 충분히 많이 뽑고 난 뒤에는 초기 상태에 관계 없이 $p$에 기반한 표본을 수집할 수 있다고 합니다.

2차원의 가우시안 분포를 만들어 내기 위해 깁스 샘플링을 적용한 예시 그림은 다음과 같습니다.

깁스 샘플링의 변형들

변수가 $(a,b,c)$ 세 개인 데이터에 대해 깁스 샘플링을 수행한다면 $b,c$를 고정시킨 채로 $a$를, $a,c$를 고정시킨 채로 $b$를, $a,b$를 고정시킨 채로 $c$를 차례대로 뽑아야 합니다.

Block Gibbs sampling 기법은 그룹으로 묶어 뽑는 기법입니다. $c$를 고정시킨 채로 $a,c$를 뽑고 $a,b$를 고정시킨 채로 $c$를 뽑는 방식입니다.

Collapsed Gibbs sampling 기법은 불필요한 일부 변수를 샘플링에서 생략하는 기법입니다. $b$가 그런 변수라 가정하면 $c$를 고정시킨 상태에서 $a$를 뽑고, $a$를 고정시킨 상태에서 $c$를 뽑습니다.

파이썬 구현

이해를 돕기 위해 ‘밑바닥부터 시작하는 데이터 과학’에 나온 간단한 예시를 들겠습니다. 깁스 샘플링은 임의의 $x$ 또는 $y$값에서 출발해서 $x$에 대한 $y$의 조건부확률과 $y$에 대한 $x$의 조건부확률 사이를 오가며 반복적으로 값을 선택하는 방법입니다. 이 과정을 여러번 반복하여 얻은 $x$와 $y$는 결합확률분포에서 얻은 샘플과 유사한 값이 됩니다.

우선 문제를 정의해 보겠습니다. 주사위 두 개를 던진다고 해 봅시다. 첫번째 주사위의 눈을 $x$, 두 주사위의 눈을 합한 값을 $y$라고 하면 $x$와 $y$의 결합확률분포 함수는 다음과 같습니다.

import random

def roll_a_die():
    # 주사위 눈은 1~6
    # 각 눈이 선택될 확률은 동일(uniform)
    return random.choice(range(1,7))

def direct_sample():
    d1 = roll_a_die()
    d2 = roll_a_die()
    return d1, d1+d2

$x$에 대한 $y$의 조건부확률과 $y$에 대한 $x$의 조건부확률 함수는 다음과 같습니다.

def random_y_given_x(x):
    # x값을 알고 있다는 전제 하에
    # y값이 선택될 확률
    # y는 x+1, x+2, x+3
    # x+4, x+5, x+6 가운데 하나
    return x + roll_a_die()

def random_x_given_y(y):
    # y값을 알고 있다는 전제 하에
    # x값이 선택될 확률
    # 첫째 둘째 주사위 값의 합이 7이거나
    # 7보다 작다면
    if y <= 7:
        # 첫번째 주사위의 눈은 1~6
        # 각 눈이 선택될 확률은 동일
        return random.randrange(1, y)
    # 만약 총합이 7보다 크다면
    else:
        # 첫번째 주사위의 눈은
        # y-6, y-5,..., 6
        # 각 눈이 선택될 확률은 동일
        return random.randrange(y-6, 7)

깁스 샘플링 함수는 다음과 같습니다.

def gibbs_sample(num_iters=100):
    # 초기값이 무엇이든 상관없음
    x, y = 1, 2
    for _ in range(num_iters):
        x = random_x_given_y(y)
        y = random_y_given_x(x)
    return x, y

깁스 샘플 수를 늘려서 결합확률분포 direct_sample로부터 뽑은 결과와 비교하면 유사한 결과가 나오는걸 확인할 수 있습니다. 다시 말해 결합확률분포를 모를 때, 이미 알고 있는 일부 조건부 확률분포에 깁스 샘플링을 적용하여 해당 결합확률분포의 표본을 얻어낼 수 있다는 것입니다.

Comments