상세 컨텐츠

본문 제목

그림으로 깁스샘플링 이해하기

그냥 공부

by 적분 ∫2tdt=t²+c 2020. 1. 3. 00:13

본문

깁스 샘플링에 대해서는 토픽 모델링 공부를 하면서 귀에 피가 나도록 들었었는데요, 사실 그 실체를 정확하게 파악하는데에는 꽤 오랜 시간이 걸렸었습니다. 막연하게만 이해하고 있던 깁스 샘플링을 정확하게 깨닫게 된 건 어떤 2장 짜리 논문에서 본 그림 덕분이었는데, 이번 포스팅에서는 그걸 공유해보고자 합니다.

참고한 논문은 다음과 같습니다.


Breslaw, J. A. (1994). Random sampling from a truncated multivariate normal distribution. Applied Mathematics Letters7(1), 1-6.


다변수 정규 분포

다변수 정규 분포(Multivariate Normal Distribution)에 대해서는 상관 토픽 모델 포스팅에서 살짝 설명했었는데요, 우리가 흔히 보던 종 모양의 정규분포 곡선을 2차원 이상으로 확장한 것이라고 보면 됩니다. 변수가 2개 이상이 되면, 단일 변수의 분산뿐만 아니라 변수들 간의 분산(공분산)이라는 개념도 생겨나게 됩니다. 이는 두 변수가 얼마나 같이 움직이는지를 보여주는 척도가 됩니다. 아래의 그래프는 2차원 정규분포의 밀도를 표현한 것인데, X와 Y의 평균은 둘 다 0, X와 Y의 분산도 둘다 1이고, X와 Y의 공분산은 0.8인 경우입니다. 괜시리 어렵게 수식으로 쓰면 다음과 같겠네요.

위 그림은 [-2, 2] 범위를 0.1 단위로 격자를 나눈 뒤, 앞서 말한 2차원 정규분포에서 점을 1만개를 뽑아 어느 격자에 속하는지를 센 뒤 그 비율에 따라 색깔을 달리하여 칠한 것입니다. 당연히 중심인 0,0 부분이 제일 진하고 주변으로 갈수록 옅어지는 모양이 됩니다. 사실 [-2, 2] 범위 바깥에서 뽑힌 점들도 있지만, 그 수가 너무 적기 때문에 그림에서는 표현하지 않았습니다.

우리가 이 포스팅에서 관심 가지는 것은 '저런 분포가 주어졌을때 어떻게 거기서 표본을 뽑아낼 것인가'입니다. 다행히도 정규분포에서 임의의 표본을 뽑아내는 방법은 박스 뮬러 변환(Box-Muller Transform)이라는 효과적인 기법이 개발되어 있어서, 이를 이용해 누구나 정규분포에서 손쉽게 표본들을 샘플링할 수 있습니다. 다변수 정규분포 역시 본질적으로 정규분포를 여러 차원에 대해 반복하는 것으로 샘플링할 수 있습니다.

이제 한 발 더 나아가봅시다. 위의 빨간 네모는 -2<X<-1, 0<Y<1인 영역을 표현하는데요, 정규분포 중 저 안에 들어가는 녀석들에 대해서만 표본을 뽑는 문제를 생각해볼 수 있겠습니다.


절단 정규 분포와 기각 샘플링

절단 정규 분포(Truncated Normal Distribution)는 정규 분포에서 특정영역만 잘라내어 만든 분포를 가리킵니다. 다르게 말하면 원래는 X의 범위가 무한했던 정규분포에서 X가 가질 수 있는 범위를 제한시킨 분포라고 할 수 있겠네요. 종종 복잡한 확률 과정을 계산하다보면 정규 분포의 정의역을 여러 가지로 제한하게 되는데 이때 사용하게 되는 분포입니다. 이번 예제에서는 위의 2변수 정규분포에서 정의역을 -2<X<-1, 0<Y<1로 제한(빨간색 박스)한 절단 정규 분포를 생각해보겠습니다. 

