/////
🏨

Chapter 9. Mixture Models and EM

Created
2021/01/28 09:41
어떠한 observed variable과 latent variable에 대한 joint distribution이 주어졌을 때, observed variable에 대한 distribution은 joint distribution을 latent variable에 대해 marginalization함으로써 얻어낼 수 있다. 이는 곧 관찰하게 될 것처럼 observed variable에 대한 marginal distribution을 joint distribution으로 풀어내 보다 간단하게 표현할 수 있도록 한다.
그리하여 이번 장에서는 mixture model을 latent variable의 측면에서 해석해 보고, 이러한 latent variable을 도입하여 maximum likelihood maximization 문제를 보다 쉽게 풀 수 있도록 하는 EM 알고리즘에 대해 배운다.

9.1. K-means Clustering

먼저 mixture model을 통해 data clustering을 할 수 있다. 그리고 clustering 알고리즘 중 흔히들 아는 K-means clustering 알고리즘이 이러한 EM 알고리즘과 상당히 밀접한 관계를 가진다. 그런 측면에서 이번 섹션에서는 간단히 K-means clustering을 짚고 넘어간다.
KK 개의 cluster가 존재한다고 하자. 또 각 cluster kk 를 대표하는 prototype인 μk\mu_k 가 존재한다고 하자. 그러면 우리의 목적은 주어진 데이터 xix_i 가 속하는 cluster과 cluster prototype μk\mu_k 를 찾는 것이다.
이제 nn 번째 데이터 xnx_nkk 번째 cluster에 속함을 나타내는 binary indicator variable rnkr_{nk} 를 정의하자. 그러면 K-means clustering 알고리즘의 objective function(혹은 distortion measure)를 다음과 같이 정의할 수 있다.
J=n=1Nk=1Krnkxnμk2J = \sum_{n=1}^N\sum_{k=1}^K r_{nk} ||x_n - \mu_k||^2
위에서 언급했듯 우리의 목적은 JJ 를 최소화하는 {rnk}\{r_{nk}\} 의 값과 {μk}\{\mu_k\} 의 값을 찾는 것인데, 이를 각각 iterative하게 최적화하는 것이 바로 K-means 알고리즘이다. 즉 {μk}\{\mu_k\} 의 초깃값을 설정하고, 이를 고정하여 {rnk}\{r_{nk}\} 를 최적화하고, 다시 이를 고정하고 {μk}\{\mu_k\} 를 최적화하는 것을 JJ 가 수렴할 때까지 반복하는 것이다. 각 데이터들의 assignment를 결정하는 단계를 E(xpectation) step이라고 부르고, 주어진 expectation에 대해 목적 함수를 최적화하는 파라미터 값({μk}\{\mu_k\})를 찾는 단계를 M(aximization) step이라고 부른다.
이 때 nn 을 포함하는 term들은 독립적이므로 각각의 nn 에 대한 rnkr_{nk} 를 독립적으로 최소화할 수 있고, 이에 따른 xnx_n 의 assignment는 μk\mu_kxnx_n 사이의 거리를 최소화하는 kk 이다. 또한 JJμk\mu_k 에 대해 미분한 후 0으로 만들면 각각의 μk\mu_k 가 해당 cluster에 속하는 데이터들의 평균값임을 알 수 있다. 이 때 JJ 는 매 iteration마다 감소하기 때문에 convergence가 보장된다.
한편 위에서는 dissimilarity를 표현하기 위해 Euclidean distance를 사용했는데, 이러한 dissimilarity를 따로 정의한 K-means 알고리즘의 확장판을 K-medoids 알고리즘이라고 부른다.

9.1.1 Image segmentation and compression

Image compression 태스크를 수행하기 위해 K-means 알고리즘을 사용하는 예시를 보여준다. 각각의 cluster prototype을 RGB 값으로 두고, 각각의 pixel을 clustering하기 위한 K-means 알고리즘을 수행하면 된다. Prototype의 값과 각 pixel의 assignment만 전송함으로써 전송에 필요한 정보량을 크게 줄일 수 있다. 다만 이는 lossless data compression이 아닌 lossy data compression임에 유의해야 한다.

9.2. Mixtures of Gaussians

