機器學習實戰ByMatlab（4）：二分K-means算法
原文出處： Liu_LongPo的專欄（@Liu_LongPo）
前面我們在是實現K-means算法的時候，提到了它本身存在的缺陷：
1.可能收斂到局部最小值
2.在大規模數據集上收斂較慢
對于上一篇博文最后說的，當陷入局部最小值的時候，處理方法就是多運行幾次K-means算法，然后選擇畸變函數J較小的作為最佳聚類結果。這樣的說法顯然不能讓我們接受，我們追求的應該是一次就能給出接近最優的聚類結果。
其實K-means的缺點的根本原因就是：對K個質心的初始選取比較敏感。質心選取得不好很有可能就會陷入局部最小值。
基于以上情況，有人提出了二分K-means算法來解決這種情況，也就是弱化初始質心的選取對最終聚類效果的影響。
二分K-means算法
在介紹二分K-means算法之前我們先說明一個定義：SSE（Sum of Squared Error）,也就是誤差平方和，它是用來度量聚類效果的一個指標。其實SSE也就是我們在K-means算法中所說的畸變函數：
SSE計算的就是一個cluster中的每個點到質心的平方差，它可以度量聚類的好壞。顯然SSE越小，說明聚類效果越好。
二分K-means算法的主要思想：
首先將所有點作為一個簇，然后將該簇一分為二。之后選擇能最大程度降低聚類代價函數（也就是誤差平方和）的簇劃分為兩個簇。以此進行下去，直到簇的數目等于用戶給定的數目k為止。
二分k均值算法的偽代碼如下：
將所有數據點看成一個簇
當簇數目小于k時
對每一個簇
計算總誤差
在給定的簇上面進行k-均值聚類（k=2）
計算將該簇一分為二后的總誤差
選擇使得誤差最小的那個簇進行劃分操作
Matlab 實現
function bikMeans
%%
clc
clear
close all
%%
biK = 4;
biDataSet = load('testSet.txt');
[row,col] = size(biDataSet);
% 存儲質心矩陣
biCentSet = zeros(biK,col);
% 初始化設定cluster數量為1
numCluster = 1;
%第一列存儲每個點被分配的質心，第二列存儲點到質心的距離
biClusterAssume = zeros(row,2);
%初始化質心
biCentSet(1,:) = mean(biDataSet)
for i = 1:row
biClusterAssume(i,1) = numCluster;
biClusterAssume(i,2) = distEclud(biDataSet(i,:),biCentSet(1,:));
end
while numCluster < biK
minSSE = 10000;
%尋找對哪個cluster進行劃分最好，也就是尋找SSE最小的那個cluster
for j = 1:numCluster
curCluster = biDataSet(find(biClusterAssume(:,1) == j),:);
[spiltCentSet,spiltClusterAssume] = kMeans(curCluster,2);
spiltSSE = sum(spiltClusterAssume(:,2));
noSpiltSSE = sum(biClusterAssume(find(biClusterAssume(:,1)~=j),2));
curSSE = spiltSSE + noSpiltSSE;
fprintf('第%d個cluster被劃分后的誤差為：%f \n' , [j, curSSE])
if (curSSE < minSSE)
minSSE = curSSE;
bestClusterToSpilt = j;
bestClusterAssume = spiltClusterAssume;
bestCentSet = spiltCentSet;
end
end
bestClusterToSpilt
bestCentSet
%更新cluster的數目
numCluster = numCluster + 1;
bestClusterAssume(find(bestClusterAssume(:,1) == 1),1) = bestClusterToSpilt;
bestClusterAssume(find(bestClusterAssume(:,1) == 2),1) = numCluster;
% 更新和添加質心坐標
biCentSet(bestClusterToSpilt,:) = bestCentSet(1,:);
biCentSet(numCluster,:) = bestCentSet(2,:);
biCentSet
% 更新被劃分的cluster的每個點的質心分配以及誤差
biClusterAssume(find(biClusterAssume(:,1) == bestClusterToSpilt),:) = bestClusterAssume;
end
figure
%scatter(dataSet(:,1),dataSet(:,2),5)
for i = 1:biK
pointCluster = find(biClusterAssume(:,1) == i);
scatter(biDataSet(pointCluster,1),biDataSet(pointCluster,2),5)
hold on
end
%hold on
scatter(biCentSet(:,1),biCentSet(:,2),300,'+')
hold off
end
% 計算歐式距離
function dist = distEclud(vecA,vecB)
dist = sum(power((vecA-vecB),2));
end
% K-means算法
function [centSet,clusterAssment] = kMeans(dataSet,K)
[row,col] = size(dataSet);
% 存儲質心矩陣
centSet = zeros(K,col);
% 隨機初始化質心
for i= 1:col
minV = min(dataSet(:,i));
rangV = max(dataSet(:,i)) - minV;
centSet(:,i) = repmat(minV,[K,1]) + rangV*rand(K,1);
end
% 用于存儲每個點被分配的cluster以及到質心的距離
clusterAssment = zeros(row,2);
clusterChange = true;
while clusterChange
clusterChange = false;
% 計算每個點應該被分配的cluster
for i = 1:row
% 這部分可能可以優化
minDist = 10000;
minIndex = 0;
for j = 1:K
distCal = distEclud(dataSet(i,:) , centSet(j,:));
if (distCal < minDist)
minDist = distCal;
minIndex = j;
end
end
if minIndex ~= clusterAssment(i,1)
clusterChange = true;
end
clusterAssment(i,1) = minIndex;
clusterAssment(i,2) = minDist;
end
% 更新每個cluster 的質心
for j = 1:K
simpleCluster = find(clusterAssment(:,1) == j);
centSet(j,:) = mean(dataSet(simpleCluster',:));
end
end
end 算法迭代過程如下
biCentSet =
-0.1036 0.0543
0 0
0 0
0 0
第1個cluster被劃分后的誤差為：792.916857
bestClusterToSpilt =
1
bestCentSet =
-0.2897 -2.8394
0.0825 2.9480
biCentSet =
-0.2897 -2.8394
0.0825 2.9480
0 0
0 0
第1個cluster被劃分后的誤差為：409.871545
第2個cluster被劃分后的誤差為：532.999616
bestClusterToSpilt =
1
bestCentSet =
-3.3824 -2.9473
2.8029 -2.7315
biCentSet =
-3.3824 -2.9473
0.0825 2.9480
2.8029 -2.7315
0 0
第1個cluster被劃分后的誤差為：395.669052
第2個cluster被劃分后的誤差為：149.954305
第3個cluster被劃分后的誤差為：393.431098
bestClusterToSpilt =
2
bestCentSet =
2.6265 3.1087
-2.4615 2.7874
biCentSet =
-3.3824 -2.9473
2.6265 3.1087
2.8029 -2.7315
-2.4615 2.7874
最終效果圖
運用二分K-means算法進行聚類的時候，不同的初始質心聚類結果還是會稍微有點不同，因為實際上這也只是弱化隨機質心對聚類結果的影響而已，并不能消除其影響，不過最終還是能收斂到全局最小。