Description
Problem 1
Consider the p-dimensional Gaussian Mixture Model:
xi ∼
X
K
k
πkp(x|µk, Σk, Ck), for i = 1, . . . , N, (1)
where p(x|µk, Σk) = 1
πp/2|Σk|
1/2 e
− 1
2
(x−µk)
0Σ
−1
k
(x−µk)
.
In the special case where Σk = σI ∀k = 1 . . . , K, discuss the connections between the
K-means algorithm and EM for fitting model (1). Additionally, show that as σ → 0
the two methods coincide.
The Gaussian Mixture Model with EM algorithm is a kind of soft K-Means algorithm. In the E step, EM algorithm we calculate the expectation of the log likelihood
as follows:
Q(θ|θ(t)) = Ez|θ(t)
[−
1
2
X
N
i=1
X
K
k=1
[(x−µk)
0Σ
−1
k
(x−µk) +Logπk −LogΣk]×δ(zi = k|xi))]
To simply the above expression, we can obtain:
Q(θ|θ(t)) = −
1
2
X
N
i=1
X
K
k=1
[(x−µk)
0Σ
−1
k
(x−µk) +Logπk −LogΣk]×Ez|θ(t)
[δ(zi = k|xi))]
(2)
= −
1
2
X
N
i=1
X
K
k=1
[(x − µk)
0Σ
−1
k
(x − µk) + Logπk − LogΣk] × P(zi = k|xi)
Here the P(zi = k|xi) = π
(t)
i,k is the posterior
P(zi = k|xi) = π
(t)
i,k =
P(xi
|zi = k, θ(t)
)P r(zi = k|θ
(t)
)
PK
k=1 P(xi
|zi = k, θ(t)
)P r(zi = k|θ
(t)
)
1
with
P(xi
|zi = k, θ(t)
) = 1
π
p/2
|Σ
(t)
k
|
1/2
e
− 1
2
(xi−µ
(t)
k
)
0Σ
−1
k
(xi−µ
(t)
k
)
and
P r(zi = k|θ
(t)
) = π
(t)
k
.
Here the posterior P(zi = k|xi) gives us the probability for a data point xi coming
from the component Ck. Thus, the EM algorithm is implementing soft labeling. If
Σk goes to zero in the EM algorithm, the P(xi
|zi = k, θ(t)
) goes to delta distribution.
This means, we have:
P(xi
|zi = k, θ(t)
) = δ(xi − µ
(t)
k
)
In addition, the posterior becomes
P(zi = k|xi) = π
(t)
i,k = δ(xi − µ
(t)
k
)π
(t)
k
/
X
K
k=1
δ(xi − µ
(t)
k
)π
(t)
k
(3)
which effectively gives hard labeling (since they are step functions over the data
space). A visualization of the above function (simplified version) can be found in
page 512 of the text book. Thus, given Σk → 0, the E step of the EM algorithm
converges to the (hard/mutual exclusive) labeling step of the K-Means algorithm.
The M step maximize the Q(θ|θ(t)) in the E step with respect to the model
parameters, which yields
µ
(t)
k =
X
N
i=1
π
(t)
i,kxi/
X
N
i=1
π
(t)
i,k (4)
If Σk goes to zero, the denominator of the above expression will converge to 1, sine
the posteriors π
(t)
i,k will converge into the forms given by Eq.(3). Accordingly, we can
write Eq.(4) as µ
(t)
k =
PN
i=1 π
(t)
i,kxi
, which is effectively the mean of all the data points
in the component Ck. This is actually the second step of the K-Means algorithm
when the new center of each cluster is calculated.
To sum up, the GMM with EM algorithm converges to the K-Means algorithm
when the Σk (or the σ) goes to zero.
Problem 2
Download “ClusterSet1.txt” from the course webpage. Apply the k-means clustering
procedure to this data set. You may code this up from scratch or use the built in
functions in either R or Matlab. Discuss how you selected ‘K’, and why you believe
it is correct.
I used R package to do the K-Means clustering. I also used the same K-Means
function in R to calculate the Within Sum of Squares (WSS) and the Total Sum of
Squares (TSS). Then I plotted the w = W SS/T SS values with increasing value of
K. The scree plot is shown in below.
2
Figure 1: Scree plot of K-Means clustering, the optimal K is 12. The y axis is the
ratio of Within Sum of Squares (WSS) and Total Sum of Squares (TSS).
It is very obvious from the scree plot that the optimal K is 12, where the w value
reaches minimal while K is kept still small. The w values for K=1,2,3,…,20 can also be
obtained from the R package: 1.000000000, 0.617371303, 0.511107268, 0.357233846,
.253278036, 0.117447680, 0.105069464, 0.068609355, 0.041399420, 0.011358326, 0.067122556,
0.008019367, 0.007076398, 0.006339394, 0.033207078, 0.006506954, 0.005367669, 0.005656337,
0.006719339, 0.031540968. It is also obvious from reading the data that the w stabilizes at very small values once K reaches 12. The clustering result is shown in Fig.2.
(Note different clusters may be dyed with one color.)
Problem 3
Recall that the EM algorithm for fitting model (1) iterates over the following updates:
For t = 1, . . . , T
1. π
(t)
i,k = p(xi ∈ Ck|µ
(t−1)
k
, Σ
(t−1)
k
)
2. µ
(t)
k =
PN
i=1 π
(t)
i,kxi/
PN
i=1 π
(t)
i,k
3. Σ
(t)
k =
PN
i=1 π
(t)
i,k(xi − µ
(t)
k
)(xi − µ
(t)
k
)
0/
PN
i=1 π
(t)
i,k
part a
Given πi,ks, show that the M.L.E for µks are given by PN
i=1 πi,kxi/
PN
i=1 πi,k.
With the Q(θ, θ(t)) as derived in Eq.(2), the µ
(t)
k
can be obtained by solving
∂Q(θ, θ(t))
∂µk
= 0 (5)
3
Figure 2: Result of K-Means clustering. Different pairwise combinations of (v1, v2, v3)
are used to display the clustering result.
This yields
X
N
i=1
X
K
k
0=1
∂[(xi − µk
0)
TΣ
−1
k
0 (xi − µk
0)πik0]
∂µk
= 0 (6)
or
X
N
i=1
∂[(xi − µk)
TΣ
−1
k
(xi − µk)πik]
∂µk
= 0 (7)
According to vector calculus, the above equation can be simplified into
X
N
i=1
(xi − µk)πik = 0 (8)
which ultimately gives
µ
(t)
k =
X
N
i=1
π
(t)
ik xi/
X
N
i=1
π
(t)
ik (9)
4
part b
Given µks and πi,ks, show that the M.L.E for Σks are given by PN
i=1 πi,k(xi −µk)(xi −
µk)
0/
PN
i=1 πi,k (Given µks).
Similarly, with the Q(θ, θ(t)) as derived in Eq.(2), the Σ(t)
k
can be obtained by
solving
∂Q(θ, θ(t))
∂Σk
= 0 (10)
which is equivalent to
X
N
i=1
X
K
k
0=1
∂[Log|Σ
−1
k
0 | − (xi − µk
0)
TΣ
−1
k
0 (xi − µk
0)πik0]
∂Σk
= 0 (11)
or
X
N
i=1
∂[Log|Σ
−1
k
| − (xi − µk)
TΣ
−1
k
(xi − µk)πik]
∂Σk
= 0 (12)
This can be simplified into
X
N
i=1
[Σk − (xi − µk)(xi − µk)
T
]πik = 0 (13)
which ultimately gives
Σ
(t)
k =
X
N
i=1
π
(t)
ik (xi − µ
(t)
k
)(xi − µ
(t)
k
)
T
/
X
N
i=1
π
(t)
ik (14)
part c
Implement the EM algorithm and using ‘ClusterSet1.txt’, compare results to those
found in Problem 2.
For the ’ClusterSet1.txt’ data, the implementation of the EM algorithm produced
results as shown in Table 1. (See Appendix for code.)
1 2 3 4 5 6 7 8 9 10 11 12
97.00 94.00 13.00 103.00 17.00 105.00 109.00 96.00 82.00 110.00 93.00 81.00
0.10 0.10 0.01 0.10 0.02 0.10 0.11 0.10 0.08 0.11 0.09 0.08
Table 1: The number of data points in each cluster and the estimated values for πk.
(k = 1, 2, …, 12).
The number of data points in each cluster is evenly distributed from the EM algorithm. Also, we can calculate the Within Sum of Square (WSS) and Total Sum of
Square (TSS) for the EM clustering results. Then we can obtain w =WSS/TSS=0.008116935.
The K-Means with 12 clusters has w =0.008019367. Thus, the w(EM) is slightly larger
than w(K-Means), meaning the performance of K-Means is slightly better.
5
Problem 4
Using a hierarchical clustering method with ‘ClusterSet1.txt’, compare results to those
found in Problem 2 and 3. Show dendrograms, and discuss the distance function you
settled on for your link function.
In this problem I used hierarchical clustering with Ward’s method. The distance
function I used is the Within Sum of Squares of Euclidean Distances (WSSED).
W SSED =
X
K
k=1
X
i<j,xiCk,xj Ck
kxi − xjk
2
2
Here kxi − xjk2
is the Euclidean distance between data point xi and data point xj
that belong to the same component Ck. The dendrograms is shown in Fig.3.
Figure 3: Dendrogram of Hierarchical Clustering with the Ward’s method. The Hight
is the WSSED.
From the dendrogram plotted for this Hierarchical algorithm, we can easily find
the best number of clusters is 9, see the red line in Fig.3 (which crosses 9 branches
in total). After labeling each data point and calculating the w = W SS/T SS for
the clustering result with 9 clusters. We obtain w = 0.01152203, which is a little
bigger than the result of K-Means (0.008019367) and the result of EM (0.008116935).
Note, the total number of clusters in the latter two algorithms is 11 instead of 9.
Accordingly, we cannot claim this Hierarchical clustering result is worse than those
of the K-Means algorithm and the EM algorithm.
Problem 5
Download “ClusterSet2.txt” from the course webpage. Using any method you find
appropriate, determine the number of clusters and and assignment labels for each
data point.
6
Here I used the Hierarchical Clustering with the Ward’s method. The Dendrogram
below clearly indicates we should set the number of clusters to be 10, see the red line
in Fig.4 (which crosses 10 branches in total).
Figure 4: Dendrogram of Hierarchical Clustering with the Ward’s method.
The labels for x1, x2, …, x1000 given by R code are shown in Fig.5.
Figure 5: The labels for the data points x1, x2, …, x1000, a result of the Hierarchical
Clustering with the Ward’s method.
7