이제 Gaussian mixture model을 latent variable의 측면에서 재조명해본다. 그리고 이를 이용한 EM 알고리즘을 소개한다.
먼저 Gaussian mixture을 다음과 같은 여러 Gaussian들의 linear combination으로 나타낼 수 있다.
p(x)=k=1KπkN(xμk,Σk)p(x) = \sum_{k=1}^K\pi_kN(x|\mu_k, \Sigma_k)
여기서 KK-dimensional binary random variable zz 를 도입한다. 즉 각 데이터마다 KK 개의 state를 가지는 one-hot-encoded variable이다.
이제 joint distribution p(x,z)p(x, z)p(z)p(z)p(xz)p(x | z) 로 표현해보고자 하는데, marginal distribution p(z)p(z) 를 다음과 같이 정의한다.
p(zk=1)=πkp(z_k = 1) = \pi_k
그러면 p(xz)p(x | z) 이 다음과 같고,
p(xzk=1)=N(xμk,Σk)p(x | z_k = 1) = N(x | \mu_k, \Sigma_k)
Joint distribution은 다음과 같이 표현된다.
p(x)=zp(z)p(xz)=k=1KπkN(xμk,Σk)p(x) = \sum_zp(z)p(x|z) = \sum_{k=1}^K\pi_kN(x|\mu_k,\Sigma_k)
이제 latent variable을 이용하여 Gaussian mixture model을 joint distribution으로 표현해냈다. 이로써 우리는 likelihood maximization을 (곧 배울 EM 알고리즘을 통해) 훨씬 간단하게 풀어낼 수 있게 된다.
그 전에 responsibility라는 개념을 알아두자. Responsibility는 zz 에 대한 posterior probability로 정의되며, component kk 가 데이터를 설명하는 데 가지는 responsibility를 의미한다. π\pi 는 posterior probability의 역할을 한다.
γ(zk)=p(zk=1x)=p(zk=1)p(xzk=1)j=1Kp(zj=1)p(xzj=1)=πkN(xμk,Σk)j=1KπjN(xμj,Σj)\gamma(z_k) = p(z_k=1|x) \\{}\\ = \frac{p(z_k=1)p(x|z_k=1)}{\sum_{j=1}^K p(z_j=1)p(x|z_j=1)} \\{}\\ =\frac{\pi_kN(x|\mu_k,\Sigma_k)}{\sum_{j=1}^K\pi_j N(x|\mu_j,\Sigma_j)}

9.2.1 Maximum likelihood

그러면 이제 우리에게는 데이터 {x1,...,xN}\{x_1, ..., x_N\} 이 주어졌고, 이 데이터의 분포를 Gaussian mixture로 모델링한다고 해보자. 먼저 log likelihood function은 다음과 같다.
lnp(Xπ,μ,Σ)=n=1Nln{k=1KπkN(xnμk,Σk)}\ln p(X|\pi,\mu,\Sigma) = \sum_{n=1}^N \ln \{ \sum_{k=1}^K \pi_k N(x_n |\mu_k,\Sigma_k) \}
이 때 mixture Gaussian의 maximum likelihood framework는 singularity라는 아주 중요한 문제점을 가지고 있다는 점을 짚고 넘어간다. 요는 mixture 중 하나의 Gaussian이 어떤 하나의 data point에 collapse되어 버리면 likelihood가 무한대가 되어 버린다는 것인데, 자세한 내용은 교재를 참고하자.
또 하나의 문제점은 바로 Gaussian들의 순서만 바꿈으로써 K!K! 개의 동일한 distribution을 생성할 수 있다는 identifiability라는 문제인데, 실제로 적용하는 데 별 문제는 없다고 한다. 역시 자세한 내용은 교재를 참고하자.

9.2.2 EM for Gaussian Mixtures

