부록 I — 모형선택의 정보 기준

I.1 Kullback-Leibler 정보

앞에서 소개한 AIC 는 두 개의 분포에 대한 거리를 반양하는 정보 기준에 의하여 유도된 모형선택의 측도이다. 이제 AIC 가 어떻게 두 개의 분포에 대한 거리에서 유도되는지 알아보자.

정의 I.1 (Kullback-Leibler 정보) 두 개의 분포 \(F\) 와 \(G\) 가 있다고 가정하고 각 분포에 대한 확률밀도함수가 \(f\) 와 \(g\) 로 주어졌다. KL-정보(Kullback-Leibler information) 는 두 개의 분포 \(F\) 와 \(G\) 의 거리를 다음과 같이 정의하는 정보기준이다.

\[ \begin{aligned} I(g;f) & = E_G \left [ \log \left \{ \frac{g(y)}{f(y)} \right \} \right ] \\ & = \int \log \left \{ \frac{g(y)}{f(y)} \right \} g(y) dy \\ & = \int \left [ \log g(y)-\log f(y) \right ] g(y) dy \end{aligned} \]

\(\blacksquare\)

두 분포의 거리를 나타내는 KL-정보는 다음과 같은 성질을 가진다.

  • \(I(g;f) \ge 0\)

  • 만약 \(I(g;f) = 0\) 이면 \(g(y)=f(y)\) a.e.

보기 I.1 (정규분포의 거리) 두 개의 정규분포 \(F \equiv N(\mu,\sigma^2)\) 와 \(G\equiv N(\xi, \tau^2)\) 을 고려하자.

분포 \(G\) 에서 \((y-\mu)^2\) 의 기대값은 다음과 같이 주어지므로

\[ E_G (y-\mu)^2 = E_G(y-\xi+\xi -\mu)^2 =\tau^2 +(\xi-\mu)^2 \]

분포 \(G\) 를 가정하고 확률밀도함수 \(f\) 의 로그에 대한 기대값은 다음과 같다.

\[ \begin{aligned} E_G[ \log f(y) ] &= E_G \left [ -\frac{1}{2} \log(2\pi \sigma^2) -(y-\mu)^2/(2\sigma^2) \right ] \\ &= -\frac{1}{2} \log(2\pi \sigma^2)-[\tau^2 +(\xi-\mu)^2]/(2\sigma^2) \end{aligned} \]

유사한 방법으로 분포 \(G\) 에서 확률밀도함수 \(g\) 의 로그에 대한 기대값도 다음과 같이 구할 수 있다.

\[ E_G[ \log g(y) ] = -\frac{1}{2} \log(2\pi \tau^2)-\frac{1}{2} \]

위에서 구한 제곱합의 기대값을 이용하면 두 개의 정규분포 \(F\) 와 \(G\) 의 KL-정보는 다음과 같이 주어진다.

\[ I(g;f) = E_G[ \log g(y) ]-E_G[ \log f(y) ] = \frac{1}{2} \left \{ \log \frac{\sigma^2}{\tau^2} + \frac{\tau^2 +(\xi-\mu)^2}{\sigma^2} -1 \right \} \]

먼저 \({\pmb y} = \{ y_1, y_2,\dots, y_n \}\) 가 참분포(true distribution) \(G(y)\) or \(g(y)\) 에서 독립적으로 얻은 확률변수라고 하자.

이제 \(F(y)\) or \(f(y)\) 는 얻어진 자료에 대한 모형으로 사용하고자 하는 분포이며 이를 후보 모형(candidate distribution)이라고 하다. 일반적으로 고려하는 분포들을 모아놓은 집합을 후보 모형군(family of candidate distributions) 이라고 하며 \(\{ f(y | \theta) | \theta \in \Theta \}\) 라고 표기한다. 이런 후보 모형군은 분포를 결정하는 모수들을 모아 놓은 모수 집합 \(\Theta\)으로 표시하기도 한다.