위의 빨간 박스부분만 떼어서 다시 그림을 그리면 아래와 같아집니다. 원래 -2<X<-1, 0<Y<1 범위 내의 표본 수는 매우 적었지만, 이제 전체 구간을 줄여버렸으므로 상대적으로 (-1, 0) 부분의 색이 진해집니다. 이제 이 범위 내에 속하는 표본을 뽑는 작업을 생각해봅시다. 이는 생각보다 어려운 작업인데요, 왜냐하면 전체 2차원 정규분포에서 표본을 뽑으면 이 구간에 들어가지 않을 확률이 들어갈 확률보다 훨씬 높기 때문입니다. (위의 2차원 정규분포의 중심은 (0,0)에 있었다는 걸 기억해보세요. 대부분의 표본은 그 근처에서 뽑힐 겁니다.)

그래도 끈기를 가지고 수많은 표본을 뽑다보면 우연히 우리가 원하는 사각형 내에 들어오는 샘플이 있을테니깐 걔네를 열심히 모아보도록 합시다. 이 방법을 기각 샘플링(Rejection Sampling)이라고 합니다. 이 방법의 장점은 '별 다른 고민없이 기존에 사용하던 샘플링 알고리즘을 그대로 사용할 수 있다'는 것이고, 또한 '의도한 분포대로 표본을 뽑을 수가 있다는 것'입니다. 다만 단점은 '표본의 기각률이 높을 경우 실용적으로 써먹기엔 틀렸다'는 것입니다. 

지금 예시로 든 이 절단 정규분포에서 1만개의 표본을 뽑기 위해서는 원 2변수 정규분포에서 약 188만개의 표본을 뽑아야 합니다. 그래야 그중 1만개 정도가 우리가 원하는 영역에 들어온다는 것이죠. 만약에 2차원이 아니라 3차원 절단 정규분포라면? 기각률은 훨씬 더 높아질 겁니다. 보통 우리가 다루는 분포는 최소 몇십 차원에서 몇 백 차원을 넘나든다는 것을 떠올려보면, 기각 샘플링은 절단 정규분포를 표집하는데에는 적절하지 않다는 결론을 얻게 됩니다.


-보너스-

기각 샘플링을 파이썬 코드로 작성하면 다음과 같습니다.



각개격파!

그래서 다른 방법을 찾아야 합니다. 우리가 뽑을 표본은 2차원 상의 점이므로, 이를 하나씩 분리해서 처리하면 조금 더 쉽지 않을까요? 사실 앞에서 밝히지 않은 게 있는데, 다변수 정규 분포를 뽑는 과정이 정규 분포를 하나씩 뽑아서 모으는 과정이라는 겁니다. 이걸 그대로 절단 정규 분포에 적용하면 되지 않을까요?

두 변수의 공분산을 ρ이라 할때 다음 2변수 정규분포는 다음과 같이 정규분포 두 개로 분해할 수 있습니다.

일때

라고 하면


가 된다.

여기서 Z1, Z2는 표준 정규분포에서 뽑은 표본입니다. X와 Y는 Z1과 Z2의 조합으로 표현할 수 있습니다. 여기서는 먼저 Z1을 X로 잡았고, 이 Z1의 값에 따라 Y값이 영향을 받게 됩니다. 혹은 정반대로

로 정의해도 무방합니다.


이제 여기서 절단 조건을 추가해보도록 합시다. a<X<b여야 하고, c<Y<d여야 한다고 해봅시다. 그럼 먼저

,

따라서 Z1은 truncN(a, b)에서 표집하면 되고, Z2는 Z1를 표집한 다음 그 값을 이용해 truncN(~~, ~~)에서 표집하면 되겠습니다. (수식 적기 귀찮아서 ~~로 대체했습니다. 다들 생략된게 뭔지 아시죠?) 여기서 truncN은 1변수 절단 정규 분포를 가리킵니다. 다행히도 1변수 절단 정규분포의 경우 샘플링 기법이 발명되었기 때문에 손쉽게 샘플링을 실시할 수 있습니다. 만약 이를 사용할 수 없는 경우라면 정규분포에서 기각 샘플링을 해도 좋습니다. (1변수 정규분포의 경우 다변수에서처럼 기각률이 높지 않기 때문에 사용해봄직합니다.)

파이썬 코드로 구현해보면 다음과 같습니다.


자, 이렇게해서 샘플링을 하면 어떤 분포가 나올까요? 짜잔~