계속 언급했듯 latent variable을 포함하는 모델의 maximum likelihood solution을 찾을 수 있는 엘레강트하고 파워풀한(;) 알고리즘이 바로 EM(expectation-maximization) 알고리즘이다.
우선 바로 위에 구했던 XX 에 대한 marginal distribution의 likelihood를 μ\mu 에 대해 미분한 후 0으로 보내면 다음과 같은 식이 나오는데,
0=n=1NπkN(xnμk,Σk)jπjN(xnμj,Σj)Σk(xnμk)0 = -\sum_{n=1}^N\frac{\pi_kN(x_n|\mu_k,\Sigma_k)}{\sum_{j}\pi_jN(x_n|\mu_j,\Sigma_j)}\Sigma_k(x_n - \mu_k)
오른쪽 식에 곱해진 저 분수 꼴의 term이 위에서 언급했던 responsibility임을 알 수 있다.
이제 cluster kk 에 속하는 point들의 effective number Nk=n=1Nγ(znk)N_k = \sum_{n=1}^N \gamma (z_{nk}) 를 정의하면, 다음의 결과들을 얻는다.
μk=1Nkn=1Nγ(znk)xnΣk=1Nkn=1Nγ(znk)(xnμk)(xnμk)Tπk=NkN\mu_k = \frac{1}{N_k}\sum_{n=1}^N \gamma(z_{nk})x_n \\ \Sigma_k = \frac{1}{N_k} \sum_{n=1}^N\gamma(z_{nk})(x_n - \mu_k)(x_n - \mu_k)^T \\ \pi_k = \frac{N_k}{N}
하지만 위와 같은 결과들은 responsibility의 복잡한 dependency로 인해 closed-form으로 풀어낼 수 없게 된다. 이를 해결하기 위한 iterative scheme이 바로 EM 알고리즘이다.
EM 알고리즘에서는 먼저 μ,Σ,π\mu, \Sigma, \pi 의 초깃값을 설정한 다음, 1) 주어진 parameter 값에 대한 responsibility(latent variable의 posterior probability)를 구하는 E step과 2) 주어진 responsibility에 대해 likelihood function을 최대화하는 parameter value를 (바로 위에서 구한 대로) 찾는 M step을 반복하는 것이다. 이러한 알고리즘이 실제로 likelihood function의 극대화 결과를 가져옴을 증명하는 것은 이 장의 맨 마지막에서 다룬다.

9.3. An Alternative View of EM

이번 섹션에서는 EM 알고리즘에서 latent variable의 핵심 역할을 알아보기 위해 EM 알고리즘을 조금 다른 시각에서 재조명한다. 결론부터 말하자면 latent variable을 명시적으로 표현함으로써 동일한 최적화 문제를 훨씬 더 쉽게 풀어낼 수 있다는 것이다.
모든 model parameter를 θ\theta 라고 표현하면, log likelihood는 다음과 같다.
lnp(Xθ)=ln{Zp(X,Zθ)}\ln p(X|\theta) = \ln\{\sum_Zp(X,Z|\theta)\}
여기서 중요한 점은 바로 summation이 logarithm 안에 나타난다는 것이다. 이 때문에 joint distribution이 exponential family여도 likelihood가 exponential family가 아닌 계산하기 어려운 형태를 가지게 되는 것이다. (exponential family들을 곱하면 exponential family가 되지만, 더하면 그렇지 않다.) 이래서 maximum likelihood를 구하는 게 굉장히 복잡해진다.
근데 위의 Gaussian mixture 예시에서는 그냥 구한 거 아닌가? 일반적으로는 어렵다는 소리인가? ⇒ 일반적으로 어렵다는 게 맞는 것 같다. 위의 예시에서도 구할 수는 있지만 더 어렵기도 하고.
이제 {X,Z}\{X, Z\} 를 complete data set, {X}\{X\} 를 incomplete data set이라고 부르자. 그러면 complete-data log likelihood function의 maximization은 굉장히 간단해짐을 잠시 후에 확인할 수 있다.
하지만 현실 세계에서는 대부분 ZZ 는 관측할 수 있는 값이 아니다. 이름부터 latent variable이지 않은가. 그렇기 때문에 우리는 complete-data log likelihood는 알아낼 수 없고, 대신에 latent variable의 posterior distribution에 대한 complete-data log likelihood의 expectation을 사용한다. 즉 보통 우리는 EM 알고리즘의 M step에서 log likelihood의 expectation을 최대화하게 된다.
Latent variable의 posterior distribution p(ZX,θold)p(Z|X,\theta^{old}) 에 대한, complete-data log likelihood function의 expectation을 다음과 같이 표현한다.
Q(θ,θold)=Zp(ZX,θold)lnp(X,Zθ)Q(\theta, \theta^{old}) = \sum_Zp(Z|X, \theta^{old})\ln p(X, Z|\theta)
이제 위와 같은 expectation QQ 는 logarithm 안에 summation 같은 것들이 없어서 굉장히 쉽게 maximization이 가능하다.
이러한 EM 알고리즘은 maximum likelihood solution뿐 아니라 MAP solution을 찾는데도 사용될 수 있다. 단순히 prior p(θ)p(\theta) 를 설정하고 위의 Q(θ,θold)+lnp(θ)Q(\theta,\theta^{old}) + \ln p(\theta) 를 maximization하면 된다고 한다.

