MCMC 수렴 진단: Gelman–Rubin R-hat과 ESS 해석

R-hat은 체인들이 같은 분포를 탐색하는지 묻는다

고차원 분포에서 MCMC를 사용했다면, 표본을 많이 얻었다는 사실만으로 사후 추론을 신뢰할 수는 없다. 여러 체인이 충분히 섞였고 정지분포에 도달했는지 확인해야 한다. Gelman–Rubin 통계량인 R-hat은 이 판단의 중심에 놓인다.

R-hat은 체인 간 분산(Between-chain variance)과 체인 내 분산(Within-chain variance)을 비교한다. 고전식은 R-hat = sqrt(Var⁺ / W)이며, Var⁺ = ((n-1)/n) W + (1/n) B로 계산한다. 여기서 W는 체인 내 평균분산이고, B는 체인 평균 간 분산에 n을 곱한 값이다. R-hat이 1에 가까울수록 수렴 신호가 강해진다. 일반적인 권고 기준은 엄격하게 R-hat < 1.01, 완화하면 < 1.05다.

진단은 R-hat 하나로 끝나지 않는다. 통계적 지표인 R-hat·ESS, trace plot과 rank plot 같은 시각 진단, HMC의 divergences·BFMI를 함께 점검해야 한다.

현대적 R-hat이 잡아내는 문제

현재는 체인을 반으로 나누는 split R-hat을 사용해 비정상성 탐지 민감도를 높이는 방식이 권장된다. rank-normalized, folded-split R-hat은 랭크 정규화와 절대값 처리를 적용해 heavy-tail과 스케일 문제에 더 강건하다.

분포의 중심과 꼬리 영역을 나누어 보는 bulk/tail R-hat도 필요하다. 중심부가 안정적으로 보이더라도 꼬리 추정이 흔들릴 수 있기 때문이다.

진단 입력은 warmup 또는 burn-in이 끝난 뒤의 표본으로 제한한다. 체인은 최소 4개를 구성하고, 초기값은 overdispersed하게 설정한다. 재현성과 체인 간 독립성을 위해 독립 시드와 랜덤 스트림도 사용한다. 벡터나 행렬 파라미터는 전개한 뒤 모든 스칼라 요소별로 R-hat을 확인한다.

임계치를 해석할 때 함께 볼 것

R-hat < 1.01은 엄격한 통과, 1.01~1.05는 경계, >1.05는 추가 샘플링이나 모형 조정이 필요한 상태로 해석한다. 다만 R-hat이 좋더라도 ESS의 bulk/tail 값, 자기상관, divergences 부재, BFMI 정상 범위를 함께 확인해야 한다.

극소수 파라미터만 기준을 벗어난 경우에는 전체 모델을 다시 다루기보다 해당 영역을 우선 재파라미터화할 수 있다. 반대로 모든 체인이 같은 잘못된 모드에 갇히면 R-hat≈1도 가능하다. 초기값 분산을 넓히고 다른 진단을 병행해야 하는 이유다.

식별성이 약하거나 경사가 폭주하고 스케일 격차가 큰 모델도 R-hat 개선을 늦춘다. 비중심화와 표준화가 이런 경우의 대응 수단이 된다. R-hat은 다중 체인을 전제하므로 단일 체인에서는 IACT와 ESS 같은 대체 지표에 의존한다.

샘플링부터 진단 게이트까지

⤢✕

HMC 경고가 있으면 step size, target_accept, 재파라미터화, 스케일 표준화를 먼저 검토한다. 경고가 없더라도 split과 rank-normalized R-hat, bulk/tail ESS를 계산해 기준 충족 여부를 확인한다. 통과하지 못하면 샘플 수·체인 수·초기값 분산을 조정하거나 재파라미터화한 뒤 다시 샘플링한다.

진단 도구가 맡는 역할

지표	성능(계산비용)	확장성(고차원/다수체인)	일관성(이론 근거/강건성)	안정성(이상치/헤비테일)	운영 편의(툴링/자동화)
R-hat (rank-normalized, split)	낮음(O(MN))	높음	높음(분산분해, 이론 정립)	높음(랭크·folded 적용)	매우 높음(기본 제공, 임계치 게이팅 용이)
ESS (bulk/tail)	낮음	높음	높음(IACT 기반)	중간(헤비테일 tail-ESS 필요)	높음
자기상관/IACT	낮음	높음	중간(추정 편향 민감)	중간	높음
Trace/Rank Plot	중간(시각화)	중간	중간(주관 보조)	중간	중간
에너지/BFMI(HMC)	낮음	높음	높음(HMC 물리량)	중간	높음(Stan/PyMC 내장)

ArviZ, Stan/CmdStan, PyMC, NumPyro는 split/rank-normalized R-hat과 ESS를 제공한다. trace plot, rank plot, energy plot, pair plot은 수렴과 혼합을 보조적으로 확인하는 데 쓴다. 운영 환경에서는 모델 파이프라인에 진단 임계치 게이트를 넣고, 실패 시 자동 튜닝과 재실행으로 연결할 수 있다.

