Metroopolis Hastings법에 기초한 설명 및 Pythhon의 실시

MCC(마르코프 체인 몬테카로파)


어떤 분포\pi(x)에서 견본을 얻으려면\pi(x)가 고차원 분포된 상황이거나 해석하기 어려운 상황에서 직접적인 방법으로 견본을 잘 얻기 어렵다.
MCMC(마르코프 프랜차이즈 몬테카로파)는 이런 상황을 처리하는 방법으로 활용된다.이것은\pi(x)를 안정적으로 분포하는 마코프 체인을 시뮬레이션하는 알고리즘이다.이 알고리즘을 반복하여 얻은 값은\pi(x)의 견본에 가깝다.
마르코프 프랜차이즈는 P(X i| X 1,..., X{i-1})=P(X i| X{i-1}) 같은 어떤 상태가 이전 상태에만 의존해 결정되는 확률 과정을 말한다.
안정분포란 확률 과정을 거쳐 다음 상태로 옮겨도 변하지 않는 분포를 말한다.마르코프 사슬을 중복 분포(극한 분포)로 수렴하면 안정적인 분포가 될 것으로 알려져 있다.

Metroopolis Hastings법


Metroopolis Hastings 방법은\pi(x)를 안정적으로 분포하는 마르코프 체인을 시뮬레이션하는 방법 중 하나입니다.
Metropolis Hastings 방법에서 사용자는 제안서의 분포를 지정해야 한다.'건의' 는 견본이 x일 때 다음 견본 y가 어느 위치에 있어야 하는지 표시하고 그 분포는 Q (y | x) 로 표시됩니다.
간단한 예로 y는 x의 위치 부근을 무작위로 선택한다.예를 들어 정적 분포에 따라 무작위 수를 사용한다
y | x\sim N(x, 1)
화중.다시 말하면 분포 Q(y|x)는
Q(y|x) =\dfrac{1}{\sqrt{2\pi}}\exp\left( -\dfrac{1}{2}(y-x)^2\right)
의 형식.
제시된 분포 Q(y|x)를 사용하면 Metropolis Hastings 방법은 다음과 같은 알고리즘이다.
  • 초기값 X1 = x_적절히 결정 1

  • t\geq1에 대해 다음 충분한 횟수를 반복합니다

  • Q(y|x t) 샘플 y의 값(권장)

  • 계산 A=\min\left(1\dfrac{\pi(y)Q(xt|y)}{pi(xt)Q(y|xt)}\right)
  • 확률A로 제안을 접수, x{t+1}=y.확률 1-A로 제안 거절x{t+1} = x_t
  • 에 머무르다

    이루어지다


    예로 들다
    \pi(x) =\begin{cases}\exp(-x) & (x\geq 0)\\0 & (x < 0)\end{cases}
    이런 분포에서 견본을 얻다.물론 이렇게 간단한 분포라면 굳이 MCMC를 사용할 필요는 없지만 결국 간단한 예일 뿐이다.
    제시된 분포 Q(y|x)는 상기 예에서 열거한 정적 분포 랜덤 수에 따라 이동한다.정적 분포 랜덤수의 경우 Q(y|x)의 공식적인 형태로 Q(xt|y)=Q(y|xt)를 구성하기 때문에 A
    A =\min\left(1,\dfrac{\pi(y)}{\pi(x_t)}\right)
    더 쉽게 쓸 수 있어요.A를 이 식으로 바꾸는 방법을 Metroopolis법이라고 하는데 MetropolisHastings법의 전신(Hastings가 나중에 일반적인 Q(y|x)로 확장됨)이다.
    파이톤에서 이 메트로폴리스 방법을 실현한 것은 다음과 같다.
    import numpy as np
    
    def pi(x):
        if x >= 0:
            return np.exp(-x)
        else:
            return 0
    	
    xs = []
    xs.append(10) # 初期値を適当に 10 と定めた
    
    for i in range(10000):
        y = np.random.normal(xs[i], 1) # 提案値
        a = min(1, pi(y) / pi(xs[i]))
    
        u = np.random.rand() # 0 ~ 1 の一様乱数を生成
    
        if u < a: # 受諾
            xs.append(y)
        else: # 拒否
            xs.append(xs[i])
    
    그러면 이렇게 계산된 xs는\pi(x)의 샘플로 여겨지기 때문에 이 직사각형은\pi(x)에 가까운 형상이라고 볼 수 있다.이거 확인해봐.
    import matplotlib.pyplot as plt
    %matplotlib inline
    
    plt.hist(xs, bins=20)
    

    확실히 상승세를 보이는 지수 분포 형태임을 알 수 있다.
    또 첫 번째 100개 샘플과 다음 100개 샘플을 비교한다.
    plt.hist(xs[0:100], bins=20)
    

    plt.hist(xs[100:200], bins=20)
    

    후자라면 오른쪽 어깨의 하강을 똑똑히 확인할 수 있을 것 같다.이렇게 시험 횟수를 늘려서\pi(x)에 더 가까운 샘플을 알 수 있다.

    References

  • 마르코프 사슬의 안정적 분포와 극한 분포
  • 제11장 표본 추출법
  • MCMC 튜토리얼부터 다봉성 분포 처리 및 응용
  • The Metropolis Hastings Algorithm
  • 좋은 웹페이지 즐겨찾기