9.3.1 Gaussian mixtures revisited

이제 위의 latent variable view를 mixture Gaussian model에 적용해보자.
일단 X,ZX, Z 의 값이 모두 주어졌다고 가정하자. 그러면 joint likelihood function은 다음과 같다.
p(X,Zμ,Σ,π)=nkπkznkN(xnμk,Σk)znkp(X,Z|\mu,\Sigma,\pi) = \prod_n\prod_k\pi_k^{z_{nk}}N(x_n|\mu_k,\Sigma_k)^{z_{nk}}
위의 식에 로그를 취하면 다음과 같다.
lnp(X,Zμ,Σ,π)=nkznk{lnπk+lnN(xnμk,Σk)}\ln p(X,Z|\mu,\Sigma,\pi) = \sum_n\sum_k z_{nk}\{\ln\pi_k + \ln N(x_n|\mu_k,\Sigma_k)\}
위의 marginal log likelihood와는 달리 kk 에 대한 summation과 logarithm의 위치가 뒤바뀐 것을 알 수 있다. 이제 logarithm이 Gaussian에 직접적으로 조져졌으니, maximum likelihood를 다루는 게 훨씬 쉽다.
znz_{n}xnx_n 이 어디에 속했는지를 알려주는 KK-dimensional one-hot-encoded vector이니, mean과 covariance를 찾는 일은 해당 cluster에 속한 데이터들을 대상으로 single-Gaussian maximization을 하는 것과 아주 동일해진다. (Gaussian들을 다 곱하면 결국 Gaussian이기 때문이다.) πk\pi_k 는 summation constraint가 있으니 Lagrange multiplier를 이용하여 maximization하면 해당 cluster에 속한 데이터들의 비율임을 쉽게 알 수 있다.
이에 따라 complete-data log likelihood function은 exponential family가 되어 closed form에 의해 최대화될 수 있다. 하지만 우리는 가정과 달리 latent variable을 관측할 수 없기 때문에 posterior를 대신하여 사용하게 된다. Posterior는 다음과 같이 표현된다.
p(ZX,μ,Σ,π)nk[πkN(xnμk,Σk)]znkp(Z|X,\mu,\Sigma,\pi) \propto \prod_n\prod_k[\pi_kN(x_n|\mu_k,\Sigma_k)]^{z_{nk}}
여기서 {znk}\{z_{nk}\} 는 독립임을 (d-separation을 통해) 쉽게 알 수 있고, 이에 따라 znkz_{nk} 의 기댓값은 다음과 같이 구할 수 있다. 그냥 E[znk]=znkp(znkX,μ,Σ,π)E[z_{nk}] = \sum z_{nk} p(z_{nk}|X,\mu,\Sigma,\pi) 을 풀어 쓴 것임에 유의하자. 이게 또 우연찮게 responsibility가 된다.
E[znk]=znkznk[πkN(xnμk,Σk)]znkznj[πjN(xnμj,Σj)]znj=πkN(xnμk,Σk)j=1KπjN(xnμj,Σj)=γ(znk)E[z_{nk}] = \frac{\sum_{z_{nk}}z_{nk}[\pi_kN(x_n|\mu_k,\Sigma_k)]^{z_{nk}}}{\sum_{z_{nj}}[\pi_jN(x_n|\mu_j,\Sigma_j)]^{z_{nj}}} \\{}\\ = \frac{\pi_kN(x_n|\mu_k,\Sigma_k)}{\sum_{j=1}^K\pi_jN(x_n|\mu_j,\Sigma_j)} = \gamma(z_{nk})
따라서 complete-data log likelihood의 expectation은 다음과 같다. 아까 전의 log likelihood function의 식에서 znkz_{nk}γ(znk)\gamma(z_{nk}) 로 바뀌었다.
EZ[lnp(X,Zμ,Σ,π)]=nkγ(znk){lnπk+lnN(xnμk,Σk)}E_Z[\ln p(X,Z|\mu,\Sigma,\pi)] = \sum_n\sum_k \gamma(z_{nk})\{\ln\pi_k + \ln N(x_n|\mu_k,\Sigma_k)\}
이제 parameter를 초기화하고, responsibility를 구하고(E step), closed-form으로 parameter를 구하는(M step) 작업을 반복하면 다시 EM 알고리즘이다. 여기서 구하는 M step에서의 parameter 값과 9.2.2에서 구했던 mixture Gaussian의 EM 알고리즘에서의 parameter 값이 동일하다고 한다.