모양이 좀 다르죠? 원래 우리가 원했던 모양은 오른쪽 아래만 새파랗고 나머지는 하얘야하는데, 이 방법으로 추출한 경우 Y = 0쪽에 뭔가 전반적으로 분포가 쏠려있습니다. 왜 이렇게 잘못된 분포가 나오게 되었을까요? 이는 우리가 X를 먼저뽑고 그에 따라 조건부 확률로 Y를 뽑았기 때문입니다. 반대로 Y를 먼저 뽑고 X를 Y에 대한 조건부로 뽑았다면 X = -1 쪽에 쏠린 분포가 나오게 됩니다.

어디가 잘못되는지 감이 안온다면 문제를 조금 단순화시켜봅시다. X, Y ∈ {a, b, c}인 확률 변수이고,  a, b, c가 뽑힐 확률은 각각 1:2:3의 비율을 이룬다고 가정해봅시다. 그런데 X와 Y를 뽑는 사건은 독립이 아니라서 X와 Y를 동시에 뽑는 경우는 단순히 X를 뽑을 확률과 Y를 뽑을 확률을 곱해서 구할 수 없고 다음과 같이 계산된다고 가정해봅시다.

X=a

X=b

X=c

Y=a0.1
Y=b0.10.2
Y=c0.10.20.3


자 만약 이 분포를 근사하기 위해 X를 먼저 뽑고, 그 다음 거기서 Y를 뽑는다면 어떻게 될까요? 먼저 X가 뽑힐 확률은 다음과 같이 됩니다.

X=a

X=b

X=c

0.1670.3330.5

보면 벌써부터 틀렸습니다. X와 Y의 결합분포에서는 X=a가 되는 경우는 0.1밖에 없으나 Y를 무시하고 X를 뽑으면 X=a가 뽑힐 확률이 0.167로 0.1보다 훨씬 큽니다. X=b인 경우도 마찬가지이죠. X=a와 X=b가 원 분포보다 더 많이 뽑히니 이 상태에서 Y가 뽑힐 확률을 고려해 X와 Y의 결합분포를 계산해보았자 틀린 결과가 나오게 됩니다.

이게 바로 위 그림에서 발생한 현상과 동일한 현상입니다. 원래는 X=-2 근방에서 뽑힐 확률이 더 적어야하지만, 우리는 Y를 무시하고 X에 대해서 우선 샘플링한 뒤 Y를 샘플링하게 되므로, X=-2 근방에서 뽑힐 확률이 증가하여 아래쪽이 두터운, 잘못된 모양의 분포를 산출하게 된 것이지요.

따라서 각개격파 방법은 다변수 절단 정규분포를 샘플링하는데에는 쓸수 없다는 사실을 깨닫게 됩니다. 그렇다면 우리가 원하는 분포를 추출하기 위해서는 기각 샘플링 삽질을 할 수 밖에 없는걸까요?


깁스 샘플링!!!

이럴 때 사용할 수 있는게 깁스 샘플링(Gibbs Sampling, 기브스 표집)입니다. 깁스 샘플링의 아이디어를 한 줄로 줄여 말하면, X와 Y가 서로 종속적인 관계에 있으니깐, X값을 정한상태에서 Y를 뽑아서 Y값을 구하고, 다시 거꾸로 이 Y값에서 X를 뽑아서 X값을 정하자는 것이죠. 수식으로 멋지게 쓰면 다음과 같겠습니다.

풀어보자면, 첫번째 X는 그냥 독립분포인 P(X)에서 뽑고, Y는 뽑은 X값을 조건부로하는 Y분포에서 뽑습니다. 그리고 각자 이전의 X 혹은 Y값을 바탕으로 조건부확률을 계산하여 다음 X혹은 Y값을 뽑는 것이죠. 포인트는 현재 샘플링하는 변수 이외의 나머지 변수들이 조건으로 고정되어 있다고 가정한다는 것입니다. (사실 초기값인 X0, Y0은 꼭 P(X)나 P(Y|X=X0)에서 뽑지 않고 아무 랜덤값을 주어도 됩니다.) 이렇게 새로운 X와 Y를 뽑는 과정을 계속 반복하다보면 여기서 나오는 X_i, Y_i값들이 원래 우리가 구하고자했던 X, Y의 결합 분포에서 추출한 표본이 된다는 겁니다. 신기하죠?

