13.2 켤레 모델과 베이지안 선형회귀

실습 · 연습문제 2개

베이지안 선형회귀의 사후분포와 예측분포 bayes-linreg-predictive

특징 행렬 Φ∈RN×D\boldsymbol\Phi \in \mathbb{R}^{N \times D}와 관측 y∈RN\mathbf{y} \in \mathbb{R}^N에 대해

y∣w∼N(Φw, σ2I),w∼N(0, α−1I)\mathbf{y} \mid \mathbf{w} \sim \mathcal{N}(\boldsymbol\Phi\mathbf{w},\ \sigma^2\mathbf{I}), \qquad \mathbf{w} \sim \mathcal{N}(\mathbf{0},\ \alpha^{-1}\mathbf{I})

인 모델을 생각한다. 두 함수를 완성한다.

  1. posterior(Phi, y, alpha, sigma): 가중치 사후분포 N(m,S)\mathcal{N}(\mathbf{m}, \mathbf{S})의 (m, S)를 반환한다.

    S=(αI+1σ2Φ⊤Φ)−1,m=1σ2SΦ⊤y\mathbf{S} = \Big(\alpha\mathbf{I} + \frac{1}{\sigma^2}\boldsymbol\Phi^\top\boldsymbol\Phi\Big)^{-1}, \qquad \mathbf{m} = \frac{1}{\sigma^2}\mathbf{S}\boldsymbol\Phi^\top\mathbf{y}
  2. predictive(Phi_new, m, S, sigma): 새 입력들의 특징 Phi_new (M,D)(M, D)에서 예측분포의 평균 (M,)(M,)과 분산 (M,)(M,)을 반환한다. 분산은 관측 잡음 σ2\sigma^2과 가중치 불확실성 ϕ⊤Sϕ\boldsymbol\phi^\top\mathbf{S}\boldsymbol\phi의 합이다.

인자 shape
Phi (N, D)
y (N,)
Phi_new (M, D)
반환 m, S (D,), (D, D)
반환 평균, 분산 (M,), (M,)

테스트는 다음을 확인한다.

  • 사후평균이 정규화 강도 λ=σ2α\lambda = \sigma^2\alpha인 ridge 해와 같다.
  • 데이터가 없으면(N=0N = 0) 사후분포가 사전분포와 같다.
  • 예측 분산은 σ2\sigma^2보다 작아지지 않고, 데이터에서 먼 입력에서 더 크다.
  • 가중치를 사후분포에서 뽑아 만든 예측의 분산(Monte Carlo)과 일치한다.

힌트: 예측 분산을 입력마다 반복문으로 계산해도 되지만, np.einsum("md,de,me->m", Phi_new, S, Phi_new)나 np.sum((Phi_new @ S) * Phi_new, axis=1)로 한 번에 계산할 수 있다.

노출 시간이 다른 계수 데이터의 Gamma–Poisson 갱신 gamma-poisson-update

검출기로 계수를 여러 번 측정했다. ii번째 측정은 노출 시간 tit_i [s] 동안 yiy_i개를 셌고, yi∼Poisson(λti)y_i \sim \mathrm{Poisson}(\lambda t_i)이다. 계수율 λ\lambda [1/s]에 사전분포 Gamma(a0,b0)\mathrm{Gamma}(a_0, b_0)(shape a0a_0, rate b0b_0, 밀도 ∝λa0−1e−b0λ\propto \lambda^{a_0-1}e^{-b_0\lambda})를 둔다. 세 함수를 완성한다.

  1. posterior(a0, b0, counts, exposure_s): 사후분포 Gamma(a,b)\mathrm{Gamma}(a, b)의 (a, b)를 반환한다.
  2. credible_interval(a, b, level=0.95): 사후분포의 등꼬리 신용구간 (lo, hi). 양쪽 꼬리에 (1−level)/2(1 - \text{level})/2씩 남긴다.
  3. prob_rate_above(a, b, threshold): 사후확률 Pr⁡[λ>threshold∣D]\Pr[\lambda > \text{threshold} \mid \mathcal{D}].

힌트: 우도를 λ\lambda의 함수로 정리하면 λ∑yie−λ∑ti\lambda^{\sum y_i} e^{-\lambda \sum t_i} 꼴이다. scipy.stats.gamma는 rate 대신 scale = 1/b를 받는다. 노출 시간이 긴 측정일수록 사후분포에 더 크게 기여해야 하며, 테스트가 이를 확인한다.