9.3.2 Relation to K-means

K-means 알고리즘은 데이터를 하나의 cluster에만 할당하는 hard assignment를 수행했다. 이와 비교하여 EM 알고리즘은 soft assignment를 수행한다고 볼 수 있다. 여기서는 EM 알고리즘으로부터 hard assignment를 수행하는 K-means 알고리즘을 도출해본다.
자세한 내용은 책을 참고하면 된다. 핵심 내용은 covariance matrix를 ϵI\epsilon I 로 정의하고, responsibility에서 극한 ϵ0\epsilon \rightarrow 0 을 취하면 responsibility 값이 hard assignment가 된다는 것이다. 이 때의 expected complete-data log likelihood 최적화 문제는 K-means 알고리즘의 목적 함수 최적화 문제와 동일함을 확인할 수 있다.

9.3.3 Mixtures of Bernoulli distributions

Bernoulli distribution들의 linear combination인 mixture Bernoulli model에서도 같은 방식으로 latent variable을 도입하여 EM 알고리즘을 수행할 수 있음을 설명한다. Mixture Gaussian과 아주 동일한 과정이니 생략한다. 교재를 참고하자.

9.3.4 EM for Bayesian linear regression

3장에서 다루었던 evidence approximation for Bayesian linear regression 문제를 EM 알고리즘으로 풀어낸다. Model evidence를 구할 때 모델 parameter ww 에 대하여 marginalization을 수행하는데, 바로 이 ww 를 latent variable로 간주하여 EM 알고리즘을 적용하면 된다. 역시 교재를 참고해보자.

9.4. The EM Algorithm in General