이제 후보 모형군에 속하는 임의의 분포 \(f(y) = f(y|\theta)\) 와 참모형 \(g(y)\) 의 KL-정보는 다음과 같다.

\[ I(g;f) = E_G[ \log g(y) ]-E_G[ \log f(y) ] \]

위에서 정의한 \(I(g;f)\)의 값이 작을수록 좋은 것이며 이는 고려한 후보 분포 \(f(y)\) 가 참모형 \(g(y)\) 에 더 가깝다는 의미이기 떄문이다.

I.2 가능도 함수

이제 자료에 대한 후보 분포를 \(F(y)\) 또는 \(f(y)\) 라고 하고 로그 가능도 함수의 기대값을 고려하자. 이 때 기대값은 참분포에서 계산된 기대값이다.

\[ E_G[ \log f(y) ] = \int \log f(y) g(y) dy \tag{I.1}\]

실제 자료를 분석하는 경우 참분포를 알 수 없기 때문에 로그 가능도 함수의 기대값 식 I.1 을 구하는 것은 불가능하다. 따라서 참분포에 대한 추정을 하여 구해야 하는데 참분포 \(G\)에 대한 추정량은 자료를 이용하여 구할 수 있는 경험 분포함수(emprical distribution)를 이용할 수 있다.

표본 자료 \({\pmb y} = \{ y_1, y_2,\dots, y_n \}\) 를 이용하여 얻은 참분포 \(G\)에 대한 경험적 추정 분포는 다음과 같다.

\[ \hat G(y) = \frac{1}{n} \sum_{i=1}^n I(y \le y_i) \tag{I.2}\]

이제 자료에서 얻은 경험분포를 사용하여 로그 가능도 함수의 기대값 식 I.1 에 대한 추정량을 구하면 다음과 같다.

\[ E_{\hat G} [ \log f(y) ] = \int \log f(y) d \hat G(y) = \frac{1}{n} \sum_{i=1}^n \log f(y_i) \tag{I.3}\]

대수의 법칙(the law of large numbers)에 의하여 표본의 개수 \(n\) 이 커지먄 다음이 성립한다.

\[ \frac{1}{n} \sum_{i=1}^n \log f(y_i) \rightarrow_{a.e.} E_{ G} [ \log f(y) ] \]

표본에 대한 모수적 확률 모형들을 모아놓은 집합, 모형공간 \(\{f(y| \pmb \theta) | \pmb \theta \in \pmb \Theta \subset R^p \}\). 을 고려하자. 표본자료 \({\pmb y} = \{ y_1, y_2,\dots, y_n \}\)로 부터 얻은 로그가능도함수는 다음과 같이 주어진다.

일단 여기서는 참분포 \(g(y)\) 가 모수적 확률 모형 집합에 속하는 분포라고 가정하자.

\[ g(x) =f(y| \pmb \theta_0) \text{ for some } \pmb \theta_0 \in \pmb \Theta \]

\[ \ell(\pmb \theta) = \sum_{i=1}^n \log f(x_i | \pmb \theta ) \]

최대가능도 추정량 \(\hat {\pmb \theta} = \hat {\pmb \theta}({\pmb y})\)은 다음과 같이 로그가능도함수를 최대로 하는 추정량이다.

\[ \hat {\pmb \theta} = arg \max_{\pmb \theta \in \pmb \Theta} \ell(\pmb \theta) \]

이제 \(\pmb \theta_0\) 을 다음에 주어진 방정식의 근이라고 하자.

\[ \int \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta}f(y | \pmb \theta ) dy =0 \]

위의 식에서 다음과 같이 로그 확률함수를 모수벡터로 미분한 양을 스코어함수(score function) \(u(\pmb \theta; \pmb y)\)이라고 부른다.

\[ u(\pmb \theta; \pmb y) = \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta} \]

따라서 \(\pmb \theta_0\) 을 다음에 주어진 방정식의 근이다.

\[ E_{\theta} [u(\pmb \theta; \pmb y) ] =0 \]

