机器学习实战ByMatlab（三）K-means算法

机器学习实战ByMatlab（三）K-means算法

http://blog.csdn.net/llp1992/article/details/45095935

K-means算法属于无监督学习聚类算法，其计算步骤还是挺简单的，思想也挺容易理解，而且还可以在思想中体会到EM算法的思想。

K-means 算法的优缺点：

1.优点：容易实现 
2.缺点：可能收敛到局部最小值，在大规模数据集上收敛较慢

使用数据类型：数值型数据

以往的回归算法、朴素贝叶斯、SVM等都是有类别标签y的，因此属于有监督学习，而K-means聚类算法只有x，没有y

在聚类问题中，我们的训练样本是

其中每个Xi都是n维实数。

样本数据中没有了y，K-means算法是将样本聚类成k个簇，具体算法如下： 
1、随机选取K个聚类质心点，记为 

2、重复以下过程直到收敛

{ 
对每个样例 i ，计算其应该属于的类： 

对每个类 j ，重新计算质心： 

}

其中K是我们事先给定的聚类数目，Ci 表示样本 i 与K个聚类中最近的那个类，Ci的值是1到K中的一个，质心uj代表我们对属于同一个类的样本中心的猜测。解释起来就是，

第一步：天空上的我们随机抽取K个星星作为星团的质心，然后对于每一个星星 i，我们计算它到每一个质心uj的距离，选取其中距离最短的星团作为Ci，这样第一步每个星星都有了自己所属于的星团；

第二步：对每个星团Ci，我们重新计算它的质心uj（计算方法为对属于该星团的所有点的坐标求平均）不断重复第一步和第二步直到质心变化很小或者是不变。

然后问题来了，怎么样才算质心变化很小或者是不变？或者说怎么判定？答案就是畸变函数（distortion function），定义如下：

J函数表示每个样本点到其质心的距离平方和，K-means的收敛就是要将 J 调整到最小，假设当前 J 值没有达到最小值，那么可以先固定每个类的质心 uj ，调整每个样例的类别 Ci 来时 J 函数减少。同样，固定 Ci ，调整每个类的质心 uj也可以是 J 减少。这两个过程就是内循环中使 J 单调变小的过程。当 J 减小到最小的时候， u 和 c 也同时收敛。（该过程跟EM算法其实还是挺像的）理论上可能出现多组 u 和 c 使 J 取得最小值，但这种情况实际上很少见。

由于畸变函数 J 是非凸函数，所以我们不能保证取得的最小值一定是全局最小值，这说明k-means算法质心的初始位置的选取会影响到最后最小值的获取。不过一般情况下，k-means算法达到的局部最优已经满足要求。如果不幸代码陷入局部最优，我们可以选取不同的初始值跑多几遍 k-means 算法，然后选取其中最小的 J 对应的 u 和 c 输出。

另一种收敛判断：

实际我们编写代码的时候，还可以通过判断“每个点被分配的质心是否改变”这个条件来判断聚类是否已经收敛

而上面所说的畸变函数则可以用来评估收敛的效果，具体将会在下面的实例中体现。

Matlab 实现

<code class="language-matlab hljs  has-numbering" style="display: block; padding: 0px; color: inherit; box-sizing: border-box; font-family: 'Source Code Pro', monospace;font-size:undefined; white-space: pre; border-radius: 0px; word-wrap: normal; background: transparent;">
<span class="hljs-function" style="box-sizing: border-box;"><span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">function</span> <span class="hljs-title" style="box-sizing: border-box;">kMeans</span></span>
clc
clear
K = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">4</span>;
dataSet = load(<span class="hljs-string" style="color: rgb(0, 136, 0); box-sizing: border-box;">'testSet.txt'</span>);
<span class="hljs-matrix" style="box-sizing: border-box;">[row,col]</span> = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">size</span>(dataSet);
<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 存储质心矩阵</span>
centSet = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">zeros</span>(K,col);
<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 随机初始化质心</span>
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">for</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>= <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>:col
    minV = min(dataSet(:,<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>));
    rangV = max(dataSet(:,<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>)) - minV;
    centSet(:,<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>) = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">repmat</span>(minV,<span class="hljs-matrix" style="box-sizing: border-box;">[K,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>]</span>) + rangV*<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">rand</span>(K,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>);
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>

