스무딩·모수적 부트스트랩으로 신뢰구간 안정화하기

복원추출만으로 충분하지 않은 조건

부트스트래핑은 관측값을 복원추출해 통계량의 경험적 분포를 만드는 비모수적 추정 기법이다. 다만 표본이 작거나 분포의 경계값을 다룰 때, 또는 왜도와 첨도가 큰 데이터에서는 기본 방식의 제약이 드러난다. 재표본이 이산적으로 남고, 최댓값·분위수처럼 경계에 민감한 통계량에서는 커버리지가 낮아질 수 있다. 소표본에서는 변동을 과소 또는 과대 추정할 가능성도 있다.

이때 선택할 수 있는 확장 방식이 Smoothed Bootstrap과 Parametric Bootstrap이다. 전자는 재표본의 연속성을 확보하는 쪽에, 후자는 데이터에 맞춘 모수적 구조를 활용하는 쪽에 무게를 둔다.

재표본 생성 접근 비교

Smoothed Bootstrap은 복원추출한 관측치마다 커널 기반의 확률적 잡음을 더한다. 이산성을 완화하고 경계값·분위수·밀도 추정의 편향과 분산을 개선하며, 신뢰구간의 커버리지를 안정화하는 것이 목적이다. 이 방식에서는 가우시안 같은 커널 형태와 대역폭(bandwidth) 선택이 핵심이다.

Parametric Bootstrap은 먼저 모수적 분포나 모형을 가정하고 추정한 뒤, 적합된 모형에서 재표본을 생성한다. 모수적 구조가 데이터와 잘 맞으면 효율을 높이고 소표본에서 분산과 커버리지를 개선할 수 있다. 반대로 모형이 어긋나면 왜곡 위험이 있으므로 적합도 검정, 잔차 진단, QQ-plot 점검이 전제된다.

관점	Classical Bootstrap	Smoothed Bootstrap	Parametric Bootstrap
왜도·경계에서의 성능	커버리지 저하 가능성	커버리지·분산 개선, 분위수·경계에 유리	모형 정합 시 최고 효율, 불일치 시 왜곡 위험
계산량	낮음~중간	중간(노이즈 생성 비용)	중간(시뮬레이션 비용, 모형 복잡도에 비례)
가정 민감도	가정 적음	커널·대역폭 민감	모형 가정 민감, 미스스펙 시 편향 큼
이상치·희소 데이터	이상치 영향 큼	적절한 h로 완화 가능	모형이 이상치에 취약할 수 있음
구현과 튜닝	구현 용이	대역폭 선택 필요	모형 선택·적합·진단 필요

대역폭과 모형을 검증하는 흐름

스무딩의 대역폭 h는 Silverman 규칙, IQR 기반 방법, 교차검증 최소화 기준으로 정할 수 있다. h가 너무 작으면 이산성이 남고 분산이 커질 수 있으며, 너무 크면 지나친 평활로 편향이 증가한다. 초기값은 규칙 기반으로 잡고, 커버리지와 표본 외 검증을 포함한 부트스트랩 진단으로 미세 조정하는 흐름이 적합하다.

모수적 방식에서는 데이터의 지원(support), 왜도·첨도, 도메인 지식을 함께 본다. 생존·대기시간 데이터에는 Gamma 또는 Weibull을 고려할 수 있다. KS/AD 검정, AIC/BIC, 잔차와 QQ-plot으로 적합 상태를 확인하고, 맞지 않으면 Generalized Gamma처럼 더 유연한 분포를 택하거나 로그·Box-Cox 변환, 단계적 모형 업데이트를 검토한다.

신뢰구간은 Percentile, Studentized, BCa(편향·가변성 보정) 방식으로 산출할 수 있다. 경계값이나 강한 왜도가 문제라면 BCa 또는 Studentized를 우선 고려하고, 단순 Percentile 방식은 왜곡 가능성을 염두에 둔다. 반복 수 B는 2,000~10,000을 권장하며, 높은 신뢰수준이나 꼬리확률을 추정할 때는 B를 높인다.

⤢✕

데이터 특성에 따른 선택

분위수·최댓값·최솟값을 추정할 때는 Smoothed Bootstrap으로 이산성을 낮추고 커버리지를 개선할 수 있다. 결함률이나 지연시간처럼 이산형 데이터의 지표를 연속 근사해야 하는 경우에도 스무딩 후 추정의 안정화를 고려할 수 있다.

생존·대기시간, 수명 시험, 로그스케일의 양수 데이터에는 Gamma·Weibull·Lognormal 모형을 사용하는 Parametric Bootstrap이 효율을 높일 수 있다. 회귀, GLM, 시계열에서는 모형별 분포 가정을 둔 모수적 잔차 부트스트랩으로 예측 구간을 산출한다.

