Kwangmin Kim - § 11.1~11.3 — Gibbs·Metropolis·Hastings 와 그 결합

1 개요 — 세 절의 공통 질문

Ch.11 § 11.1~11.3 은 MCMC 알고리즘의 수학적 정합성에 답한다.

절	질문
11.1	Gibbs 가 왜 목표 분포에서 샘플하는가
11.2	Metropolis-Hastings 가 왜 목표 분포에서 샘플하는가
11.3	두 알고리즘을 어떻게 조합하며, 그 결합이 여전히 수렴하는가

Overview (02-11-0) 가 “MCMC 가 뭐고 왜 혁명인가” 였다면, 이 포스트는 “이 알고리즘이 옳다는 걸 어떻게 증명하는가”. 증명을 이해하면 알고리즘 변형·결합을 안전하게 할 수 있다.

2 § 11.1 — Gibbs Sampler 정밀 분석

2.1 아이디어 복습

$\theta = (\theta_1, \ldots, \theta_d)$ 를 블록으로 분할. 각 반복에서 각 블록을 나머지 조건부 추출:

\[ \theta_j^t \sim p(\theta_j \mid \theta_{-j}^{t-1}, y) \]

$\theta_{-j}^{t-1}$ 은 “이번 반복에서 이미 업데이트된 것 + 아직 $t-1$ 값인 것”:

\[ \theta_{-j}^{t-1} = (\theta_1^t, \ldots, \theta_{j-1}^t, \theta_{j+1}^{t-1}, \ldots, \theta_d^{t-1}) \]

2.2 정상 분포 유지 — 1 블록 업데이트

핵심 정리: 각 조건부 추출이 결합 분포의 정상성을 유지.

증명 (2 블록 case). $\theta = (\theta_1, \theta_2)$, 정상 분포 $p(\theta_1, \theta_2 \mid y)$. 블록 1 업데이트:

\[ \theta_1^{\mathrm{new}} \sim p(\theta_1 \mid \theta_2, y), \quad \theta_2^{\mathrm{new}} = \theta_2 \]

$(\theta_1, \theta_2) \sim p(\theta_1, \theta_2 \mid y)$ 라고 하자. 업데이트 후 결합 분포:

\[ p(\theta_1^{\mathrm{new}}, \theta_2 \mid y) = \int p(\theta_1^{\mathrm{new}} \mid \theta_2, y) p(\theta_1, \theta_2 \mid y) d\theta_1 \]

$\theta_1^{\mathrm{new}}$ 추출이 $\theta_1$ 에 무관하므로:

\[ = p(\theta_1^{\mathrm{new}} \mid \theta_2, y) \int p(\theta_1, \theta_2 \mid y) d\theta_1 = p(\theta_1^{\mathrm{new}} \mid \theta_2, y) \cdot p(\theta_2 \mid y) = p(\theta_1^{\mathrm{new}}, \theta_2 \mid y) \]

업데이트 후 결합 분포가 같다. 정상성 유지.

$d$ 블록으로 일반화: 각 블록 업데이트가 정상성 유지 → 전체 순환 업데이트도 정상성 유지.

직관 — 왜 조건부 추출이 “자동 정상 분포” 인가

조건부 분포에서 추출하면 마치 “목표에서 처음부터 추출한 것처럼” 보인다. 조건부 정보 $\theta_2$ 는 그대로 유지 → 나머지 $\theta_1$ 만 조건부 정확 재샘플.

더 형식적 관점: $p(\theta_1, \theta_2) = p(\theta_2) p(\theta_1 \mid \theta_2)$ 분해에서 $\theta_2$ 주변 분포는 건드리지 않고, $\theta_1 \mid \theta_2$ 를 정확한 조건부 분포로 재샘플. 주변 + 조건부 = 결합의 두 요소 모두 보존.

이것이 Gibbs 의 우아함. 수용·기각 없이 항상 목표로 샘플링.

2.3 수렴 조건 — Irreducibility

Gibbs 는 정상성 유지하지만 수렴은 별개 문제. 체인이 시작 분포에서 정상 분포로 수렴하려면:

Irreducible: 어떤 $(\theta_1, \theta_2)$ 에서 어떤 $(\theta_1', \theta_2')$ 로 유한 시간 내 도달 가능.
Aperiodic: 거의 자동 성립.

반례 — 병리적 사후:

$p(\theta_1, \theta_2)$ 가 두 섬 모양:

     θ_2
      |
  ##  |
  ##  |
      |  ##
      |  ##
------+------ θ_1

각 섬 내부에서는 Gibbs 이동 가능하지만, 섬 사이는 건널 수 없다 (조건부 추출로 도약 불가). Irreducibility 깨짐.

실제 사례: 다봉 (multimodal) 혼합 모형. Gibbs 가 한 mode 에 갇힘.

2.4 이변량 정규 $\rho = 0.8$ 예제

Gelman BDA Ch.11 의 고전 예제. $\theta = (\theta_1, \theta_2)$, 목표 분포:

\[ \theta \mid y \sim \mathrm{N}\!\left(\begin{pmatrix} 0 \\ 0 \end{pmatrix}, \begin{pmatrix} 1 & 0.8 \\ 0.8 & 1 \end{pmatrix}\right) \]

조건부:

\[ \theta_1 \mid \theta_2 \sim \mathrm{N}(0.8 \theta_2, 1 - 0.8^2) = \mathrm{N}(0.8 \theta_2, 0.36) \]

\[ \theta_2 \mid \theta_1 \sim \mathrm{N}(0.8 \theta_1, 0.36) \]

Gibbs:

$\theta_1^{t} \sim \mathrm{N}(0.8 \theta_2^{t-1}, 0.36)$.
$\theta_2^{t} \sim \mathrm{N}(0.8 \theta_1^t, 0.36)$.

2.5 왜 $\rho$ 가 크면 수렴이 느린가

Gibbs 는 한 단계에서 평균 $0.8 \theta$ 쪽으로 이동, 분산 $0.36$ 의 무작위 성분.

$\rho \to 1$: 조건부 평균이 $\theta$ 에 가까움 (이동 폭 작음), 분산이 0 에 가까움 (흔들림 작음).

결과: 각 업데이트가 거의 움직이지 않음 → 전체 공간 탐색에 많은 반복 필요.

정량: 자기상관 시간 $\tau \propto 1/(1 - \rho^2)$. $\rho = 0.99$ 면 $\tau \approx 50$ — 한 효과적 독립 표본 얻으려면 50 반복.

직관 — “축 정렬 이동” 의 기하학

사후 분포가 $\rho = 0.8$ 로 기울어진 타원. Gibbs 는 수평 (θ_1 변경) 또는 수직 (θ_2 변경) 만 가능.

타원의 장축 을 따르려면 수평 + 수직 교대 이동 — 지그재그. 장축 방향으로 매 단계 이동 거리 = 조건부 분산 = $\sqrt{0.36} = 0.6$. 타원 길이 (장축 $\approx \sqrt{1 + 0.8} \cdot 2 \approx 2.7$) 가로지르려면 많은 단계.

해결 1: 재매개변수화 (§ 11.3). $\phi_1 = \theta_1 + \theta_2$ 같은 새 좌표로 변환.

해결 2: HMC (Ch.12). 임의 방향 점프 가능.

3 § 11.2 — Metropolis & Hastings 정밀 분석

3.1 Metropolis 알고리즘 재정의

대칭 제안 분포 $J$ (즉 $J(a \mid b) = J(b \mid a)$).

제안 $\theta^* \sim J(\theta^* \mid \theta^{t-1})$.
$r = p(\theta^* \mid y) / p(\theta^{t-1} \mid y)$.
수용: 확률 $\min(r, 1)$ 로 $\theta^t = \theta^*$, 아니면 $\theta^t = \theta^{t-1}$.

3.2 정상 분포 증명 — Detailed Balance