EM 알고리즘은 latent variable을 가지는 확률 모델의 MLE solution을 찾아내는 일반적인 알고리즘이라고 한다. 위에서 여러 가지 예시를 살펴보았으니, 이제 EM 알고리즘의 보다 일반적인 형태를 살펴보고, 이러한 EM 알고리즘이 실제로 likelihood function을 최대화하게 된다는 것을 증명해보자.
마찬가지로 observed variable XX 와 hidden(latent) variable ZZ 가 주어졌다고 하자. Joint distribution p(X,Zθ)p(X, Z|\theta) 는 parameter θ\theta 를 가진다. 우리의 목표는 다음과 같은 likelihood function을 최대화하는 것이다.
p(Xθ)=Zp(X,Zθ)p(X|\theta) = \sum_Zp(X,Z|\theta)
여기서 위에서 본 것처럼 p(Xθ)p(X|\theta) 를 직접 최적화하는 것은 어렵고, 이에 비해 p(X,Zθ)p(X,Z|\theta) 를 최적화하는 문제는 훨씬 쉬운 상황이라고 하자. 이제 latent variable에 대한 distribution q(Z)q(Z) 를 정의하면, 다음의 decomposition이 성립한다.
lnp(Xθ)=L(q,θ)+KL(qp)L(q,θ)=Zq(Z)ln{p(X,Zθ)q(Z)}KL(qp)=Zq(Z)ln{p(ZX,θ)q(Z)}\ln p(X|\theta) = L(q,\theta) + KL(q||p) \\{}\\ L(q,\theta) = \sum_Zq(Z)\ln\{\frac{p(X,Z|\theta)}{q(Z)}\} \\ KL(q||p) = -\sum_Zq(Z)\ln\{\frac{p(Z|X,\theta)}{q(Z)}\}
우선 KL divergence가 0 이상의 실수인 것은 아주 잘 알려진 사실이며, 이에 따라 L(q,θ)L(q,\theta) 는 log likelihood의 lower bound를 형성한다. 이제 이를 이용하여 EM 알고리즘이 실제로 log likelihood를 최대화함을 증명할 수 있다.
먼저 E step에서, θ\theta 를 고정한 상태에서 q(Z)q(Z) 에 대하여 lower bound L(q,θ)L(q,\theta) 를 최대화한다. 그런데 lnp(Xθ)\ln p(X|\theta)q(Z)q(Z) 와 독립적이기 때문에 값이 변하지 않고, 이에 따라 q(Z)q(Z) 에 대해 L(q,θ)L(q,\theta)를 최대화하기 위해서는 KL(qp)KL(q||p) 의 값이 0이 되는 qq 를 찾게 된다.
이제 M step에서, 이번에는 q(Z)q(Z) 를 고정한 상태에서 θ\theta 에 대하여 lower bound L(q,θ)L(q,\theta) 를 최대화한다. 이는 lower bound L(q,θ)L(q,\theta) 의 값을 증가시킬 것이고, 결과적으로 log likelihood의 값을 증가시키게 된다. 이 때 qq 는 변하기 이전 parameter에 대해 최적화되었으므로, 새로운 pp 에 대한 KL divergence는 0보다 큰 값을 가지게 된다. 이에 따라 log likelihood의 증가분은 lower bound보다 크다. 그리고 KL divergence가 0보다 큰 값을 가지게 되었으므로 다음 E step에서 또 lower bound가 증가할 수 있게 된다.
이 때 q(Z)=p(ZX,θold)q(Z) = p(Z|X,\theta^{old}) 를(즉 KL divergence가 0이 되는 q(Z)q(Z) 를) 위의 L(q,θ)L(q,\theta) 의 정의에 대입하면 다음과 같은 식을 얻는데,
L(q,θ)=Zp(ZX,θold)lnp(X,Zθ)Zp(ZX,θold)lnp(ZX,θold)=Q(θ,θold)+constL(q,\theta) = \sum_Zp(Z|X,\theta^{old})\ln p(X,Z|\theta) - \sum_Z p(Z|X,\theta^{old}) \ln p(Z|X,\theta^{old}) \\ = Q(\theta,\theta^{old}) + \text{const}
이는 E step 이후의 lower bound를 뜻한다. 이는 다시 M step에서의 lower bound maximization은 expected complete-data log likelihood의 maximization과 동일하며, 따라서 expected complete-data log likelihood를 최대화하는 것이 log likelihood를 최대화하는 것과 동일함을 뜻한다!
EM 알고리즘을 이용하여 MLE solution이 아니라 MAP solution을 찾을 수도 있다. 먼저 prior p(θ)p(\theta) 를 정의하면 log of posterior는 다음과 같다.
lnp(θX)=lnp(θ,X)lnp(X)\ln p(\theta|X) = \ln p(\theta,X) - \ln p(X)
이제 위의 decomposition을 이용하면 다음과 같은데,
lnp(θX)=L(q,θ)+KL(qp)+lnp(θ)lnp(X)L(q,θ)+lnp(θ)lnP(X)\ln p(\theta|X) = L(q,\theta) + KL(q||p) + \ln p(\theta) - \ln p(X) \\ \geq L(q,\theta) + \ln p(\theta) - \ln P(X)
위와 같은 lower bound(맨 오른쪽 식)를 대상으로 qq 에 대해 최적화하는 E step과 θ\theta 에 대해 최적화하는 M step을 반복하면 되는 것이다.
이외에도 E step 혹은 M step의 연산이 intractable 할 경우에 각각의 step을 (SGD 같은) nonlinear optimization 알고리즘을 통해 최적화하는 GEM(Generalized EM), lower bound에 대한 partial optimization을 수행하는 EM, EM의 sequential scheme (혹은 batch version) 등을 소개하며 장을 마무리한다.
이제... 오랜 PRML 스터디의 끝이 조금씩 보이기 시작한다...
E.O.D.