스무딩이나 모수화를 적용하면 커버리지가 목표치에 가까워지고 소표본의 분산이 줄어들 수 있다. 모형이 정합할 경우 동일 정확도에 필요한 B를 줄일 가능성도 있다. 진단과 튜닝 절차를 포함하면 재현 가능한 파이프라인을 만들 수 있어 배포와 모니터링에도 유리하다.

Python으로 분위수 신뢰구간 계산

아래 예시는 90% 분위수를 통계량으로 사용한다. Python 3.10+, numpy>=1.24, scipy>=1.10 환경을 전제로 하며, x를 실제 데이터로 교체해 실행할 수 있다.

# Python 3.10+, numpy 1.24+, scipy 1.10+
import numpy as np
from scipy import stats

rng = np.random.default_rng(42)

def silverman_bandwidth(x: np.ndarray) -> float:
    x = np.asarray(x)
    n = len(x)
    if n < 2:
        raise ValueError("표본 크기 n>=2 필요")
    sd = np.std(x, ddof=1)
    iqr = np.subtract(*np.percentile(x, [75, 25]))
    sigma = min(sd, iqr / 1.34) if iqr > 0 else sd
    h = 1.06 * sigma * n ** (-1/5)
    # 하한선 보호
    return max(h, 1e-8)

def percentile_ci(samples: np.ndarray, alpha: float = 0.1):
    lower = np.quantile(samples, alpha/2)
    upper = np.quantile(samples, 1 - alpha/2)
    return lower, upper

def smoothed_bootstrap_ci(x, stat_fn, B=5000, alpha=0.1, h=None):
    x = np.asarray(x)
    n = len(x)
    if h is None:
        h = silverman_bandwidth(x)
    stats_boot = np.empty(B)
    for b in range(B):
        idx = rng.integers(0, n, size=n)
        xb = x[idx] + rng.normal(0.0, h, size=n)
        stats_boot[b] = stat_fn(xb)
    return percentile_ci(stats_boot, alpha), stats_boot

def fit_distribution(x):
    x = np.asarray(x)
    # 간단한 규칙: 양수·우측왜도면 Gamma, 아니면 Normal
    skew = stats.skew(x, bias=False)
    if np.all(x > 0) and skew > 0.5:
        # Gamma with loc=0
        a, loc, scale = stats.gamma.fit(x, floc=0)
        dist = stats.gamma(a=a, loc=loc, scale=scale)
        name = "Gamma"
    else:
        mu, sd = np.mean(x), np.std(x, ddof=1)
        dist = stats.norm(loc=mu, scale=sd)
        name = "Normal"
    return name, dist

def parametric_bootstrap_ci(x, stat_fn, B=5000, alpha=0.1):
    name, dist = fit_distribution(x)
    n = len(x)
    stats_boot = np.empty(B)
    for b in range(B):
        xb = dist.rvs(size=n, random_state=rng)
        stats_boot[b] = stat_fn(xb)
    # 간단 적합도 점검(가능 시 KS)
    try:
        ks_p = stats.kstest(x, dist.cdf)[1]
    except Exception:
        ks_p = np.nan
    return (percentile_ci(stats_boot, alpha), stats_boot, name, ks_p)

# 예시 데이터(우측왜도): 실제 데이터로 교체 권장
x = stats.lognorm(s=0.6, scale=np.exp(2.0)).rvs(size=150, random_state=rng)

# 관심 통계량: 90% 분위수
q = 0.9
stat_fn = lambda s: np.quantile(s, q)

# Smoothed Bootstrap
(ci_sm, boot_sm) = smoothed_bootstrap_ci(x, stat_fn, B=4000, alpha=0.1)[0:2]
print(f"Smoothed Bootstrap {int(q*100)}% 분위수 CI: {ci_sm}")

# Parametric Bootstrap
(ci_pa, boot_pa, model_name, ks_p) = parametric_bootstrap_ci(x, stat_fn, B=4000, alpha=0.1)
print(f"Parametric Bootstrap({model_name}) {int(q*100)}% 분위수 CI: {ci_pa}  | KS p-value={ks_p:.3f}")

# 점추정치
point = stat_fn(x)
print(f"표본 {int(q*100)}% 분위수 점추정: {point:.4f}")

난수 시드를 고정하고 B를 늘리면 결과 재현성을 확보할 수 있다. 스무딩 대역폭은 h × {0.5, 1, 2} 비교로 민감도를 확인한다. 모수적 부트스트랩은 적합도 p-value 하나로 선택하지 않고, 잔차·QQ-plot·도메인 제약을 함께 반영해 판단한다.