PRML [2.3] 가우시안 분포
앞에서도 살펴봤지만, 가우시안 분포(정규 분포)는 연속 변수를 모델하는 분포로 자주 사용된다.
식 2.42는 단일 변수에 대한 가우시안 분포의 형태이다.
\[N(x|\mu, \sigma^2) = \dfrac{1}{(2\pi\sigma^2)^{1/2}}\exp\left\{-\frac{1}{2\sigma^2}(x-\mu)^2\right\} \qquad{(2.42)}\]- $\mu$ : 평균
- $\sigma^2$ : 분산
식 2.43은 다변량 가우안 분포의 형태이다.
\[N({\bf x}|{\pmb \mu}, {\bf \Sigma}) = \dfrac{1}{(2\pi)^{D/2}|{\bf \Sigma}|^{1/2}}\exp\left\{-\frac{1}{2}({\bf x}-{\pmb \mu})^T{\bf \Sigma}^{-1}({\bf x}-{\pmb \mu})\right\} \qquad{(2.43)}\]- $\mu$ : D 차원 평균 벡터
- $\Sigma$ : D x D 공분산 행렬
- $|\Sigma|$ : $\Sigma$ 의 행렬식
중심 극한 정리 (central limit theorem)
라플라스의 중심 극한 정리에 따르면, 여러 개의 확률 변수들의 합에 해당하는 확률 변수는, 합새지는 확률 변수의 개수가 증가함에 따라 점점 가우시안 분포에 수렴한다. 사실, 이건 시행이 있는 모든 확률에 해당하기 때문에 이미 직관적으로 알고 있는 내용과 비슷하다. 주사위 엄청 많이 던지면 한 면이 나올 확률이 1/6에 수렴하겠지.
그림에서 볼 수 있는 것과 같이, $N$이 커지면 가우시안 분포의 모습을 띤다.
따라서, $N$개의 이진 확률 변수 $x$의 관찰값의 합인 $m$에 대한 분포인 이항 분포 역시 $N \rightarrow \infty$ 이 됨에 따라 가우시안의 형태를 띤다는 것을 알 수 있다.
가우시안 분포의 기하학적 형태
$\bf x$에 대한 가우시안 분포의 함수적 종속성은 지수상에서 나타난다. 가우시안에 로그를 취해서 계산하는 경우가 많기 때문에 결국은 지수부를 다루게 되는 것이다.
\[\Delta^2 = ({\bf x}-{\pmb \mu})^T{\bf \Sigma}^{-1}({\bf x}-{\pmb \mu}) \qquad{(2.44)}\]- $\bf x$ : 관측된 데이터 포인트
- $\mu$ : 기댓값 벡터
- $\rm T$ : 전치
- ${\bf \Sigma}^{-1}$ : 공분산의 역행렬
- $\Delta$ : $\mu$ 로부터 $\bf x$까지의 마할라노비스 거리
마할라노비스 거리 (Mahalanobis distance)는 다변량 데이터에 한 점이 데이터 집합의 평균으로부터 얼마나 떨어져있는지 나타내는 척도이다.
$\bf \Sigma$ 가 항등 행렬(모든 변수가 서로 독립적이고 분산이 1인 경우)인 경우에는 이마할라노비스 거리가 유클리디안 거리와 같아진다. $\bf \Sigma$ 가 항등 행렬이면 $\bf \Sigma^{-1}$도 항등 행렬이다. 이때, 식은 아래와 같이 간소화된다.
$\Delta^2 = ({\bf x}-{\pmb \mu})^T({\bf x}-{\pmb \mu})=$ $\bf x$와 $\pmb \mu$ 사이의 유클리디안 거리의 제곱 (각 차원에 대한 단순한 차이의 제곱합)
이 이차식이 상수가 되는 $\bf x$ 공간의 표면에서는 가우시안 분포도 상수가 된다. 이 이차식이 상수가 되는 경우는 모든 점이 동일한 확률 밀도 값을 가질 때이고, 이 공간은 사실상 가우시안 분포의 등고선 (contour, 평균 $\mu$를 중심으로 하는 타원형(다차원에서는 하이퍼타원))이다.
공분산 행렬은 대칭 행렬이며, 이를 고유벡터(eigenvector)와 고유값(eigenvalue)으로 분석될 수 있다. 고유벡터와 고유값은 공분산 행렬의 선형 변환 후에도 변하지 않는 성질을 가지며, 이를 통해 가우시안 분포의 기하학적 형태를 이해할 수 있다.
가우시안 분포의 수학적 변환
고유벡터와 고유값의 활용
가우시안 분포를 분석할 때, 공분산 행렬을 고유벡터의 선형 결합으로 표현할 수 있으며, 이를 통해 분포의 축을 변환하여 차원 간 독립적인 정규 분포를 얻을 수 있다.
이 변환은 가우시안 분포의 원점을 평균($μ$)으로 옮기고 고유 벡터를 축으로 하는 회전 변환이다.
- 빨간 선은 2차원 공간 $\textbf x=(x_1, x_2)$ 상에서의 상수 가우시안 확률 분포의 타원형 표면이다.
- 이 타원의 축들은 공분산 행렬의 고유 벡터 $u_i$에 의해 정의된다.
- 각 축은 고윳값 $\lambda_i$에 대응된다.
야코비안과 좌표 변환
$\bf x$ 좌표계에서 $\bf y$ 좌표계로의 전환은 야코비안(Jacobian) 행렬을 이용한다.
이는 선형 변환 시 발생하는 공간의 부피 변화율을 반영하며, 이를 통해 가우시안 분포를 새로운 좌표계로 표현할 수 있다.
\[p({\bf y}) = p(x)|{\bf J}| = \prod_{j=1}^{D}\dfrac{1}{(2\pi\lambda_j)^{1/2}}\exp\left\{-\dfrac{y_j^2}{2\lambda_j}\right\} \qquad{(2.56)}\]- $p({\bf y})$ : 변환된 좌표계에서의 확률 밀도 함수
- $p(x)$ : 원래 좌표계에서의 확률 밀도 함수
- $|{\bf J}|$ : 좌표 변환의 부피 변화율
식 (2.56)은 가우시안 분포를 $\bf x$ 좌표계에서 $\bf y$ 좌표계로 야코비안 좌표 변환한 결과이다.
이 변환 과정에서 각 차원 $j$에 대해 독립적인 정규 분포의 곱으로 표현되는 가우시안 분포를 얻는다. 이 결과는 $\bf y$ 좌표계에서 각 차원이 서로 독립적인 정규 분포를 따름을 나타낸다. 이는 고유 벡터로의 좌표 변환을 통해 가우시안 분포의 축이 정렬되어 각 차원이 독립적으로 분리될 수 있음을 의미한다.
가우시안 분포의 적률
가우시안 분포의 평균($μ$)과 분산($Σ$)은 모멘트를 통해 계산된다. 식 2.56을 변수 변환 후, 적분하고, 함수의 좌우 대칭성을 사용하여 평행 이동하면 가우시안 분포에서 기댓값이 분포의 중심인 $\mu$와 같음을 다시 한 번 알 수 있다.
\[E[x]=\mu\]이차 모멘트값(두 변수의 곱의 기댓값)은 공분산, 분산과 밀접한 관계가 있기 때문에 분포의 형태와 분산을 이해하는 데 도움이 된다. 이차 모멘트 $E[{\bf x}{\bf x}^T]$ 식을 변수 변환 후 적분하고, $z$를 고유 벡터 $u_i$로 표현하면 이차 모멘트를 계산할 수 있다.
\[E[{\bf x}{\bf x}^T]={\pmb \mu}{\pmb \mu}^T + \Sigma \qquad{(2.62)}\]공분산은 $cov[{\bf x}]=E[({\bf x}-E[{\bf x}])({\bf x}-E[{\bf x}])^T]$ 으로 정의된다. $E[x]=\mu$를 이용하여 계산하면 공분산행렬이 된다.
\[cov[{\bf x}]=\Sigma \qquad{(2.64)}\]가우시안 분포의 한계
$D$차원 데이터에 대해서 $D(D+3)/2$개의 매개변수를 갖게 된다. 일반적으로 대칭인 공분산 행렬 $\Sigma$는 $D(D+1)/2$ 개의 독립적인 매개변수를 갖는다. 여기서 $D$ 값이 커지면 행렬과 역행렬을 계산하는 것이 느려질 수 있다.
이에 대한 해결 방안 중 하나는 제한된 형태의 공분산 행렬을 사용하는 것이다.
위 그림은 공분산 행렬의 형태에 따른 이차원 가우시안 분포에서의 상수 확률 밀도의 경로를 나타낸 것이다.
- (a) : 일반적인 공분산 행렬
- (b) : 대각 행렬
- 공분한 행렬로 대각 행렬(diagonal matrix)만 사용하면 $2D$ 개의 매개변수만 고려하면 된다.
- (c) : 항등 행렬의 상수배
- 등방성(isotropic) 공분산을 사용하면 $D+1$ 개의 매개변수가 있고, 상수 밀도의 경로가 구의 형태를 띤다.
이렇게, 자유도를 줄임으로써 역행렬 계산을 더 빠르게 할 수 있다. 그러나, 확률 밀도의 형태를 제약시키며, 데이터의 상관관계를 표현하는 데 실패할 수도 있다.
또 다른 가우시안 분포의 한계점은 가우시안 분포가 단봉(unimodal)이라서 다봉(multimodal) 분포에 대한 적절한 근사치를 제공할 수가 없다는 것이다.
이를 해결하기 위해 잠재 변수(latent variable), 숨은 변수(hidden variable), 비관측 변수(unobserved variable) 을 사용할 수 있다.
혹은, 이산 잠재 변수를 통해 가우시안 혼합 모델을 사용할 수도 있다.
2.3.1 조건부 가우시안 분포
만약 두 변수 집합이 결합적으로 가우시안 분포를 보인다면, 하나의 변수 집합에 대한 다른 변수 집합의 조건부 분포 역시 가우시안 분포를 보인다는 성질이 있다. 또한, 각 변수 집합의 주변 분포 역시 가우시안 분포이다.
D차원 벡터 $\textbf x$ 가 $N({\bf x}| {\pmb \mu}, \Sigma)$의 가우시안 분포를 보일 때, $\textbf x$를 두 부분 집합 ${\textbf x}_a$와 ${\textbf x}_b$로 나누자. 그러면 $\textbf x$는 다음과 같이 표현할 수 있다.
${\bf x}=\dbinom{ {\bf x}_a }{ {\bf x}_b} \qquad{(2.65) }$
그리고 평균과 공분산 행렬은 다음과 같이 정의되겠다.
\[\mu=\dbinom{ {\pmb \mu}_a }{ {\pmb \mu}_b } \qquad{(2.66)}\] \[\Sigma=\left(\begin{array}{cc}\Sigma_{aa} & \Sigma_{ab}\\ \Sigma_{ba} & \Sigma_{bb}\end{array}\right) \qquad{(2.67)}\]공분산 행렬의 역행렬을 다음과 같이 정밀도 행렬(precision matrix)로 정의해서 사용하는 것이 편리할 수 있다.
\[\Lambda=\Sigma^{-1} \qquad{(2.68)}\] \[\Lambda=\left(\begin{array} {}\Lambda_{aa} & \Lambda_{ab} \\ \Lambda_{ba} & \Lambda_{bb} \end{array}\right) \qquad(2.69)\]- 대칭 행렬의 역행렬도 대칭
- $\Lambda_{aa}$와 $\Lambda_{bb}$는 대칭
- $\Lambda_{ab}^{\rm T}=\Lambda_{ba}$
- 단, $\Lambda_{aa}$는 $\Sigma_{aa}$의 역행렬이 아님
다변량 가우시안 분포의 지수부를 이차형식으로 전개하면 조건부 분포 $p(\textbf x_a | \textbf x_b)$가 가우시안 분포를 따른다는 것을 확인할 수 있다.
조건부 분포의 평균과 공분산을 구하기 위해, 이차형식의 계수를 분석한다. 이차형식의 계수는 공분산의 역행렬(정밀도 행렬, precision matrix)과 관련이 있으며, 이를 통해 조건부 공분산 $Σ_{a|b}$를 도출한다. 그리고, 1차 계수를 분석하여 조건부 평균 $μ_{a|b}$를 구한다.
그 다음, 계산된 정확도 행렬(역공분산)을 공분산 행렬로 변환하기 위해 역행렬 관계를 사용한다. Schur 보완(Schur complement)을 이용하여 복잡한 행렬의 역을 계산하는 방법을 적용하면
최종적으로, 조건부 분포 $p(\textbf x_a | \textbf x_b)$에 대한 평균과 분산은 다음과 같다.
\[{\pmb \mu}_{a|b} = {\pmb \mu}_a + \Sigma_{aa}\Sigma_{bb}^{-1}({\bf x}_b-{\pmb \mu}_b) \qquad{(2.81)}\] \[\Sigma_{a|b} = \Sigma_{aa} - \Sigma_{ab}\Sigma_{bb}^{-1}\Sigma_{ba} \qquad{(2.82)}\]조건부 분포 $p(\textbf x_a | \textbf x_b)$를 표현할 때 계산과정에서 분할 정밀 행렬을 사용하는 것이 분할 공분산 행렬을 사용하는 것보다 더 단순한 형태를 띤다.
조건부 분포 $p(\textbf x_a | \textbf x_b)$의 평균은 $\textbf x_b$에 대한 일차식이며, 공분산은 $\textbf x_b$에 대해 독립적이다.
이것이 선형 가우시안 모델(linear Gaussian)이다.
2.3.2 주변 가우시안 분포
결합 분포 $p(\textbf x_a, \textbf x_b)$가 가우시안 분포이면 조건부 분포 $p(\textbf x_a | \textbf x_b)$도 가우시안 분포임을 확인했다.
\[p({\bf x}_a) = \int p({\bf x}_a, {\bf x}_b)d{\bf x}_b \qquad{(2.83)}\]주변 분포(Marginal Distribution)는 결합 확률 분포(Joint Probability Distribution)에서 특정 변수들을 집중해서 보고 나머지 변수들을 전체적으로 고려하여 얻어진 분포이다. 즉, 결합 확률 분포에서 한 변수 또는 몇몇 변수들에 대한 전체 확률 분포를 알고 싶을 때, 나머지 변수들에 대해서는 관심이 없어 그 변수들을 모두 합산하거나 적분하여 제거함으로써 얻는다.
이러한 특징을 생각해보면, 식 2.83의 주변 분포 역시 가우시안 분포이다. 결합 분포 지수상의 이차식에 초점을 맞춰서 주변 분포 $p(\textbf x_a )$의 평균과 공분산을 구하면 다음과 같다.
\[E[{\bf x}_a] = {\pmb \mu}_a \qquad{(2.92)}\] \[cov[{\bf x}_a] = \Sigma_{aa} \qquad{(2.93)}\]조건부 분포에서는 분할 정밀 행렬을 사용할 때 평균과 공분산이 단순하게 표현됐었다. 주변 분포는 분할 공분산 행렬을 사용할 때 가장 단순하게 표현된다.
분할 가우시안 (Partitioned Gaussians) 정리
$\Lambda \equiv \Sigma^{-1}$이고, 결합 가우시안 분포 $N({\bf x}|{\pmb \mu}, \Sigma)$가 주어졌을 때,
\[{\bf x}=\dbinom{ {\bf x}_a }{ {\bf x}_b },\ \ \ \ {\pmb \mu}=\dbinom{ {\pmb \mu}_a }{ {\pmb \mu}_b } \qquad{(2.94)}\] \[\Sigma=\begin{pmatrix}\Sigma_{aa} & \Sigma_{ab}\\ \Sigma_{ba} & \Sigma_{bb}\end{pmatrix}\ , \ \ \ \ \Lambda = \begin{pmatrix}\Lambda_{aa} & \Lambda_{ab}\\ \Lambda_{ba} & \Lambda_{bb}\end{pmatrix} \qquad{(2.95)}\]로 표현할 수 있고, 조건부 분포는 다음과 같다.
\[p({\bf x}_a|{\bf x}_b) = N({\bf x}_a\;|\;{\pmb \mu}_{a|b}, \Lambda_{aa}^{-1}) \qquad{(2.96)}\] \[{\pmb \mu}_{a|b} = {\pmb \mu}_a - \Lambda_{aa}^{-1}\Lambda_{ab}({\bf x}_b-{\pmb \mu}_b) \qquad{(2.97)}\]주변 분포는 다음과 같다.
\[p({\bf x}_a) = N({\bf x}_a\;|\;{\pmb \mu}_a, \Sigma_{aa}) \qquad{(2.98)}\]다음은 두 개의 변수에 대한 다변량 가우시안 분포의 조건부 분포와 주변 분포의 예시이다.
왼쪽은 두 변수에 대한 가우시안 분포 $p(x_a, x_b)$의 경로를 나타낸 것이다. 오른쪽의 빨간 곡선이 조건부 분포, 파란색 곡선이 주변 분포이다.
2.3.3 가우시안 변수에 대한 베이지안 정리
지금까지 가우시안 분포 $p(\textbf x)$에 대해 벡터를 $\textbf x=(\textbf x_a, \textbf x_b)$로 나누어 조건부 분포 $p({\bf x}_a|{\bf x}_b)$와 주변 분포 $p({\bf x}_a)$가 가우시안이 된다는 것을 확인했다., 또한, 조건부 분포의 평균이 $\textbf x_b$에 대하여 선형임도 확인했다.
선형 가우시안의 예시 중 하나가 가우시안 주변 분포 $p(\textbf x)$와 가우시안 조건부 분포 $p(\textbf y|\textbf x)$의 평균이 $\textbf x$에 대해 선형이며 공분산이 $\textbf x$에 대해 독립적인 상황이다.
이러한 조건을 가정했을 때,
주변 분포 $p(\textbf y)$와 조건부 분포 $p(\textbf x|\textbf y)$를 구하고 싶다면 어떻게 해야 하는가?
베이즈 정리를 써야 한다는 감이 세게 온다.
일단, $\bf x$와 $\bf y$의 결합 확률 분포를 가우시안 분포로 설정하고, 로그를 취하여 이차형식의 형태로 전개한다.
그리고, $\bf x$와 $\bf y$에 대한 이차형식을 분리하여 정확도 행렬 $\bf R$을 도출하고, 이를 통해 공분산 행렬을 구한다.
\[{\bf R} = \left(\begin{array}{cc}\Lambda+{\bf A}^T{\bf L}{\bf A} & -{\bf A}^T{\bf L}\\-{\bf L}{\bf A} & {\bf L}\end{array}\right) \qquad{(2.104)}\] \[cov[{\bf z}]={\bf R}^{-1} = \left(\begin{array}{cc}\Lambda^{-1} & \Lambda^{-1}{\bf A}^T \\ {\bf A}\Lambda^{-1} & {\bf L}^{-1}+{\bf A}\Lambda^{-1}{\bf A}^T \end{array}\right) \qquad{(2.105)}\] \[E[{\bf z}] = \dbinom{ {\pmb \mu} }{ {\bf A} {\pmb \mu} - {\bf b}} \qquad{(2.108)}\]p(y)
이를 바탕으로 주변 확률 분포 $p(y)$를 구해보자. 결합 분포에서 $\bf x$에 대해 적분하여 $\bf y$의 주변 분포를 구한다. $\bf y$의 평균과 공분산을 계산하여, $\bf y$만을 고려한 가우시안 분포의 특성을 파악할 수 있다.
\[E[{\bf y}] = {\bf A}{\pmb \mu} + {\bf b} \qquad{(2.109)}\] \[cov[{\bf y}] = {\bf L}^T + {\bf A}\Lambda^{-1}{\bf A}^T \qquad{(2.110)}\]여기서 $\rm A=I$ 인 경우, 이 결괏값은 두 가우시안 분포의 콘볼루션(convolution)에 해당한다. 평균은 두 가우시안의 평균의 합이고, 공분산은 두 가우시안의 공분산의 합이다.
p(x∣y)
조건부 확률 분포 $p(x∣y)$의 표현식을 구해보자.
주어진 $\bf y$에 대한 $\bf x$의 조건부 평균과 공분산을 계산하면 조건부 분포에 대한 평균과 공분산은 쉽게 나온다.
\[E[{\bf x}|{\bf y}] = (\Lambda+{\bf A}^T{\bf L}{\bf A})^{-1}\{ {\bf A}^T{\bf L}({\bf y}-{\bf b})+\Lambda{\pmb \mu}\} \qquad{(2.111)}\] \[cov[{\bf x}|{\bf y}] = (\Lambda+{\bf A}^T{\bf L}{\bf A})^{-1} \qquad{(2.112)}\]주변 가우시안 분포와 조건부 가우시안 분포
$\textbf x$에 대한 주변 가우시안 분포와 $\textbf x$가 주어졌을 때의 $\textbf y$의 조건부 분포가 다음과 같이 주어졌다고 하자.
\[p({\bf x}) = N({\bf x}\;|\;{\pmb \mu}, \Lambda^{-1}) \qquad{(2.113)}\] \[p({\bf y}|{\bf x}) = N({\bf y}\;|\;{\bf A}{\bf x}+{\bf b}, {\bf L}^{-1}) \qquad{(2.114)}\]베이지안을 사용하여 다음과 같은 주변 분포와 조건부 분포를 도출해 낼 수 있다.
\[p({\bf y}) = N({\bf y}\;|\;{\bf A}{\pmb \mu}+{\bf b}, {\bf L}^{-1}+{\bf A}\Lambda^{-1}{\bf A}^T) \qquad{(2.115)}\] \[p({\bf x}|{\bf y}) = N({\bf x}\;|\;\Sigma\{ {\bf A}^T {\bf L}({\bf y}-{\bf b})+\Lambda{\pmb \mu}\}, \Sigma) \qquad{(2.116)}\]단, $\Sigma$는 다음과 같다.
\[\Sigma = (\Lambda + {\bf A}^T{\bf L}{\bf A})^{-1} \qquad{(2.117)}\]2.3.4 가우시안 분포의 최대 가능도
데이터 집합 ${\bf X}=({\bf x}_1,…,{\bf x}_n)^{\rm T}$의 관측값 {${\bf x}_n$} 들이 다변량 가우시안 분포로부터 독립적으로 추출되었다고 하자. (i.i.d)
이때, 원 분산의 매개변수들을 MLE를 이용하여 추정할 수 있다.
\[\ln p({\bf X}|{\pmb \mu}, \Sigma) = -\frac{ND}{2}\ln(2\pi) - \frac{N}{2}\ln|\Sigma|-\frac{1}{2}\sum_{n=1}^{N}({\bf x}_n-{\pmb \mu})^T\Sigma^{-1}({\bf x}_n-{\pmb \mu}) \qquad{(2.118)}\]로그 가능도 함수는 위와 같다. 이 식5을 정리하면 가능도 함수값은 다음 두 값에만 종속되어 있음을 알 수 있다.
\[\sum_{n=1}^{N}{\bf x}_n \ , \ \ \ \ \sum_{n=1}^{N}{\bf x}_n{\bf x}_n^T \qquad{(2.119)}\]위 두 값이 충분 통계량이다.
로그 가능도의 $\mu$에 대한 미분값을 다음과 같다.
\[\dfrac{\partial}{\partial {\pmb \mu}}\ln p({\bf X}|{\pmb \mu}, \Sigma) = \sum_{n=1}^{N}\Sigma^{-1}({\bf x}_n-{\pmb \mu}) \qquad{(2.120)}\]이 미분값을 0으로 놓으면 최대 가능도 추정값의 해를 구할 수 있다.
\[{\pmb \mu}_{ML} = \frac{1}{N}\sum_{n=1}^{N}{\bf x}_n \qquad{(2.121)}\]식 2.121은 관측된 데이터 포인트들의 평균값이다. 공분산에 대해서도 해를 찾으면 다음과 같다.
\[\Sigma_{ML} = \frac{1}{N}\sum_{n=1}{N}({\bf x}-{\pmb \mu}_{ML})({\bf x}-{\pmb \mu}_{ML})^T \qquad{(2.122)}\]이때의 평균과 공분산의 최대 가능도 추정치의 기댓값은 다음과 같다.
\[E[{\pmb \mu}_{ML}] = {\pmb \mu} \qquad{(2.123)}\] \[E[\Sigma_{ML}] = \frac{N-1}{N}\Sigma \qquad{(2.124)}\]평균은 MLE 추정치와 실제 평균이 같은 것을 확인할 수 있다. 그러나, 공분산은 추정치의 기댓값이 실제 공분산보다 작게 편향되어 있다.
다음과 같이 이 편향성을 수정할 수 있다.
\[\tilde{\Sigma} = \frac{1}{N-1}\sum_{n=1}^{N}({\bf x}_n - {\pmb \mu}_{ML})({\bf x}_n - {\pmb \mu}_{ML})^T \qquad{(2.125)}\]2.3.5 순차 추정
순차적인 방법론은 데이터 포인트들을 하나씩 처리하고 바로 바로 버린다. 실시간으로 계속 학습을 해야 하거나 데이터 집합의 크기가 너무 커저 다 저장하고 있을 수 없을 때 유용하다.
식 1.121은 $N$개의 관측값에 대하여 추정한 것이다. 순차 추정의 개념을 이해하기 위해 $N - 1$개의 데이터 포인트들은 이미 학습이 되었고, 마지막 $N$ 번째 데이터 포인트를 반영하는 과정을 살펴보자.
\[\begin{align} \boldsymbol{\mu}_{ML}^{(N)} &= \frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_n \notag \\ &= \frac{1}{N}\mathbf{x}_N + \frac{1}{N}\sum_{n=1}^{N-1}\mathbf{x}_n \notag \\ &= \frac{1}{N}\mathbf{x}_N + \frac{N-1}{N}\boldsymbol{\mu}_{ML}^{(N-1)} \notag \\ &= \boldsymbol{\mu}_{ML}^{(N-1)}+\frac{1}{N}(\mathbf{x}_N-\boldsymbol{\mu}_{ML}^{(N-1)}) \notag \qquad{(2.126)} \end{align}\]$N - 1$개의 데이터 포인트를 관찰한 후에 원 분포의 평균을 $\mu_{ML}^{(N-1)}$로 추정을 했다. 마지막 데이터 포인트 $\textbf x_N$ 을 관측한 후, 예전 관측값을 1/N에 비례하는 만큼 오차 방향 $ (\textbf x_N - \pmb \mu_{ML}^{(N-1)}) $ 으로 이동시켰다. 이렇게 새로운 추정값 $ \mu_{ML}^{(N)} $ 을 구했다.
로빈스 몬로(Robbins-Monro) 알고리즘
이와 같은 순차적인 학습을 일반화 한 것이 로빈스 몬로 알고리즘이다.
$\theta$가 주어졌을 때의 $z$의 조건부 기댓값으로 결정 함수 $f(\theta)$를 정의할 수 있다. 이렇게 정의되는 함수가 회귀 함수(regression function)이다.
\[f(\theta)\equiv E[z|\theta] = \int zp(z|\theta)dz \qquad{(2.127)}\]이 함수를 도식화하면 위와 같고, 여기서 회귀 함수 $f(\theta)$의 근 $\theta^\star$을 구하고 싶다.
로민스 몬로 알고리즘에서 업데이트는 다음 식에 따라 이루어진다.
\[\theta^{(N)} = \theta^{(N-1)} - a_{N-1} z(\theta^{N-1}) \qquad{(2.129)}\]계수 $a_N$의 조건
- $\lim_{N\rightarrow\infty}a_N=0$ : 시간이 지남에 따라 업데이트의 그키가 0으로 수렴하여 $\theta$가 특정 값에 수렴하도록 한다.
- $\sum_{N=1}^{\infty}a_N=\infty$ : 충분히 많은 업데이트를 통해 $\theta$가 실제 해에 도달할 수 있도록 한다.
- $\sum_{N-1}^{\infty}a_N^2<\infty$ : 업데이트 과정에서의 노이즈가 무한히 커지지 않도록 한다.
가능도 함수의 해와 회귀 함수의 근 찾기
$f(\theta) = -\log p(x|\theta)$ 의 최솟값을 찾는 것은 $f’(\theta)=0$이 되는 $\theta$ 값을 찾는 것과 같으며, 이는 $f(\theta)$의 근을 찾는 문제로 해석할 수 있다. 따라서, 최대 가능도 해를 찾는 것은 회귀 함수의 근을 찾는 것과 같다.
MLE에 로빈스 몬로 방법론을 적용하면 다음과 같다.
\[\theta^{(N)} = \theta^{(N-1)} - a_{N-1}\frac{\partial}{\partial\theta^{(N-1)}}\left[-\ln p(x_N|\theta^{(N-1)})\right] \qquad{(2.135)}\]가능도 함수를 $\mu_{ML}$로 미분하여 변수 $z$를 구할 수 있다. 데이터 포인트 $x$와 현재 추정된 평균 $\mu_{ML}$ 사이의 차이를 분산 $\sigma^2$로 정규화한 값이다. 로빈스 몬로에서 가우시안 분포에 대한 파라미터 $\mu_{ML}$의 업데이트를 할 때 사용된다.
\[z=\frac{\partial}{\partial\mu_{ML}}[-\ln p(x|\mu_{ML}, \sigma^2)]=-\frac{1}{\sigma^2}(x-\mu_{ML}) \qquad{(2.136)}\]$z$ 를 구하면 이 값을 대입하여 $E[x|\theta]=-\frac{1}{\sigma^2}(\mu-\mu_{ML})$를 구할 수 있다. ($\theta$는 모델의 파라미터)
$E[z|\theta]$는 실제로 $\mu$와 $\mu_{ML}$ 사이의 차이에 비례하며, 이 차이는 회귀선을 따라 직선 함수로 나타낼 수 있다. 즉, $\mu_{ML}$이 $\mu$에 가까워질수록 $E[z|\theta]$는 0에 가까워지고, 이는 최적의 추정치 $\mu_{ML}$가 $\mu$에 수렴함을 의미한다.
$z$의 분포는 다음 그림과 같이 $-\frac{1}{\sigma^2}(\mu-\mu_{ML})$를 평균으로 하는 가우시안 분포가 된다.
가우시안 분포에 대해 $\theta$가 평균 $\mu_{ML}$에 해당하는 경우, 회귀 함수는 빨간색 직선의 형태를 띠게 된다. 이 회귀 함수의 근은 실제 평균 $\mu$에 해당한다.
2.3.6 가우시안 분포에서의 베이지안 추론
MLE를 토대로 매개변수 $\mu$와 $\Sigma$에 대한 점 추정값을 구했다. 이제 이 매개변수들에 대한 사전 분포를 도입하여 베이지안 방법론을 사용하자.
가우시안 분포의 분산을 알고 있는 상황에서 평균을 추정
분산 $\sigma^2$을 알고 있고, N개의 관찰값 ${\bf x} = (x_1,…,x_N)^{\rm T}$이 주어졌을 때 평균 $\mu$를 추정하고 싶다.
$\mu$의 가능도 함수($\mu$가 주어졌을 때 N개의 관찰값의 확률)는 다음과 같다.
\[p({\bf x}|\mu) = \prod_{n=1}^{N}p(x_n|\mu) = \dfrac{1}{(2\pi\sigma^2)^{N/2}}\exp\left\{-\frac{1}{2\sigma^2}\sum_{n=1}^{N}(x_n-\mu)^2\right\} \qquad{(2.137)}\]가능도 함수는 $\mu$에 대한 이차식 지수 형태를 띤다. 따라서, 사전 분포 $p(\mu)$로 가우시안 분포를 선택하면 사전분포가 이 가능도 함수의 켤레 분포가 될 것이다.
\[p(\mu) = N(\mu|\mu_0, \sigma_0^2) \qquad{(2.138)}\]사전 분포를 위와 같이 설정하면 사후 분포는 다음과 같다.
\[p(u|{\bf x}) \propto p({\bf x}|\mu)p(\mu) \qquad{(2.139)}\]이 식을 조금 조작하면 사후 분포는 다음과 같다.
\[p(\mu|{\bf x}) = N(\mu|\mu_N, \sigma_N^2) \qquad{(2.140)}\]- $\mu_N = \frac{\sigma^2}{N\sigma_0^2+\sigma^2}\mu_0 + \frac{N\sigma_0^2}{N\sigma_0^2+\sigma^2}\mu_{ML} \qquad{(2.141)}$
- 사후 분포의 평균이 사전 평균 $\mu_0$와 최대 가능도 해 $\mu_{ML}$의 절충값에 해당한다.
- $\frac{1}{\sigma_N^2}=\frac{1}{\sigma_0^2}+\frac{N}{\sigma^2} \qquad{(2.142)}$
- 분산의 사후 분포식
$\mu_{ML}$은 $\mu$의 최대 가능도 해이자, 표본의 평균으로 주어진다.
\[\mu_{ML} = \frac{1}{N}\sum_{n=1}^{N}x_n \qquad{(2.143)}\]사후 분산의 정밀도는 사전 분산의 정밀도에 관측된 데이터 하나당 하나씩의 관측 데이터 정밀도를 합한 것과 같다. 따라서, 데이터 포인트의 수가 증가하면 정밀도는 꾸준히 증가한다. = 사후 분포의 분산이 꾸준히 감소한다.
분산을 알고 있다는 가정하에 가우시안 분포의 평균 $\mu$에 대한 베이지안 추론을 나타낸 것이다. $N$이 커질 수록 분산은 작아진다.
앞에서 $N - 1$개의 데이터 포인트들을 관찰한 시점에서 마지막 N 번째 데이터 포인트의 기여도를 살펴봤다. 여기서도 마지막 포인트를 떼어 놓으면 사후 분포는 다음과 같이 나타낼 수 있다.
\[p(\mu|{\bf x}) \propto \left[p(\mu)\prod_{n=1}^{N-1}p(x_n|\mu)\right]p(x_N|\mu) \qquad{(2.144)}\]가우시안 분포의 평균을 알고 있는 상황에서 분산을 추정
위에서는 평균을 추정하고 싶어서 $\mu$에 대한 가능도 함수에서 시작했었다. 지금은 분산을 알고 싶다. 분산 대신 정밀도로 계산하는 것이 더 편리하다. $\lambda$의 가능도 함수는 다음과 같다.
\[p({\bf x}|\lambda) = \prod_{n=1}^{N} N(x_n|\mu, \lambda^{-1}) \propto \lambda^{N/2} \exp\left\{ -\frac{\lambda}{2} \sum_{n=1}^{N}(x_n-\mu)^2\right\} \qquad{(2.145)}\]따라서 켤레 사전 분포는 $\lambda$의 거듭제곱과 $\lambda$의 선형 함수의 지수 함수를 곱한 것에 비례하는 형태를 띠어야 한다. 감마 분포 (gamma distribution)은 다음과 같다.
\[Gam(\lambda|a,b)=\frac{1}{\Gamma(a)}b^a\lambda^{a-1}\exp(-b\lambda) \qquad{(2.146)}\]- $\Gamma(x) = \int_{0}^{\infty}u^{x-1}e^{-u}du$
감마 분포의 평균과 분산은 다음과 같다.
\[E[\lambda] = \frac{a}{b} \qquad{(2.147)}\] \[var[\lambda] = \frac{a}{b^2} \qquad{(2.148)}\]사전 분포 $Gam(\lambda|a_0,b_0)$ 에 가능도 함수 식을 곱하면 사후 분포를 얻을 수 있다.
\[p({\bf x}|\lambda) \propto \lambda^{a_0-1}\lambda^{N/2} \exp \left\{-b_0\lambda-\frac{\lambda}{2}\sum_{n=1}^{N}(x_n-\mu)^2\right\} \qquad{(2.149)}\]처음에 감마 분포를 켤레 사전 분포로 사용했던 것이기 때문에, 사후 분포 역시 구조가 똑같다. 여기서 $a$와 $b$는 다음과 같다.
\[a_N = a_0 + \frac{N}{2} \qquad{(2.150)}\] \[b_N = b_0 + \frac{1}{2}\sum_{n=1}^{N}(x_n-\mu)^2 = b_0 + \frac{N}{2}\sigma_{ML}^2 \qquad{(2.151)}\]식 2.150을 보면 $N$개의 데이터 포인트를 관측하면 $a$ 값을 $N/2$ 만큼 증가시킨다. 따라서, 사전 분포의 매개변수 $a_0$은 $2a_0$개의 유효 사전 관측값에 해당한다.
식 2.151에서는 $N$개의 데이터 포인트들이 $b$ 값에 $N\sigma^2_{ML}/2$ 만큼 기여하고 있다. 따라서, $b_0$은 $2b_0/2a_0=b_0/a_0$의 분산을 갖는 $2a_0$개의 유효 사전 관측값에 해당한다.
정밀도 대신에 분산을 이용할 때는 역감마 분포 (inverse gamma distribution)을 사용한다.
가우시안 분포의 평균과 정밀도를 둘 다 모를 때
일단 켤레 사전 분포를 찾긴 해야 하니 가능도 함수에서 $\mu$와 $\lambda$에 대한 의존도를 살펴보자.
\[\begin{align} p({\bf x}|\mu, \lambda) &= \prod_{n=1}^{N}\left(\frac{\lambda}{2\pi}\right)^{1/2} \exp\left\{-\frac{\lambda}{2}(x_n-\mu)^2\right\} \notag \\ &\propto \left[\lambda^{1/2}\exp\left(-\frac{\lambda\mu^2}{2}\right)\right]^N\exp\left\{\lambda\mu\sum_{n=1}^{N}x_n-\frac{\lambda}{2}\sum_{n=1}^{N}x_n^2\right\} \notag \qquad{(2.152)} \end{align}\]가능도 함수와 같은 의존도 구조를 갖는 사전 분포는 다음과 같이 설정할 수 있다. 상수 $c$, $d$, $\beta$ 를 좀 도입해야 한다.
\[\begin{align} p(\mu, \lambda) &\propto \left[\lambda^{1/2}\exp\left(-\frac{\lambda\mu^2}{2}\right)\right]^\beta\exp\{c\lambda\mu-d\lambda\} \notag \\ &= \exp\left\{-\frac{\beta\lambda}{2}(\mu-c/\beta)^2\right\}\lambda^{\beta/2}\exp\left\{-(d-\frac{c^2}{2\beta})\lambda\right\} \notag \qquad{(2.153)} \end{align}\]$p(\mu, \lambda) = p(\mu|\lambda)p(\lambda)$이므로 다음과 같은 정규화된 사전 분포를 이끌어낼 수 있다.
\[p(\mu, \lambda) = N(\mu|\mu_0, (\beta\lambda)^{-1})Gam(\lambda|a,b) \qquad{(2.154)}\]- $\mu_0 = \frac{c}{\beta}$
- $a=(1+\beta)/2$
- $b=d-\frac{c^2}{2\beta}$
식 2.154를 정규 감마 (normal gamma), 가우시안 감마 (Gaussian gamma) 분포라고 한다. 이 분포는 단순히 $\mu$와 $\lambda$의 사전 분포의 곱이 아니다. $\mu$와 $\lambda$에 대해서 독립적인 사전 분포를 선택해도 사후 분포에서는 연결성이 있다. 아래는 정규 감마 분포를 도식화 한 것이다.
다변량 가우시안 분포의 추정
정밀도가 알려져 있을 때
$D$차원 변수 $\bf x$에 대한 다변량 가우시안 분포 $N({\bf x}|{\pmb \mu}, \Lambda^{-1})$의 경우, 정밀도가 알려져 있으면 평균 $\mu$에 대한 켤레 사전 분포는 다시 가우시안 분포이다.
평균이 알려져 있을 때
켤레 사전 분포는 위샤트 분포 (Wishart distribution)이며, 다음과 같다.
\[W(\Lambda|{\bf W}, v) = B|\Lambda|^{(v-D-1/2)}\exp\left(-\frac{1}{2}Tr({\bf W}^{-1}\Lambda)\right) \qquad{(2.155)}\]- $v$ : 분포의 자유도 (degree of freedom)
- W : $D \times D$ 척도 행렬
- $Tr()$ : 행렬의 대각합
- $B({\bf W}, v) = |{\bf W}|^{-v/2}\left(2^{vD/2}\pi^{D(D-1)/4}\prod_{i=1}^{D}\left(\frac{v+1-i}{2}\right)\right)^{-1} \qquad{(2.156)}$
정밀도 대신 공분산 행렬을 이용하여 정의할 때는 켤레 사준 분포는 역 위샤트 분포 (inverse Wishart distribution)이다.
평균과 정밀도를 둘 다 모를 때
단변량과 비슷한 추론 과정을 통해 다음의 사전 분포를 구할 수 있다.
\[p({\pmb \mu}, \Lambda|\mu_0, \beta, {\bf W}, v) = N({\pmb \mu}|{\pmb \mu}_0, (\beta\Lambda)^{-1})W(\Lambda|{\bf W}, v) \qquad(2.157)\]이 식을 정규 위샤트 (normal Wishart), 가우시안 위샤트 (Gaussian Wishart)라고 한다.
2.3.7 스튜던트 t 분포
사전 분포로 감마 사전 분포 $Gam(\tau|a,b)$를 갖는 단변량 가우시안 분포 $N(x|\mu, \tau^{-1})$가 주어졌을 때, 정밀도를 적분해서 없애면 $x$에 대한 주변 분포를 구할 수 있다.
\[\begin{align} p(x|\mu,a,b) &= \int_{0}^{\infty}N(x|\mu, \tau^{-1})Gam(\tau|a,b)d\tau \notag \\ &= \int_{0}^{\infty}\frac{b^a e^{(-b\tau)}\tau^{(a-1)}}{\Gamma(a)}\left(\frac{\tau}{2\pi}\right)^{1/2}\exp\left\{-\frac{\tau}{2}(x-\mu)^2\right\}d\tau \notag \\ &= \frac{b^a}{\Gamma(a)}\left(\frac{1}{2\pi}\right)^{1/2}\left[b+\frac{(x-\mu)^2}{2}\right]^{-a-1/2}\Gamma(a+1/2) \qquad{(2.158)} \notag \end{align}\]$z = \tau[b+(x-\mu)^2/2]$로 변수 변환을 하면 위 식을 다음과 같이 나타낼 수 있다.
\[St(x|\mu,\lambda,v) = \frac{\Gamma(v/2+1/2)}{\Gamma(v/2)}\left(\frac{\lambda}{\pi v}\right)^{1/2}\left[1+\frac{\lambda(x-\mu)^2}{v}\right]^{-v/2-1/2} \qquad{(2.159)}\]- $v$ : 자유도 (degree of freedom)
- $\lambda$ : t 분포의 정밀도 (분산의 역은 아니다)
이 분포를 스튜던트 t 분포 (student’s t-distribution)라고 한다.
그림에서 볼 수 있듯이, $v \rightarrow \infty$ 일 때 분포는 가우시안 분포가 되고, $v=1$일 때는 코시 분포(cauchy distribution)이 된다.
스튜던트 t 분포는 같은 평균과 다른 정밀도를 가진 무한히 많은 가우시안 분포들을 합산함으로써 구할 수 있다. 무한한 숫자의 가우시안 분포가 혼합된 것이다.
t 분포의 최대 가능도 해는 EM(expectation maximization, 기댓값 최대화) 알고리즘을 통해 찾을 수 있다.
그림에서도 볼 수 있듯이, $v$ 값이 충분히 크지 못할 때는 이상값(outlier, 특이값)에 해당하는 데이터 포인트들의 존재에 대해 덜 예민하다. 이 성질을 강건성(robustness)라고 한다. 가우시안보다 더 긴 꼬리를 가지게 되는 것이다.
막대 그래프는 가우시안 분포로부터 추출한 데이터 포인트들이다. (a)에서 빨간 곡선은 최대 가능도를 이용하여 근사한 t 분포이다. 녹색 곡선은 그 뒤에 겹쳐서 거의 보이지 않는다. 두 분포가 비슷한 해를 보인다는 것이다. (b)에서 이상값을 3개 추가했더니 가우시안 분포는 이상값의 영향을 받아서 왜곡되었지만, t 분포는 영향을 덜 받은 것을 확인할 수 있다. 이상값에 덜 예민한 것이다.
최소 제곱법 (Least Squares Method)은 회귀 문제를 푸는 데 널리 사용되는 방법이다. 잔차의 제곱합을 최소화하여, 데이터를 가장 잘 설명하는 선형 회귀선을 찾는다. 그런데, 이 최소 제곱법은 가우시안 분포를 가정한다. 오차가 평균 0과 일정한 분산을 가지며, 모든 데이터 포인트에 대해 동일하게 적용된다는 것이다. 또한, 이상치에 민감하다. t 분포와 같은 Robust Regression Model (강건한 회귀 모델)은 이런 이상치의 영향을 줄이기 위해 설계되어, 더 넓은 범위의 오차 분포를 허용한다. 이러한 분포를 바탕으로 회귀 모델을 수립하면 더 강건한 모델을 얻을 수 있다.
다변량 스튜턴트 t 분포
식 2.158에 $v=2a$, $\lambda=a/b$, $\eta=\tau b/a$ 와 같이 매개변수를 설정하고 대입하면 t 분포를 다음과 같이 적을 수 있다.
\[St(x|\mu,\lambda,v)=\int_{0}^{\infty}N(x|\mu,(\eta\lambda)^{-1})Gam(\eta|v/2,v/2)d\eta \qquad{(2.160)}\]이를 다변량 가우시안 분포 $N(\textbf x | \pmb \mu, \pmb \Lambda, v)$에 대해 일반화하면, 다변량 스튜턴트 t 분포는 다음과 같다.
\[St({\bf x}|{\pmb \mu}, \pmb \Lambda, v) = \int_{0}^{\infty}N({\bf x}|{\pmb \mu}, (\eta\pmb \Lambda)^{-1})Gam(\eta|v/2,v/2)d\eta \qquad{(2.161)}\]이 적분식을 계산하면 다음과 같다.
\[St({\bf x}|{\pmb \mu}, \pmb\Lambda, v) = \frac{\Gamma(D/2+v/2)}{\Gamma(v/2)}\frac{|\pmb\Lambda|^{1/2}}{(\pi v)^{D/2}}\left[1+\frac{\nabla^2}{v}\right]^{-D/2-v/2} \qquad{(2.162)}\]- $D$ : $\bf x$의 차원수
- $\nabla^2 = ({\bf x}-{\bf \mu})^T\Lambda({\bf x}-{\bf \mu})$ = 마할라노비스 거리의 제곱
다변량 스튜던트 t 분포의 형태이며, 다음과 같은 성질을 만족한다. 이 성질들은 단변량 형태에도 만족된다.
\[\begin{array}{lc} E[{\bf x}] = {\bf \mu}, & \;if\;v>1 & \qquad{(2.164)}\\ cov[{\bf x}]=\frac{v}{(v-2)}\Lambda^{-1}, & \;if\;v>2 & \qquad{(2.165)}\\ mode[{\bf x}]={\bf \mu} & & \qquad{(2.166)} \end{array}\]2.3.8 주기적 변수
풍향이나 시간(달력)과 같은 주기적(periodic) 변수는 가우시안 분포를 연속 변수의 밀도 모델로 사용하는 것이 적합하지 않은 경우가 종종 발생한다.
이런 주기적 변수들은 극좌표 $0 \leq \theta \leq 2\pi$를 이용하여 나타낼 수 있다. 그러나, 어떤 방향 하나를 원점으로 삼고 가우시안 분포 같은 종래의 분포를 적용하면, 어느 값을 원점으로 잡는지에 대해 종속적인 결과를 보인다.
좌표축의 선택에 종속적이지 않도록 평균값을 측정하기 위해 관측값을 단위 원 위의 포인트로 잡으면 다음 그림과 같이 크기가 1인 2차원 벡터들로 표현될 수 있다.
이 벡터들의 평균은 다음과 같다.
\[\bar{ {\bf x}} = \frac{1}{N}\sum_{n=1}^{N}{\bf x}_n \qquad{(2.167)}\]이 평균값에 해당하는 각도 $\bar \theta$를 구하면 된다. 관측값은 $\textbf x_n=(\cos \theta_n, \sin \theta_n)$으로 표현하고, 표본 평균은 $\bar {\textbf x}_n=({\bar r}\cos \bar \theta_n, {\bar r}\sin \bar \theta_n)$으로 표현할 수 있다. 이를 위 식에 대입하고 $x_1$과 $x_2$ 성분이 같다고 놓으면 다음과 같이 얻을 수 있다.
\[\bar{ {\bf x}}_1=\bar{r}\cos\bar{\theta}=\frac{1}{N}\sum_{n=1}^{N}\cos\theta_n \\\bar{ {\bf x}}_2=\bar{r}\sin\bar{\theta}=\frac{1}{N}\sum*{n=1}^{N}\sin\theta_n \qquad{(2.168)}\]$\tan \theta=\sin \theta /\cos \theta$ 이므로, 둘을 나누면 $\bar \theta$를 구할 수 있다.
\[\bar{\theta}=\tan^{-1}\left\{\frac{\sum_n\sin\theta_n}{\sum_n\cos\theta_n}\right\} \qquad{(2.169)}\]폰 미제스 분포 (von Mises distribution)
폰 미제스 분포는 가우시안 분포를 주기적 변수에 적용할 수 있도록 일반화한 것이다.
$2\pi$를 주기고 갖는 분포 $p(\theta)$는 다음의 조건을 만족시킨다.
- $p(\theta) \ge 0$ : 확률 분포는 0보다 크거나 같다.
- $\int_{0}^{2\pi}p(\theta)d\theta=1$ : 적분하면 1
- $p(\theta+2\pi)=p(\theta)$ : 주기적이다.
두 개의 변수 $\textbf x=(x_1, x_2)$에 대한 가우시안 분포를 고려해 보자. 평균 $\pmb \mu=(\mu_1, \mu_2)$ 이며, 공분산 행렬 $\pmb\Sigma=\sigma^2 \rm\pmb I$ 를 가진다 ( $\rm\pmb I$는 2x2 단위 행렬(identity matrix)/항등 행렬 ). 이 분포는 다음과 같이 적을 수 있다.
\[p(x_1,x_2)=\frac{1}{2\pi\sigma^2}\exp\left\{-\frac{(x_1-\mu_1)^2+(x_2-\mu_2)^2}{2\sigma^2}\right\} \qquad{(2.173)}\]이 분포는 그림에서 파란색으로 나타내어져 있다. $p(\textbf x)$ 가 상수일 때의 컨투어(윤곽선)은 고정된 반지름의 원주 위에 존재한다. 이 분포는 빨간색으로 그려진 단위 윈에 대해서 조건부이다.
이 주기적 분포를 극좌표호 변환하면 다음과 같다.
\[x_1=r\cos\theta\ , \ \ \ \ x_2=r\sin\theta \qquad{(2.174)}\] \[\mu_1=r_0\cos\theta_0 \ ,\ \ \ \ \mu_2=r_0\sin\theta_0 \qquad{(2.175)}\]이 변환값들을 이차원 가우시안 분포식(식 2.173)에 대입하고 단위 원 조건(r = 1)을 추가하면 가우시안 분포의 지수부는 다음과 같다.
\[\begin{align} -\frac{1}{2\sigma^2}\{(r\cos\theta-r_0\cos\theta_0)^2+(r\sin\theta-r_0\sin\theta_0)^2\} &\notag \\ =-\frac{1}{2\sigma^2}\{1+r_0^2-2r_0\cos\theta\cos\theta_0-2r_0\sin\theta\sin\theta_0\} &\notag \\ =\frac{r_0}{\sigma^2}(\theta-\theta_0)+\text{const} \qquad{(2.176)} &\notag \end{align}\]- const : $\theta$로부터 독립적인 항들의 집합
$m=r_0/\sigma^2$ 으로 두면 단위 원을 가정하면 다음과 같다.
\[p(\theta|\theta_0, m)=\frac{1}{2\pi I_0(m)}\exp\{m\cos(\theta-\theta_0)\} \qquad{(2.179)}\]- $\theta_0$ : 분포의 평균
- $m$ : 가우시안 분포의 정밀도. 집중 매개변수 (concentration parameter)
- $I_0(m)=\frac{1}{2\pi}\int_{0}^{2\pi}\exp{m\cos\theta}d\theta$
- 0차 1종 베젤 함수 (Bessel function)
- 위 식의 정규화 계수
위 식을 폰 미제스 분포(von Mises distribution), 원형 정규 분포 (circular normal)라고 한다.
왼쪽 그림에서 볼 수 있듯이, m 값이 커지면 가우시안 분포에 가까워진다. 오른쪽은 왼쪽 데카르트 좌표계를 극좌표계에 그린 것이다.
$I_0(m)$은 이렇게 생겼다고 한다.
폰 미제스 분포의 MLE
폰 미제스 분포의 매개변수 $\theta_0$과 $m$의 최대 가능도 추정값(MLE)를 구해보자.
로그 가능도 함수는 다음과 같다.
\[\ln p(D|\theta_0, m) = -N\ln(2\pi)-N\ln I_0(m) + m\sum_{n=1}^{N}\cos(\theta-\theta_0) \qquad{(2.181)}\]$\theta_0$에 대해 미분한 값을 0으로 놓으면 다음과 같다.
\[\sum_{n=1}^{N}\sin(\theta_n-\theta_0) = 0 \qquad{(2.182)}\]$\sin(A-B)=\cos{B}\cos{A}-\cos{A}\sin{B}$ 이므로 다음과 같이 쓸 수 있다.
\[\theta_0^{ML} = \tan^{-1}\left\{\frac{\sum_n\sin\theta_n}{\sum_n\cos\theta_n}\right\} \qquad{(2.184)}\]이 값이 관찰값들의 평균과 같다.
이제 $m$의 MLE를 구해보자. 로그 가능도 함수를 미분하고 $I_0’(m)=I_1(m)$로 놓고, $A(m)=\frac{I_1(m)}{I_0(m)}$ 로 정의하면 다음과 같이 적을 수 있다.
\[A(m_{ML})=\frac{1}{N}\sum_{n=1}^{N}\cos(\theta_n-\theta_0^{ML}) \qquad{(2.185)}\]코사인 합차공식을 이용하면 아래와 같이 적을 수 있다.
\[A(m_{ML})=\left(\frac{1}{N}\sum_{n=1}^{N}\cos\theta_0^{ML}\right)+\left(\frac{1}{N}\sum_{n=1}^{N}\sin\theta_n\right)\sin\theta_0^{ML} \qquad{(2.187)}\]여기서 우변을 계산하고, $A(m)=\frac{I_1(m)}{I_0(m)}$ 에 따라서 $m$을 계산하면 그 값이 $m_{ML}$이다.
폰 미제스 분포의 한계점은 단봉 분포라는 것이다. 따라서, 단봉 분포라는 한계점을 갖고 있는 가우시안과 같이 다봉 데이터를 잘 설명하지 못한다. 그래도 주기적인 특성을 가진 데이터는 모델링 잘 한다. ^^
해결책으로, 폰 미제스 분포들을 혼합해서 주기적 변수들의 다봉성을 다룰 수 있다.
2.3.9 가우시안 분포의 혼합
이전에 말했다시피, 가우시안은 단봉 분포라서 실제 데이터를 모델링하는 데에 한계가 있다.
따라서, 우리는 여러 개의 가우시안 분포를 혼합하여 혼합 분포 (mixture distribution)을 사용한다.
위 그림은 미국 옐로스톤 국립 공원의 간헐 온천 분화 측정치이다. 왼쪽은 이 데이터를 한 개의 가우시안 분포로 근사한 것이고, 오른쪽은 두 개를 선형 결합하여 근사한 것이다.
이 그림에서 파란색이 세 개의 가우시안 분포를 나타내고 있고, 빨간 색은 이 세 분포의 합이다. 이렇게 사용하면 복잡한 데이터에 더 잘 근사할 수 있겠다.
이를 식으로 이해하면 다음과 같다. $K$개의 사우시안 밀도의 중첩이다.
\[p({\bf x}) = \sum_{k=1}^{K}\pi_k N({\bf x} | {\pmb \mu}_k, \Sigma_k) \qquad{(2.188)}\]- $N({\bf x} | {\pmb \mu}_k, \Sigma_k)$ : 각각의 가우시안 밀도. 혼합의 성분(component)
- ${\pmb \mu}_k$ : 평균
- $\Sigma_k$ : 공분산
- $\pi_k$ : 혼합 계수 (mixing coefficient)
이를 가우시안 혼합 분포 (mixture of Gaussians)라고 부흔다.
3개의 성분을 가진 분포이다.(a)는 각 성분의 상수 밀도의 경로이고, (b)는 혼합 분포의 주변 확률 밀도 $p(\textbf x)$의 경로이고, (c)는 $p(\textbf x)$의 표면이다.
매개변수 $\pi_k$ 는 각 가우시안 분포가 전체 모델에서 차지하는 비중/가중치를 나타낸다. 해당 가우시안 분포에 데이터 포인트가 속할 확률을 나타낸다고 보면 되겠다. $\pi_k$ 를 조정함으로써 모델에서 특정 성분의 영향력을 증가시키거나 감소시킬 수 있다. $\pi_k$ 는 아래의 조건을 만족시켜야 한다.
\[\sum_{k=1}^{K}\pi_k=1 \qquad{(2.189)}\\ 0 \le \pi_k \le 1 \qquad{(2.190)}\]따라서 밀도 $N({\bf x} | {\pmb \mu}_k, \Sigma_k)=p(\textbf x|k)$는 $k$가 주어졌을 때의 $\bf x$의 확률이다. 아래와 같이 써도 똑같다.
\[p({\bf x}) = \sum_{k=1}^{K}p(k)p({\bf x}|k) \qquad{(2.191)}\]여기서 사후확률 $p(\textbf x|k)$를 책임값 (responsibilities)이라고 한다. 이름은 직관적으로 이해가 된다. 베이지안 정리에 따라서 이 사후 확률을 다음과 같이 표현할 수 있다.
\[\gamma_k({\bf x})\equiv p(k|{\bf x})=\frac{p(k)p({\bf x}|k)}{\sum_l p(l)p({\bf x}|l)}=\frac{\pi_kN({\bf x}|{\pmb \mu}_k, \pmb \Sigma_k)}{\sum_l \pi_l N({\bf x}|{\pmb \mu}_l, \pmb \Sigma_l)} \qquad{(2.192)}\]이 가우시안 혼합 분포의 형태는 매개변수 $\pmb \pi, \pmb \mu, \pmb \Sigma$에 의해 결정된다. 이 매개변수들은 MLE를 써서 찾을 수 있다. 아래는 가우시안 혼합 모델의 로그 가능도 함수이다.
\[\ln p({\bf X}|{\pmb \pi}, {\pmb \mu}, \pmb\Sigma) = \sum_{n=1}^{N}\ln\left\{\sum_{k=1}^{K} \pi_k N({\bf x}_n|{\pmb \mu}_k,\pmb \Sigma_k) \right\} \qquad{(2.193)}\]주어진 매개변수 $\pmb \pi, \pmb \mu, \pmb \Sigma$ 하에서 관측된 데이터 $\bf X$가 나타날 확률을 나타낸다. 이 확률을 최대화할 때의 매개변수 값을 찾고 싶은 것이다. 하지만, 로그 안에 $k$에 대한 합산이 포함되어 있어 해가 결정된 형태로 나오지 않는다. 해를 구할 때는 반복적인 최적화가 필요하다.
기댓값 최대화 (expectation maximization) 알고리즘을 통해 로그 가능도를 최대화 하는 매개변수 값을 찾는 과정을 로그 가능도의 증가가 더 이상 유의미하지 않을 때까지 (어떠한 값에 수렴할 때까지) 반복하여 데이터에 대한 모델의 적합도를 점진적으로 개선해야 한다.