모델 유형에 따라 달라지는 대응

계층 베이지안 모델에서는 그룹 수준의 스케일 차이로 체인 혼합이 저하되고 R-hat>1.05가 지속될 수 있다. 비중심화(non-centered) 재파라미터화와 공변량 표준화, target_accept 0.9~0.95 상향을 적용하면 split R-hat<1.01, bulk-ESS 증가, 예측 신뢰구간 안정화로 이어질 수 있다.

시계열 상태공간 모형은 긴 상관구조로 IACT가 증가하고 ESS가 감소하며 tail R-hat이 경계에 놓일 수 있다. 상태 변수를 차분하거나 표준화하고, 관측 오차와 상태 전이 스케일 사전을 강화하며, 더 긴 warmup을 사용한다. tail-ESS가 개선되면 꼬리 추정이 안정돼 극단값 예측 오류가 줄어든다.

멀티모달 잠재공간 모형에서는 체인별 모드 할당이 어긋나 R-hat이 상승한다. 초기값을 모드별로 명시적으로 분산시키고, 병렬 tempering 또는 가벼운 tempering 프라이어 도입을 검토한다. 모드 교차 탐색이 늘어나면 R-hat과 rank plot도 개선될 수 있다.

진단 자동화는 불량 러닝을 일찍 차단해 재샘플링 횟수를 20~40% 감소시키는 효과가 보고된다. 불필요하게 긴 체인 연장을 막아 총 샘플 대비 유효샘플 비율(ESS/N)을 높이고, 과적합·비정상 체인을 배제해 추론 일관성도 높인다. 임계치 기반 품질 게이트는 모델 배포 프로세스의 표준화에도 연결된다.

PyMC와 ArviZ로 R-hat·ESS 게이트 만들기

전제조건은 Python 3.10+, pymc 5.x, arviz 0.16+, numpy 1.23+이며 CPU 병렬 체인 4개를 권장한다.

설치는 다음과 같다.

• pip install pymc arviz

# 환경: Python 3.11, pymc==5.10+, arviz==0.16+
import pymc as pm
import arviz as az
import numpy as np

rng = np.random.default_rng(42)
y = rng.normal(loc=0.0, scale=1.0, size=200)

with pm.Model() as model:
    mu = pm.Normal("mu", 0, 5)
    sigma = pm.HalfNormal("sigma", 2)
    obs = pm.Normal("obs", mu=mu, sigma=sigma, observed=y)

    idata = pm.sample(
        draws=2000,
        tune=1000,
        chains=4,
        cores=4,
        target_accept=0.9,   # HMC 안정화
        random_seed=[1, 2, 3, 4],
        init="jitter+adapt_diag"  # 초기값 분산 확대
    )

# R-hat 및 ESS 요약
summary = az.summary(idata, var_names=["mu", "sigma"], round_to=2)
print(summary[["mean", "sd", "r_hat", "ess_bulk", "ess_tail"]])

# 임계치 점검(자동 게이트)
rhat = az.rhat(idata)
essb = az.ess(idata, method="bulk")
esst = az.ess(idata, method="tail")

def check_diagnostics(rhat, essb, esst, min_ess=1000, max_rhat=1.01):
    rhat_bad = (rhat.to_array() > max_rhat).values.any()
    essb_bad = (essb.to_array() < min_ess).values.any()
    esst_bad = (esst.to_array() < min_ess/2).values.any()
    return not (rhat_bad or essb_bad or esst_bad)

ok = check_diagnostics(rhat, essb, esst)
print("Diagnostics pass:", ok)

이 예시의 실무 기준은 r_hat < 1.01, ess_bulk ≥ 1000, ess_tail ≥ 500이다. 실패하면 draws와 chains를 늘리고, target_accept를 0.9→0.95로 조정하거나 non-centered 재파라미터화와 스케일 표준화를 검토한다.

신뢰성과 비용 사이의 선택

최소 4체인을 병렬로 두고 split/rank-normalized R-hat과 bulk/tail ESS를 함께 확인하는 구성이 기본이다. 초기값은 overdispersed하게 잡고, 데이터와 파라미터를 표준화하며, 비식별성을 줄이기 위한 사전 정규화를 적용한다. divergences, BFMI<0.3~0.4, max_treedepth 히트도 R-hat과 분리하지 않고 해석한다.

R-hat<1.01 같은 엄격한 임계치는 신뢰성을 높이는 대신 샘플 수와 시간 비용을 높인다. 체인을 4→8로 늘리면 혼합 검출력은 높아지지만 메모리와 CPU 비용도 증가한다. 강한 사전이나 재파라미터화는 수렴성을 높일 수 있으나 모델 해석성을 낮출 가능성이 있다.