깁스샘플링은 변수가 더 많을 때에도 마찬가지로 잘 작동하기 때문에 변수가 굉장히 많고 그 조건들이 까다로워 전체의 결합 확률분포를 계산해 거기서 샘플을 추출하기 어려울 때 아주 유용하게 사용할 수 있습니다. 이번 포스팅에서 예시로 들고 있는 다변수 절단 정규 분포가 딱 그에 맞는 상황입니다. 조건부 확률에서는 쉽게 샘플링할 수 있지만, 결합 확률에서는 샘플링이 어려운 경우 말입니다.


이 과정을 파이썬 코드로 정리하면 다음과 같습니다.

짜잔! 그러면 이렇게 처음에 우리가 그렸던 그림과 거의 동일한 분포를 얻을 수 있습니다. 이 그림은 깁스 샘플링을 1회 돌려서 나온 X, Y값을 사용한 경우이구요, 더 많은 횟수를 돌린 경우에도 위와 거의 동일한 모양의 분포가 나옵니다. 보통 깁스 샘플링 사용시 초기 몇 십 회 이내의 표본은 편향되어 있을 가능성이 높아서 사용하지 않습니다만 (이렇게 앞부분의 샘플들을 생성하는 과정을 Burn-in iteration이라고 합니다. 본격적으로 샘플을 뽑기에 앞서 불을 때는 과정이라고 생각하면 되겠죠), 예시로 든 2변수 절단 분포의 경우 차원수가 낮기 때문에 1회만으로도 원 분포와 거의 유사한 분포가 나왔습니다.


그리고 위 코드에서는 매번 깁스 샘플링 결과를 얻기 위해 초기값으로 X, Y를 설정하고 burn-in을 거쳐 최종 샘플을 뽑는 과정을 반복하고 있는데, 굳이 그럴 필요 없이 매 iteration 마다 생성되는 표본들을 연속적으로 사용해도 동일한 분포를 얻을 수 있습니다. 단, 이 경우 이전의 표본과 바로 다음의 표본은 서로 직접적으로 종속관계에 있으므로 대게 1~2개씩 건너뛰면서 사용하는 경우가 많습니다. (10회째, 12회째, 14회째 ... 결과를 사용하는 것처럼)


수치로 비교해보기

그림으로만 보면 우리가 근사한 분포가 얼마나 원 분포와 얼마나 가까운지 확인하기 어렵습니다. 따라서 두 확률분포의 거리를 측정하는데 쓰이는 KL Divergence 척도를 이용해 유사도를 확인해보았습니다.

각개격파깁스 샘플링 1회3회5회
KL Div0.21490.01390.01300.0105

KL Divergence는 두 분포가 완전히 같은 경우 0, 다를수록 큰 양수의 값을 가집니다. 처음 시도했던 각개격파 방법의 경우 0.2149로 매우 큰 오차를 보였지만, 깁스 샘플링의 경우 0.01 정도로 거의 원 분포와 동일한 근사 분포를 보이고 있음을 확인할 수 있습니다. 특히 그 차이는 샘플링 횟수가 늘어날수록 줄어들고 있구요.

정리해보자면, 변수가 많고 복잡한 확률분포에서 표본을 표집하는 것은 때때로 어려운데, 이를 단순한 확률 분포 상의 표집의 반복으로 바꿔주는 장치가 깁스 샘플링이라 할 수 있겠습니다. 본 포스팅에서는 1변수 절단 정규분포를 이용해서 다변수 절단 정규분포를 표집하는 것을 실험해보았구요. 현업에서 사용되는 다양한 베이지안 모델은 이보다 더 복잡한 경우가 많으니, 당연히 이 곳에서도 매우 유용하게 쓰이고 있습니다. 이번 기회에 그 개념을 확실히 잡아서 다음번에 깁스 샘플링을 만나도 당황하지 않으시면 좋겠네요!

관련글 더보기