전이 밀도 $T(\theta' \mid \theta)$ 가 detailed balance 만족:

\[ p(\theta \mid y) T(\theta' \mid \theta) = p(\theta' \mid y) T(\theta \mid \theta') \]

이 조건이 충분 조건으로 정상 분포 보장.

증명. Detailed balance 가정 하에:

\[ \int p(\theta \mid y) T(\theta' \mid \theta) d\theta = \int p(\theta' \mid y) T(\theta \mid \theta') d\theta = p(\theta' \mid y) \underbrace{\int T(\theta \mid \theta') d\theta}_{= 1} = p(\theta' \mid y) \]

즉 $p(\theta \mid y)$ 가 $T$ 의 정상 분포. $\square$

3.3 Metropolis 가 Detailed Balance 만족 — 증명

Metropolis 의 전이 밀도:

\[ T(\theta' \mid \theta) = J(\theta' \mid \theta) \min\!\left(1, \frac{p(\theta' \mid y)}{p(\theta \mid y)}\right) + \delta(\theta' - \theta) \cdot (1 - A(\theta)) \]

$A(\theta)$ 는 $\theta$ 에서의 평균 수용 확률. 두 번째 항 ($\delta$) 은 “머무는” 확률.

Detailed balance 확인 — 두 경우 분리. $\theta \ne \theta'$ 인 경우 ($\delta$ 항 0).

WLOG $p(\theta' \mid y) \ge p(\theta \mid y)$. 그러면:

$\theta \to \theta'$: $\min(1, r_{\theta \to \theta'}) = 1$. $T(\theta' \mid \theta) = J(\theta' \mid \theta)$.
$\theta' \to \theta$: $\min(1, r_{\theta' \to \theta}) = p(\theta \mid y)/p(\theta' \mid y)$. $T(\theta \mid \theta') = J(\theta \mid \theta') \cdot p(\theta \mid y)/p(\theta' \mid y)$.

양변:

\[ \text{LHS} = p(\theta \mid y) J(\theta' \mid \theta) \]

\[ \text{RHS} = p(\theta' \mid y) \cdot J(\theta \mid \theta') \cdot \frac{p(\theta \mid y)}{p(\theta' \mid y)} = p(\theta \mid y) J(\theta \mid \theta') \]

대칭성 $J(\theta' \mid \theta) = J(\theta \mid \theta')$ 이므로 LHS = RHS. Detailed balance 성립. $\square$

$\theta = \theta'$ 인 경우도 정상성에 영향 없음 ($\delta$ 의 대칭성).

3.4 Metropolis-Hastings 확장

대칭 가정 제거. 수용 비율:

\[ r = \frac{p(\theta^* \mid y) / J(\theta^* \mid \theta^{t-1})}{p(\theta^{t-1} \mid y) / J(\theta^{t-1} \mid \theta^*)} \tag{11.2} \]

3.5 MH Detailed Balance 증명

같은 라벨링. $\theta \to \theta'$ 의 수용 확률 $= 1$ 이 되도록 $\theta' = \arg\max$ 으로 선택 가정.

\[ \text{LHS} = p(\theta \mid y) J(\theta' \mid \theta) \cdot 1 \]

$\theta' \to \theta$ 수용 확률:

\[ \min\!\left(1, \frac{p(\theta)/J(\theta \mid \theta')}{p(\theta')/J(\theta' \mid \theta)}\right) = \frac{p(\theta) J(\theta' \mid \theta)}{p(\theta') J(\theta \mid \theta')} \]

(라벨링에 의해 이 비율 $\le 1$).

\[ \text{RHS} = p(\theta' \mid y) \cdot J(\theta \mid \theta') \cdot \frac{p(\theta) J(\theta' \mid \theta)}{p(\theta') J(\theta \mid \theta')} = p(\theta) J(\theta' \mid \theta) \]

LHS = RHS. $\square$

Hastings 비대칭 보정의 의미: $J$ 의 방향성 편향을 밀도비에 역수로 삽입해 균형 회복.

직관 — 왜 $J(\theta^{t-1} \mid \theta^*) / J(\theta^* \mid \theta^{t-1})$ 비율인가

$J$ 가 비대칭: 예를 들어 “위로 올라가기 쉬움, 내려오기 어려움” 인 제안. 이런 제안 하에서 Metropolis 의 단순 수용 규칙은 위쪽에 과다 샘플 — 균형 깨짐.

Hastings 의 보정: 되돌아오기 확률이 적으면 가는 확률도 그만큼 적게 수용. 수식적으로 $r$ 분모에 $1/J(\theta^{t-1} \mid \theta^*)$ — 되돌아오기 어려우면 분모 커짐 → $r$ 작아짐 → 수용 확률 낮아짐.

결과: 이동의 “쉬움/어려움” 비대칭을 자동 상쇄. 임의의 제안 분포 $J$ 에도 정확한 사후 분포 회복.

3.6 Ideal Jumping Rule

“이상적” $J$: $J(\theta^* \mid \theta) = p(\theta^* \mid y)$ — 사후 자체에서 제안.

이 경우 $r$:

\[ r = \frac{p(\theta^* \mid y) / p(\theta^* \mid y)}{p(\theta^{t-1} \mid y) / p(\theta^{t-1} \mid y)} = 1 \]

항상 수용. 연속 샘플이 독립 (이전과 무관). 하지만 $p$ 에서 추출 가능하면 MCMC 필요 없음 — 모순. 이상 기준선만.

실무: $J$ 를 $p$ 에 근접하게 설계. Laplace 근사 + Hastings 보정.

3.7 Random Walk Metropolis

가장 흔한 선택: $J(\theta^* \mid \theta) = \mathrm{N}(\theta, c^2 \Sigma)$. 대칭. $c, \Sigma$ 가 튜닝 대상.

$\Sigma$ 선택: 사후의 covariance 근사. Laplace 또는 이전 MCMC run 에서 추정.

$c$ 선택: Roberts-Gelman-Gilks (1997) 이론 — 고차원 정규 목표에서 $c^* = 2.38/\sqrt{d}$ 로 수용률 $\approx 23.4\%$ 달성.

실무:

수용률 측정.
너무 낮으면 ($< 10\%$) $c$ 축소.
너무 높으면 ($> 50\%$) $c$ 확대.
목표 $20{-}40\%$.

4 § 11.3 — Gibbs 와 Metropolis 의 결합

4.1 Gibbs 가 MH 의 특수 케이스 — 증명

Gibbs 에서 $\theta_j$ 업데이트 시 제안 분포를 조건부 사후로 선택:

\[ J_j(\theta_j^* \mid \theta_{-j}^{t-1}) = p(\theta_j^* \mid \theta_{-j}^{t-1}, y) \]

(새 값 $\theta_j^*$ 가 현재 $\theta_j^{t-1}$ 에 의존 안 함 — 조건부에서 직접 추출).

MH 수용 비율 계산:

\[ r = \frac{p(\theta_j^*, \theta_{-j}^{t-1} \mid y) / J_j(\theta_j^* \mid \theta_{-j}^{t-1})}{p(\theta_j^{t-1}, \theta_{-j}^{t-1} \mid y) / J_j(\theta_j^{t-1} \mid \theta_{-j}^{t-1})} \]

여기 주의: 역방향 제안 $J_j(\theta_j^{t-1} \mid \theta_{-j}^{t-1})$ 은 조건부가 동일하므로 $p(\theta_j^{t-1} \mid \theta_{-j}^{t-1}, y)$.

대입:

\[ r = \frac{p(\theta_j^*, \theta_{-j}^{t-1} \mid y) / p(\theta_j^* \mid \theta_{-j}^{t-1}, y)}{p(\theta_j^{t-1}, \theta_{-j}^{t-1} \mid y) / p(\theta_j^{t-1} \mid \theta_{-j}^{t-1}, y)} \]

분자 = $p(\theta_{-j}^{t-1} \mid y)$ ($\theta_j^*$ 제거). 분모 = $p(\theta_{-j}^{t-1} \mid y)$. 따라서 $r = 1$. 항상 수용.

결론: Gibbs 에서 조건부 추출은 “MH 에서 수용 확률 1 의 제안”. MH 의 특수 케이스.

직관 — 왜 Gibbs 가 “수용 확률 100%” 인가

MH 의 수용은 “제안이 좋은 방향인가” 를 판단. 제안이 이미 조건부 사후에서 나왔으면, 이미 “그 조건부 상 최선 분포” 에서 추출된 것 → 거부할 이유 없음.

이것이 Gibbs 의 효율성: 수용/기각 오버헤드 없음. 반면 조건부 샘플링이 가능해야 하는 제약이 있음.

조건부 공액 없으면: Metropolis-within-Gibbs — 해당 블록만 Metropolis, 나머지는 Gibbs.

4.2 Metropolis-within-Gibbs

현실적 혼합 패턴. 블록 1 은 조건부 공액 → Gibbs, 블록 2 는 비공액 → Metropolis.

절차:

$\theta_1^t \sim p(\theta_1 \mid \theta_2^{t-1}, y)$ (Gibbs).
$\theta_2^*$ 제안 from $J(\cdot \mid \theta_2^{t-1})$. Metropolis 수용/기각으로 $\theta_2^t$ 결정.

정당성: 각 블록 업데이트가 조건부 정상 분포 유지 → 전체 체인이 결합 정상 분포 유지.

수학적으로: MH 의 부분 MH 구성도 MH — 수용 규칙이 정확한 비율 유지하면 detailed balance 성립.

4.3 블록 업데이트 전략

4.4 변수 묶기 (Blocking)

$\theta = (\theta_1, \theta_2, \ldots, \theta_d)$ 를 스칼라별로 업데이트하면 상관 큰 변수 쌍에서 느림.

해결: 상관된 변수를 한 블록으로 묶어 동시 업데이트.

예: $(\mu, \tau)$ 가 강한 상관 → 블록 $(\mu, \tau)$ 공동 Metropolis.

효과: 블록 내부 상관 상쇄, 체인 혼합 개선.

단점: 블록 크기가 클수록 수용률 낮아짐 → 튜닝 더 어려움. 실무 절충.

4.5 재매개변수화 — 상관 제거

$\theta$ 가 사후에서 상관되면, 선형 변환 $\phi = A\theta$ 로 상관 제거 가능. 목표: $\phi$ 의 사후 공분산 대각.

예: 이변량 정규 $\rho = 0.8$.

\[ \phi_1 = \frac{\theta_1 + \theta_2}{\sqrt{2}}, \quad \phi_2 = \frac{\theta_1 - \theta_2}{\sqrt{2}} \]

$\phi_1, \phi_2$ 의 공분산: 계산하면 0, 분산 $(1 + 0.8)$ 와 $(1 - 0.8)$.

결과: $\phi$ 좌표에서 Gibbs 는 독립 변수 업데이트 — 한 번에 하나씩이지만 타원 장축·단축 정렬.

4.6 회전의 일반화

$= $ 사후 공분산. 고유값 분해 $\Sigma = Q \Lambda Q^\top$. 변환 $\phi = Q^\top \theta$ 하면 $\mathrm{Cov}(\phi) = \Lambda$ — 대각.

Gibbs 효율: $\phi$ 축이 정확히 사후 타원 축 정렬. 각 Gibbs 업데이트가 “축 방향” 이동 = “사후 타원 축 방향” 이동. 지그재그 없음.

실무 한계: $\Sigma$ 미지 → 반복 추정 후 변환.

4.7 Non-centered Parameterization

계층 모형의 funnel 문제에 대응. 전통:

\[ \theta_j \sim \mathrm{N}(\mu, \tau^2) \]

Non-centered:

\[ \eta_j \sim \mathrm{N}(0, 1), \quad \theta_j = \mu + \tau \eta_j \]

새 변수 $\eta_j$ 가 $(\mu, \tau)$ 와 독립. MCMC 에서 $\eta_j$ 를 샘플 → $\theta_j$ 는 결정적 변환.

효과: $\tau$ 가 작아져도 $\eta_j$ 는 표준 정규 (funnel 제거). Ch.12 에서 HMC 맥락에서 상세.

5 이변량 정규 $\rho = 0.8$ 에서 알고리즘 비교

5.1 Gibbs

자기상관 시간: $\tau_{\mathrm{ac}} \approx 1/(1 - \rho^2) = 1/0.36 \approx 2.8$. 상대적 효율 $\approx 36\%$.

5.2 Metropolis with $\mathrm{N}(0, 0.2^2 I)$ 점프

Gelman et al. Figure 11.1 의 설정. 작은 점프 → 수용률 높지만 탐색 느림.

자기상관 시간 $\tau_{\mathrm{ac}}$ 측정 — 대략 50~100 (매우 느림). 의도적으로 튜닝 부족 예시.

5.3 Metropolis with 적응 점프

$J = \mathrm{N}(\theta^{t-1}, c^2 \Sigma)$ where $\Sigma$ = 추정 사후 공분산, $c = 2.38/\sqrt{2}$.

자기상관 시간 $\tau_{\mathrm{ac}} \approx 5\sim10$. Gibbs 와 비교 가능, 튜닝 필요.

5.4 재매개변수화 후 Gibbs

$\phi = Q^\top \theta$ 변환. 사후 대각 → Gibbs 완벽 작동.

자기상관 시간 $\tau_{\mathrm{ac}} \approx 1$. 실질적 독립 샘플.

5.5 순서 비교

방법	자기상관 시간	$n_{\mathrm{eff}}/n$
Naive Gibbs ($\theta$ 좌표)	2.8	36%
Metropolis 작은 점프	50~100	1~2%
Metropolis 적응 점프	5~10	10~20%
재매개변수화 Gibbs	1	100%

교훈: 좌표 선택 + 알고리즘 선택이 곱으로 작용. 재매개변수화의 위력.

6 Bioassay 로지스틱 회귀 — Gibbs 실패 & Metropolis 성공

Ch.3.7 의 bioassay 예제. 모형:

\[ y_i \mid n_i, \theta_i \sim \mathrm{Bin}(n_i, \theta_i), \quad \theta_i = \mathrm{logit}^{-1}(\alpha + \beta x_i) \]

6.1 Gibbs 의 곤란

조건부 $p(\alpha \mid \beta, y)$ 는 로지스틱 우도 $\times$ prior. 표준 분포 아님 — 직접 추출 불가.

→ Gibbs 적용 불가. Metropolis 또는 Metropolis-within-Gibbs 필요.

6.2 Metropolis 적용

전체 블록 Metropolis:

$(\alpha^*, \beta^*) \sim \mathrm{N}(\alpha^{t-1}, \beta^{t-1}; c^2 \Sigma)$ 제안.
$r = p(\alpha^*, \beta^* \mid y) / p(\alpha^{t-1}, \beta^{t-1} \mid y)$.
$\min(r, 1)$ 로 수용.

$\Sigma$ 초기값: Laplace 근사의 covariance. 수용률 모니터하며 $c$ 조정.

6.3 One-at-a-time Metropolis-within-Gibbs

$\alpha^* \sim \mathrm{N}(\alpha^{t-1}, c_\alpha^2)$. $r_\alpha$ 로 수용/기각.
$\beta^* \sim \mathrm{N}(\beta^{t-1}, c_\beta^2)$. $r_\beta$ 로 수용/기각.

장점: 각 차원 독립 튜닝. 단점: $\alpha, \beta$ 강한 상관이면 느림.

Bioassay 에선 $\alpha, \beta$ 가 실제로 상관 ($\alpha + \beta x$ 관계). 전체 블록 Metropolis 가 유리.

7 코드 — 이변량 정규 Gibbs & Metropolis

7.1 Gibbs Sampler

import numpy as np
import matplotlib.pyplot as plt

rng = np.random.default_rng(11)
n_iter = 5000
rho = 0.8

# Gibbs
theta = np.zeros((n_iter, 2))
theta[0] = [2.5, 2.5]  # 시작점

for t in range(1, n_iter):
    # theta_1 | theta_2 ~ N(rho * theta_2, 1 - rho^2)
    theta[t, 0] = rng.normal(rho * theta[t-1, 1], np.sqrt(1 - rho**2))
    # theta_2 | theta_1 (new) ~ N(rho * theta_1, 1 - rho^2)
    theta[t, 1] = rng.normal(rho * theta[t, 0], np.sqrt(1 - rho**2))

# Burn-in 50% 제거
theta_gibbs = theta[n_iter//2:]
print(f"Gibbs 결과:")
print(f"  평균: {theta_gibbs.mean(axis=0).round(3)}")
print(f"  공분산:\n{np.cov(theta_gibbs.T).round(3)}")

7.2 Metropolis (작은 점프, 의도적 비효율)

def log_target(theta, rho):
    # 이변량 정규 로그 밀도 (상수 제외)
    t1, t2 = theta
    return -0.5 / (1 - rho**2) * (t1**2 - 2*rho*t1*t2 + t2**2)

theta_m = np.zeros((n_iter, 2))
theta_m[0] = [2.5, 2.5]
scale = 0.2
accept = 0

for t in range(1, n_iter):
    # 제안
    proposal = theta_m[t-1] + rng.normal(0, scale, size=2)
    log_r = log_target(proposal, rho) - log_target(theta_m[t-1], rho)
    if np.log(rng.uniform()) < log_r:
        theta_m[t] = proposal
        accept += 1
    else:
        theta_m[t] = theta_m[t-1]

print(f"\nMetropolis (scale={scale}):")
print(f"  수용률: {accept/n_iter:.3f}")
theta_metro = theta_m[n_iter//2:]
print(f"  평균: {theta_metro.mean(axis=0).round(3)}")
print(f"  공분산:\n{np.cov(theta_metro.T).round(3)}")

7.3 자기상관 시간 비교

from scipy.signal import correlate

def autocorr_time(x, max_lag=100):
    x = x - x.mean()
    n = len(x)
    r = correlate(x, x, mode="full")[n-1:n+max_lag]
    r = r / r[0]
    # 자기상관이 0.05 아래로 떨어지는 시차까지 합
    tau = 1 + 2 * r[1:np.argmax(r < 0.05)].sum()
    return tau

tau_gibbs = autocorr_time(theta_gibbs[:, 0])
tau_metro = autocorr_time(theta_metro[:, 0])

print(f"\n자기상관 시간 비교:")
print(f"  Gibbs:      tau = {tau_gibbs:.1f}, n_eff = {len(theta_gibbs)/tau_gibbs:.0f}")
print(f"  Metropolis: tau = {tau_metro:.1f}, n_eff = {len(theta_metro)/tau_metro:.0f}")

기대 결과: Gibbs $\tau \approx 3$ (이론값), Metropolis $\tau \approx 50{-}100$ (의도적 작은 점프). Gibbs 가 약 20 배 효율적 (이 특수 상황에서).

7.4 재매개변수화

# 선형 변환: phi_1 = (theta_1 + theta_2)/sqrt(2), phi_2 = (theta_1 - theta_2)/sqrt(2)
Q = np.array([[1, 1], [1, -1]]) / np.sqrt(2)
phi_samples = theta_gibbs @ Q.T

print(f"\n재매개변수화 후:")
print(f"  phi 평균: {phi_samples.mean(axis=0).round(3)}")
print(f"  phi 공분산:\n{np.cov(phi_samples.T).round(3)}")
# 기대: diagonal ~ [1.8, 0.2]

$\phi$ 좌표에서 공분산이 대각 ($\rho_\phi \approx 0$) → Gibbs 가 $\phi$ 에서 훨씬 효율적일 것.

7.5 Bioassay Metropolis 실전

# Bioassay 데이터
x_bio = np.array([-0.86, -0.30, -0.05, 0.73])
n_bio = np.array([5, 5, 5, 5])
y_bio = np.array([0, 1, 3, 5])

def log_posterior_bioassay(params, x, n, y):
    alpha, beta = params
    eta = alpha + beta * x
    theta = 1 / (1 + np.exp(-eta))
    # 로그 이항 우도
    ll = (y * np.log(theta + 1e-12) + (n - y) * np.log(1 - theta + 1e-12)).sum()
    # 균등 prior
    return ll

# Metropolis
n_iter = 10_000
params = np.zeros((n_iter, 2))
params[0] = [0.8, 8]  # MLE 근사 시작점
scale_cov = np.array([[1, 0.5], [0.5, 3]])  # Laplace 추정
L = np.linalg.cholesky(scale_cov)
accept_bio = 0

for t in range(1, n_iter):
    proposal = params[t-1] + L @ rng.normal(0, 1, size=2)
    log_r = log_posterior_bioassay(proposal, x_bio, n_bio, y_bio) - log_posterior_bioassay(params[t-1], x_bio, n_bio, y_bio)
    if np.log(rng.uniform()) < log_r:
        params[t] = proposal
        accept_bio += 1
    else:
        params[t] = params[t-1]

print(f"\nBioassay Metropolis:")
print(f"  수용률: {accept_bio/n_iter:.3f}")
print(f"  alpha 평균: {params[n_iter//2:, 0].mean():.2f}")
print(f"  beta 평균: {params[n_iter//2:, 1].mean():.2f}")
print(f"  LD50 (-alpha/beta) 평균: {(-params[n_iter//2:, 0] / params[n_iter//2:, 1]).mean():.3f}")

예상: $\alpha \approx 1$, $\beta \approx 10$, LD50 $\approx -0.1$. Ch.3.7 의 결과와 일치.

8 세 절을 관통하는 직관 모음

8.1 “정상성 유지” 가 모든 알고리즘의 공통 척도

Gibbs·Metropolis·MH 모두 detailed balance 로 정상성 보장. 이를 이해하면:

새 알고리즘 안전성 판단 (Hamiltonian MC, slice, NUTS 모두 detailed balance).
알고리즘 결합 시 검증 (Metropolis-within-Gibbs 가 왜 작동하는가).
변형 시 교정 (Hastings 가 비대칭 $J$ 를 어떻게 보정).

Detailed balance 가 MCMC 설계의 핵심 언어.

8.2 $J$ 선택이 수렴 속도 지배

정상 분포 도달은 보장되지만 얼마나 빨리 가 $J$ 에 달림.

$J$ 가 $p$ 에 가까움 → 수용률 높고 이동 빠름.
$J$ 가 $p$ 와 다름 → 수용률 낮거나 이동 느림.

이게 Ch.12 의 HMC 가 $p$ 의 기울기를 이용해 데이터 기반 $J$ 를 만드는 이유.

8.3 재매개변수화는 “무료 가속”

알고리즘 바꾸지 않고 좌표만 변경 → 큰 성능 향상. 계산 비용 거의 없음.

현대 베이즈 실무: 먼저 재매개변수화 시도, 그 다음 알고리즘 선택.

9 실전 체크리스트

§ 11.1~11.3 의 교훈을 실무 절차로:

조건부 공액 확인 — 있으면 Gibbs 우선.
블록 구조 파악 — 상관 변수 묶기.
비공액 블록은 Metropolis — Metropolis-within-Gibbs.
Laplace 로 초기 $\Sigma$ — Metropolis 점프 공분산 시작점.
$c = 2.38 / \sqrt{d}$ 시도 — Roberts-Gelman-Gilks 최적.
수용률 20~40% 목표 — 튜닝 지표.
자기상관 시간 측정 — 효율 진단.
사후 상관 크면 재매개변수화 — 회전 또는 스케일.
Funnel 구조면 non-centered — 계층 모형 표준.
여러 알고리즘 비교 — 같은 문제에 Gibbs, Metropolis, 재매개변수화 병렬 실행해 $n_{\mathrm{eff}}$ 비교.

10 관련 주제

선행 지식

Ch.11 Overview (02-11-0) — MCMC 기본 지도
Ch.10 심화 (02-10-1) — Rejection sampling (MH 와 구조 유사)
Ch.3 Multiparameter Models — 조건부 공액 구조
Ch.5 Hierarchical Models — Gibbs 의 자연스러운 응용

Ch.11 세부 절 (후속 작성 예정)

02-11-2-* — § 11.4~11.6 심화 (수렴 진단·ESS·8 학교)
02-11-3-* — § 11.7~11.8 심화 (문헌·연습)

후속 주제

Ch.12 Efficient MCMC — HMC, NUTS, 자동 튜닝
Ch.13 Variational Inference — MCMC 대안

관련 개념

Geman & Geman (1984) — Gibbs sampler 원저 (이미지 처리)
Gelfand & Smith (1990) — 베이즈 Gibbs 도입
Metropolis et al. (1953) — 원자력 물리에서 원 알고리즘
Hastings (1970) — 일반화
Roberts, Gelman, Gilks (1997) — 최적 수용률
Papaspiliopoulos, Roberts, Sköld (2007) — Non-centered parameterization 이론
Liu (2001), Monte Carlo Strategies — 종합 참고서
Brooks, Gelman, Jones, Meng (2011), Handbook of MCMC — 현대 표준