만약 정상적인 조건들(regularity conditions)이 만족하면 다음과 같은 결과를 얻을 수 있다.

  • 가능도 방정식 \(\partial \ell(\pmb \theta)/\partial \pmb \theta =0\) 은 점근적으로(with probability 1) 방정식의 해 \(\hat {\pmb \theta}\) 를 가진다.

  • 최대가능도추정량(MLE) \(\hat {\pmb \theta}\) 는 점근적으로 \(\pmb \theta_0\) 에 수렴한다.

  • 최대가능도추정량은 점근적으로 다음과 같은 정규분포를 따른다.

\[ \sqrt{n} ( \hat {\pmb \theta} - \pmb \theta_0) \rightarrow_d N(0, I(\pmb \theta_0)) \tag{I.4}\]

위의 식에서 \(I(\pmb \theta)\) 는 피셔정보(Fisher information matrix) 이라고 부르며 다음과 같이 정의된다.

\[ I(\pmb \theta) = \int f(y | \pmb \theta ) \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta} \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta^t} dy \]

위에서 최대가능도 추정량에 대한 모든 성질은 참분포 \(g(y)\) 가 모수적 확률 모형들의 집합 \(\{f(y| \pmb \theta) | \pmb \theta \in \pmb \Theta \subset R^p \}\) 에 속한다고 가정하였다

\[ g(x) =f(y| \pmb \theta_0)) \text{ for some } \pmb \theta_0 \]

만약 \(g(y)\) 가 우리가 고려하고 있는 모수적 모형들의 집합에 속해있지 않다면 앞에서 구한 최대가능도 추정량의 점근적 성질들은 어떻게 될까?

\[ g(x) \ne f(y| \pmb \theta) \text{ for all } \pmb \theta \]

이제 참분포 \(g(y)\) 가 모수적 확률 모형들의 집합에 속하지 않을 수도 있다고 가정하자.

또한 \(\pmb \theta_0\) 를 다음 방정식의 근이라고 하자.

\[ \int g(y) \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta} dy =0 \]

이러한 가정하에서는 다음이 성립한다.

  • 최대가능도추정량(MLE) \(\hat {\pmb \theta}\) 는 점근적으로 \(\pmb \theta_0\) 에 수렴한다.

  • 최대가능도추정량은 점근적으로 다음과 같은 정규분포를 따른다.

\[ \sqrt{n} ( \hat {\pmb \theta} - \pmb \theta_0) \rightarrow_d N(0, J^{-1}(\pmb \theta_0) I(\pmb \theta_0) J^{-1}(\pmb \theta_0) ) \]

위의 식에서 \(I(\pmb \theta)\) 와 \(J(\pmb \theta)\) 는 다음과 같이 정의되는 양이다.

\[ \begin{aligned} I(\pmb \theta) &= \int g(y) \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta} \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta^t} dy \\ J(\pmb \theta) &= - \int g(y) \frac{\partial^2 \log f(y | \pmb \theta ) } {\partial \pmb \theta \pmb \theta^t} dy \end{aligned} \]

주목할 점은 만약 참분포가 모수적 분포의 집합에 속하면, 즉 \(g(y) =f(y| \pmb \theta_0)\) 이면 \(I(\pmb \theta_0) =J(\pmb \theta_0)\) 이 성립하고 식 I.4 이 성립힌다.

I.3 AIC

이제 \({\pmb y} = \{ y_1, y_2,\dots, y_n \}\) 는 참모형 \(g(y)\)에서 얻어진 독립표본이라고 하자. 모수적 확률 분포의 집합 \(\{f(y| \pmb \theta) | \pmb \theta \in \pmb \theta \subset R^p \}\)을 고려한다. 또한 모르는 모수 \(\pmb \theta\) 는 최대가능도 추정량 \(\hat {\pmb \theta}\) 에 의하여 추정된다고 하자.