<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 用于存储每个点被分配的cluster以及到质心的距离</span>
clusterAssment = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">zeros</span>(row,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">2</span>);
clusterChange = true;
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">while</span> clusterChange
    clusterChange = false;
    <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 计算每个点应该被分配的cluster</span>
    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">for</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span> = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>:row
        <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 这部分可能可以优化</span>
        minDist = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">10000</span>;
        minIndex = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">0</span>;
        <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">for</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span> = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>:K
            distCal = distEclud(dataSet(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>,:) , centSet(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span>,:));
            <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">if</span> (distCal < minDist)
                minDist = distCal;
                minIndex = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span>;
            <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>
        <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>
        <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">if</span> minIndex ~= clusterAssment(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>)            
            clusterChange = true;
        <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>
        clusterAssment(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>) = minIndex;
        clusterAssment(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">2</span>) = minDist;
    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>

    <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 更新每个cluster 的质心</span>
    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">for</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span> = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>:K
        simpleCluster = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">find</span>(clusterAssment(:,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>) == <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span>);
        centSet(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">j</span>,:) = mean(dataSet(<span class="hljs-transposed_variable" style="box-sizing: border-box;">simpleCluster'</span>,:));
    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>

<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>
figure
<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">%scatter(dataSet(:,1),dataSet(:,2),5)</span>
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">for</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span> = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>:K
    pointCluster = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">find</span>(clusterAssment(:,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>) == <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">i</span>);
    scatter(dataSet(pointCluster,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>),dataSet(pointCluster,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">2</span>),<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">5</span>)
    hold on
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>
<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">%hold on</span>
scatter(centSet(:,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span>),centSet(:,<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">2</span>),<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">300</span>,<span class="hljs-string" style="color: rgb(0, 136, 0); box-sizing: border-box;">'+'</span>)
hold off

<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span>

<span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">% 计算欧式距离</span>
<span class="hljs-function" style="box-sizing: border-box;"><span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">function</span> <span class="hljs-title" style="box-sizing: border-box;">dist</span> = <span class="hljs-title" style="box-sizing: border-box;">distEclud</span><span class="hljs-params" style="color: rgb(102, 0, 102); box-sizing: border-box;">(vecA,vecB)</span></span>
    dist  = <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">sqrt</span>(sum(power((vecA-vecB),<span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">2</span>)));