댓글 영역

  • 프로필 사진
    2020.03.17 20:06
    중간에 설명을 보면, A,B,C가 1:2:3의 비율로 뽑힌다는 설명이 있는데 조금 헷갈려서 질문을 드려요. 위와 같은 table에서는 2:5:9의 비율로 뽑히는거일까요? (a,c)(b,b)(b,c)(c,a)(c,b)(c,c)가 각각 [0.1, 0.1, 0.2, 0.1, 0.2, 0.3] 의 확률로 뽑히는 경우인데, 결국 c가 0.9의 확률로, b가 0.5의 확률로, a가 0.2의 확률로 등장하거든요. 이렇게 되면, 2:5:9이지 않나 싶은데 헷갈리네요 ㅠ
    • 프로필 사진
      2020.03.18 01:51 신고
      앗 X랑 Y를 둘다 a, b, c로 놓았더니 헷갈리기가 쉽네요.
      X는 {a, b, c} 중 하나고, Y는 {ㄱ, ㄴ, ㄷ} 중 하나라고 바꿔서 설명드리도록 할게요. X가 a, b, c일 확률이 1:2:3이고, Y가 ㄱ, ㄴ, ㄷ일 확률도 1:2:3이라고 가정합시다. 그런데 두 사건 X와 Y는 독립이 아니라서 (X, Y)를 뽑을 확률은 (a, ㄷ), (b, ㄴ), (c, ㄱ), (b, ㄷ), (c, ㄴ), (c, ㄷ)에 대해 각각 1:1:1:2:2:3인 상황을 얘기한 거였습니다.
    • 프로필 사진
      2020.03.18 14:50
      항상 감사드립니다. 같은 학교 선배님이신데 정말 매번 큰 도움을 받고 있습니다.
    • 프로필 사진
      2020.03.18 14:57
      혹시 위의 결과물을 모두 볼 수 있는 파이썬 코드를 받을 수 있을까요?
    • 프로필 사진
      2020.03.18 19:08 신고
      아 반갑습니다. 후배님이셨군요. 도움이 될 수 있어서 아주 기쁩니다 ㅎㅎ
      위 결과물들은 본문에 있는 코드 조각들을 각각 돌리시면 확인하실 수 있습니다. 단 블로그에 올린 그림들은 위 코드 실행 결과들을 수작업으로 잘라내고 편집한거라서 코드 실행 결과가 그림과 완전히 동일하지는 않다는 점 양해바랍니다.
  • 프로필 사진
    2020.04.20 02:00
    좋은 설명 감사합니다 :) 블로그에 유익한 글이 많아서 읽던 중에 이번 포스팅을 통해 깁스 샘플링을 알게 되었는데요, 잘 이해을 하지 못해서 질문을 드리게 되었습니다. 깁스 샘플링이 어떤 주어진 다변수 확률변수 시스템에서 특정 영역에 놓인 샘플들을 추출하는 목적으로 받아드리기는 했는데.. 이부분이 아직 관련 문제를 실제로 접하지 못했다보니 피부로 와닿지가 않습니다.. 왜냐하면, 파이썬으로 코드로 본다면,
    wholedata = [ [ x1, y1] , [x2, y2], ...]
    mixitem(wholedata) # 리스트 원소의 순서를 랜덤하게 섞음
    num_sams = 100
    for loopcnt, (samx, samy) in enumerate(wholedata):
    ....if loopcnt == num_sams:
    ........break
    ....if samx < -1 and samx > -2 and samy > 0 and samy < 1:
    ........print("wannabe sample", samx, samy)
    pritn("sampling end")

    뭔가 제가 샘플링이라는 근본적인 의미나 목적을 놓치고 있는 것 같기는 합니다만.. 그래도 한번 질문드려요 ! 감사합니다 : )
    • 프로필 사진
      2020.04.21 17:25 신고
      네 구현하신 코드가 기각샘플링인걸로 보입니다. 기각 샘플링은 차원이 높아지고, 조건이 복잡해지면 if에 걸리는 샘플의 개수가 극도로 적어집니다. (2차원에서는 100개 중에 30~40개가 걸린다면 200차원에서는 몇십만개를 해도 한개도 안 걸릴수 있지요)
      그렇기 때문에 고차원의 복잡한 분포를 샘플링하는데에는 기각샘플링을 사용하기 어렵습니다.