이제 추정된 모수로 구한 확률분포 \(f(y|\hat {\pmb \theta})\) 와 참분포 \(g(y)\) 가 얼마나 차이가 있는지 관심이 있으며 이 거리를 K-L 정보 를 이용하여 구하면 다음과 같다.

\[ I(g(z);~f(z|\hat {\pmb \theta})) = E_G [ \log g(z)] - E_G[\log f(z|\hat {\pmb \theta}) ] \tag{I.5}\]

위의 식 I.5 에서 기대값 \(E_G()\) 는 참분포 \(g(z)\)에 추출한 새로운 확률변수 \(z\)에 대한 기대값이며 표본으로 부터 구한 \(\hat {\pmb \theta}=\hat {\pmb \theta}({\pmb y}_n)\) 는 표본 \(\pmb y\)의 함수로서 기대값 \(E_G()\)과 관계없이 고정된 양이다.

위의 식 I.5 에 주어진 K-L 정보에서 앞의 기대값 \(E_G [ \log g(z)]\) 은 언제나 주어진 상수이므로 참분포와 모수적 분포의 거리를 나타내는 양으로 K-L 정보에서 뒤의 기대값이 모수적 분포의 적함도를 반영하는 중요한 측도이다.

\[ E_G[\log f(z|\hat {\pmb \theta}) ] = \int \log f(z|\hat {\pmb \theta}) g(z) dz \tag{I.6}\]

\(I(g(z);~f(z|\hat {\pmb \theta})) \ge 0\) 이므로 위의 식 I.6 에 주어진 양이 크면 클수록 참분포와 거리가 작아지므로 더 좋은 분포의 추정량이라고 말할 수 있다.

여기서 중요한 점은 실제 문제에서는 참분포 \(g(y)\) 를 알 수 없으며 식 I.6 의 값을 추정하려면 참분포 \(g(y)\) 에 대한 추정량이 필요하다. 가장 간단한 추정량은 참분포의 분포함수 \(G\) 를 경험적 표본 분포함수로 추정하는 것이다. 아래는 참분포의 분포함수에 대한 단순 추정량 \(\hat G\) 이다.

\[ E_{\hat G} [\log f(z|\hat {\pmb \theta}) ] = \int \log f(z|\hat {\pmb \theta}) d \hat G(z) = \frac{1}{n} \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) \tag{I.7}\]

사실 식 I.7 는 최대값을 가지는 로그 가능도함수를 \(n\) 으로 나눈 양이다

이제 식 I.7 으로 주어진 추정량으로 식 I.6 를 추정해야 하는데 사실 두 양이 모두 표본 \(y_1, y_2, \dots, y_n\)에 의해 얻어진 \(\hat {\pmb \theta}\)의 함수이다. 따라서 두 통계량 모두 참분포 \(g(y)\) 에서 얻어진 표본 \(y_1, y_2, \dots, y_n\) 도 고려해야 한다.

\[ E_{G(y)} \left [ \frac{1}{n} \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) \right ] =_? E_{G(y)} \left [ E_{G(z)} [\log f(z|\hat {\pmb \theta}) ] \right ] \tag{I.8}\]

불행하게도 위의 두 기대값의 값이 다르기 때문에 식 I.7 에 나타난 추정량은 식 I.6 의 불편추정량이 아니다. 따라서 식 I.7 에 나타난 추정량은 다음과 같이 주어진 편이(bias)를 구해서 보종할 수 있다.

\[ b(G) = E_{G({\pmb y})} \left [ \log f({\pmb y}|\hat {\pmb \theta}) - E_G(\log f(z|\hat {\pmb \theta}(\pmb y))) \right ] \tag{I.9}\]