<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">end</span></code><ul class="pre-numbering" style="box-sizing: border-box; position: absolute; width: 50px; top: 0px; left: 0px; margin: 0px; padding: 6px 0px 40px; border-right-width: 1px; border-right-style: solid; border-right-color: rgb(221, 221, 221); list-style: none; text-align: right; background-color: rgb(238, 238, 238);"><li style="box-sizing: border-box; padding: 0px 5px;">1</li><li style="box-sizing: border-box; padding: 0px 5px;">2</li><li style="box-sizing: border-box; padding: 0px 5px;">3</li><li style="box-sizing: border-box; padding: 0px 5px;">4</li><li style="box-sizing: border-box; padding: 0px 5px;">5</li><li style="box-sizing: border-box; padding: 0px 5px;">6</li><li style="box-sizing: border-box; padding: 0px 5px;">7</li><li style="box-sizing: border-box; padding: 0px 5px;">8</li><li style="box-sizing: border-box; padding: 0px 5px;">9</li><li style="box-sizing: border-box; padding: 0px 5px;">10</li><li style="box-sizing: border-box; padding: 0px 5px;">11</li><li style="box-sizing: border-box; padding: 0px 5px;">12</li><li style="box-sizing: border-box; padding: 0px 5px;">13</li><li style="box-sizing: border-box; padding: 0px 5px;">14</li><li style="box-sizing: border-box; padding: 0px 5px;">15</li><li style="box-sizing: border-box; padding: 0px 5px;">16</li><li style="box-sizing: border-box; padding: 0px 5px;">17</li><li style="box-sizing: border-box; padding: 0px 5px;">18</li><li style="box-sizing: border-box; padding: 0px 5px;">19</li><li style="box-sizing: border-box; padding: 0px 5px;">20</li><li style="box-sizing: border-box; padding: 0px 5px;">21</li><li style="box-sizing: border-box; padding: 0px 5px;">22</li><li style="box-sizing: border-box; padding: 0px 5px;">23</li><li style="box-sizing: border-box; padding: 0px 5px;">24</li><li style="box-sizing: border-box; padding: 0px 5px;">25</li><li style="box-sizing: border-box; padding: 0px 5px;">26</li><li style="box-sizing: border-box; padding: 0px 5px;">27</li><li style="box-sizing: border-box; padding: 0px 5px;">28</li><li style="box-sizing: border-box; padding: 0px 5px;">29</li><li style="box-sizing: border-box; padding: 0px 5px;">30</li><li style="box-sizing: border-box; padding: 0px 5px;">31</li><li style="box-sizing: border-box; padding: 0px 5px;">32</li><li style="box-sizing: border-box; padding: 0px 5px;">33</li><li style="box-sizing: border-box; padding: 0px 5px;">34</li><li style="box-sizing: border-box; padding: 0px 5px;">35</li><li style="box-sizing: border-box; padding: 0px 5px;">36</li><li style="box-sizing: border-box; padding: 0px 5px;">37</li><li style="box-sizing: border-box; padding: 0px 5px;">38</li><li style="box-sizing: border-box; padding: 0px 5px;">39</li><li style="box-sizing: border-box; padding: 0px 5px;">40</li><li style="box-sizing: border-box; padding: 0px 5px;">41</li><li style="box-sizing: border-box; padding: 0px 5px;">42</li><li style="box-sizing: border-box; padding: 0px 5px;">43</li><li style="box-sizing: border-box; padding: 0px 5px;">44</li><li style="box-sizing: border-box; padding: 0px 5px;">45</li><li style="box-sizing: border-box; padding: 0px 5px;">46</li><li style="box-sizing: border-box; padding: 0px 5px;">47</li><li style="box-sizing: border-box; padding: 0px 5px;">48</li><li style="box-sizing: border-box; padding: 0px 5px;">49</li><li style="box-sizing: border-box; padding: 0px 5px;">50</li><li style="box-sizing: border-box; padding: 0px 5px;">51</li><li style="box-sizing: border-box; padding: 0px 5px;">52</li><li style="box-sizing: border-box; padding: 0px 5px;">53</li><li style="box-sizing: border-box; padding: 0px 5px;">54</li><li style="box-sizing: border-box; padding: 0px 5px;">55</li><li style="box-sizing: border-box; padding: 0px 5px;">56</li><li style="box-sizing: border-box; padding: 0px 5px;">57</li><li style="box-sizing: border-box; padding: 0px 5px;">58</li><li style="box-sizing: border-box; padding: 0px 5px;">59</li><li style="box-sizing: border-box; padding: 0px 5px;">60</li><li style="box-sizing: border-box; padding: 0px 5px;">61</li><li style="box-sizing: border-box; padding: 0px 5px;">62</li><li style="box-sizing: border-box; padding: 0px 5px;">63</li><li style="box-sizing: border-box; padding: 0px 5px;">64</li></ul>

效果如下：

这是正常分类的情况，很明显被分为了4个类，不同颜色代表不同的类，cluster的质心为 “ + ”

当然，这只是其中一种情况，很有可能我们会出现下面这种情况：

这就是上面所说的，K-means的缺点之一，随机初始点的选择可能会让算法陷入局部最优解，这时候我们只需重新运行一次程序即可。

至于每一个看似都可以正常聚类的情况呢，我们则利用上面所说的“畸变函数”来衡量聚类的效果，当然是J越小聚类效果越好。

实际使用的时候，我们只需多次运行程序，选取J最小的聚类效果。