식 I.9 에 나타난 편이에 대한 해석은 다음과 같이 말할 수 있다.

  • 실제로 추정해야 하는 측도는 \(E_{G(y)} [ E_{G(z)} [\log f(z|\hat {\pmb \theta}) ] ]\)이며 이는 표본 \(\pmb y\) 에서 추정량을 이용하여 새로운 반응변수 \(z\) 를 예측할 때의 측도이다.

  • 하지만 표본 추정량에 근거한 \(E_{G(y)} [ \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) /n ]\)은 표본 \(\pmb y\) 에서 추정량을 이용하여 다시 표본에서 얻은 반응값을 예측하는 측도이다.

  • 따라서 두 측도 사이에는 차이가 존재하며 표본 추정량에 근거한 \(E_{G(y)} [ \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) ]\)는 실제로 과대 추정된다(과적합 발생).

  • 이러한 이유로 가능도함수로 나타난 측도 \(-2 \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta})\) 를 추정된 편이로 보정해주어야 올바른 추론이다.

만약 우리가 식 I.9 에 나타난 편이 \(b(G)\)을 추정할 수 있다면 우리가 찾은 모수적 최적 모형과 참모형의 K-L거리에 근거한 모형의 적합도 \(IC({\pmb y};\hat {\pmb \theta})\) 를 다음과 같이 정의할 수 있다.

\[ \begin{aligned} IC({\pmb X};\hat G) &= -2(\text{log-likelihood of the model} - \text{bias estimator}) \\ &= -2 \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) +2 (\text{bias estimator } b(G)) \end{aligned} \]

본 강의에서는 구하지 않겠지만 식 I.9 에 나타난 편이 \(b(G)\)을 는 다음과 같이 두 행렬의 곱에 대각원소의 합으로 나타난다.

\[ b(G) = tr [ I(\pmb \theta_0) J(\pmb \theta_0)^{-1} ] \tag{I.10}\]

위의 식에서 \(I(\pmb \theta)\) 와 \(J(\pmb \theta)\)는 다음과 같이 정의된 양이다.

\[ \begin{aligned} I(\pmb \theta) &= \int g(y) \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta} \frac{\partial \log f(y | \pmb \theta ) } {\partial \pmb \theta^t} dy \\ J(\pmb \theta) &= - \int g(y) \frac{\partial^2 \log f(y | \pmb \theta ) } {\partial \pmb \theta \pmb \theta^t} dy \end{aligned} \]

만약 참모형 \(g(z)\) 이 다음과 같이 고려한 모수적 모형의 집합 \(\{f(y| \pmb \theta) | \pmb \theta \in \pmb \theta \subset R^p \}\) 에 속한다면

\[ \text{ If } g(y) =f(y|\pmb \theta_0) \quad \text{ for some } \theta_0 , \quad \text{then} \quad I(\pmb \theta) = J(\pmb \theta) \]

이런 조건에서는 식 I.10 에 주어진 편이가 다음과 같이 모수의 개수로 나타난다.

\[ b(G) = tr [ I(\pmb \theta_0) J(\pmb \theta_0)^{-1} ] = tr({\pmb I}_p) = p \]

따라서 K-L 정보기준으로 유도된 참모형과 최대가능도 추정법으로 선택된 분포의 거리로 표현되는 AIC(Akaike Information Criteria)는 다음과 같이 정의된다.

\[ AIC = -2 \sum_{i=1}^n \log f(y_i|\hat {\pmb \theta}) +2p \]

I.4 BIC

베이지안 정보 기준(BIC) 또는 슈바르츠 정보 기준(SIC, SBC, SBIC)은 모형공간에서 최적의 모형을 선택하는 기준으로, 일반적으로 BIC가 낮은 모델이 선호된다. 이는 부분적으로 (AIC)과 밀접한 관련이 있다.

먼저 BIC 를 유도하는 과정을 살펴보려면 베이지안 통계에서 나타나는 모수에 대한 사전분포(prior distribution) \(g(\theta)=p(\theta)\) 을 고려해야 한다.

또한 모형 공간을 \(\mathcal{M}=\left\{m_1, \ldots, m_M\right\}\) 이라고 하자. 또한 \(m \in \mathcal{M}\) 을 모형공간에 속하는 하나의 모형을 나타낸다.

\(\ell(\theta)\) 을 로그 가능도 함수라고 하자. 아레 식에서 \(f(y \mid \theta, m)\) 은 모형 \(m\) 에서 주어진 모수 \(\theta\) 에 대한 반응변수 \(y\) 의 조건부 확률밀도함수이다.

\[ \ell(\theta)=\log f(y \mid \theta, m) \]

다음으로 로그 가능도 함수와 사전분포를 이용하여 함수 \(g\) 와 \(h\) 를 다음과 같이 정의하자.

\[ \begin{aligned} & g(\theta)=p(\theta \mid m ) \\ & h(\theta)=\frac{1}{n} \ell (\theta) . \end{aligned} \tag{I.11}\]

베이지안 통계에 의하면 반응변수 \(y\)에 대한 주변분포(marginal distribution) \(p(y \mid m)\) 는 다음과 같이 주어진다.

\[ \begin{aligned} p(y \mid m) & =\int_{\Theta} f(y \mid \theta, m) p(\theta \mid m) \mathrm{d} \theta \\ & =\int_{\Theta} \exp [n h(\theta)] g(\theta) \mathrm{d} \theta \end{aligned} \tag{I.12}\]

This is an integral suitable for Laplace approximation which states that

이제 위의 적분에 대한 라플라스 근사를 적용하면 다음과 같이 주어진다.

\[ \int_{\Theta} \exp [n h(\theta)] g(\theta) \mathrm{d} \theta=\left(\sqrt{\frac{2 \pi}{n}}\right)^p \exp \left[n h\left(\theta_0\right)\right]\left(g\left(\theta_0\right)\left| H \left(\theta_0\right)\right|^{-1 / 2}+O(1 / n)\right) \tag{I.13}\]

위의 식에서

\(\theta_0\) 는 \(h(\theta)\) 를 최대화하는 값이고 \(H \left(\theta_0\right)\) 는 \(\theta_0\) 에서 계산된 \(h(\theta)\)수의 헤시안 행렬(Hessian matrix) 이다.

우리는 지금 최대가능도 추정법을 다루고 있으므로 위의 식에서 나타난 \(\theta_0\) 는 최대가능도 추정량 \(\hat{\theta}\) 이다.

\[ \hat{\theta}=\underset{\theta}{\arg \max } ~\ell(\theta) . \]

위의 결과에서 식 I.13 를 식 I.11 와 식 I.12 에 적용하면 다음 결과를 얻을 수 있다.

\[ p(y \mid m ) \approx\left(\sqrt{\frac{2 \pi}{n}}\right)^p f(y \mid \hat{\theta}, m) p(\hat{\theta} \mid m)| H (\hat{\theta})|^{-1 / 2} . \tag{I.14}\]

식 I.14 에 로그를 취하고 \(-2\) 다음과 같은 결과를 얻는다.

\[ -2 \log p(y \mid m ) \approx-2 \ell(\hat{\theta})+p \log n-p \log (2 \pi)-2 \log p(\hat{\theta} \mid m)+\log |J(\hat{\theta})| . \tag{I.15}\]

표본의 크기가 커지면(\(n \rightarrow \infty\)) , 식 I.15 의 마지막 3개의 항은 \(O_p(1)\) 으로 나머지 항에 비교하여 무시할 수 있다.

이제 모형 \(\mathcal{M}=\left\{m_1, \ldots, m_M\right\}\) 에서 최적의 모형을 선택하는 기준은 모형에 대한 사후분포 \(p\left(m_j \mid y\right)\) 를 최대로 하는 기준을 사용하는데 이는 우리가 근사한 주변분포 \(p\left(y \mid m_j\right)\) 에 비례하는 것을 이용하여 다음과 같이 모형의 선택 기준으로 BIC 를 정의할 수 있다. \[ \operatorname{BIC}(m)=-2 \log f(y \mid \hat{\theta}, m)+p \log n . \]