Kmeans++

Kmeans算法源码

标签（空格分隔）： 算法学习

---

（1）算法概述

首先要来了解的一个概念就是聚类，简单地说就是把相似的东西分到一组，同 Classification (分类)不同，对于一个 classifier ，通常需要你告诉它“这个东西被分为某某类”这样一些例子，理想情况下，一个 classifier 会从它得到的训练集中进行“学习”，从而具备对未知数据进行分类的能力，这种提供训练数据的过程通常叫做 supervised learning (监督学习)，而在聚类的时候，我们并不关心某一类是什么，我们需要实现的目标只是把相似的东西聚到一起，因此，一个聚类算法通常只需要知道如何计算相似 度就可以开始工作了，因此 clustering 通常并不需要使用训练数据进行学习，这在 Machine Learning 中被称作 unsupervised learning (无监督学习)。

　　我们经常接触到的聚类分析，一般都是数值聚类，一种常见的做法是同时提取 N 种特征，将它们放在一起组成一个 N 维向量，从而得到一个从原始数据集合到 N 维向量空间的映射——你总是需要显式地或者隐式地完成这样一个过程，然后基于某种规则进行分类，在该规则下，同组分类具有最大的相似性。

　　假设我们提取到原始数据的集合为(x1, x2, …, xn)，并且每个xi为d维的向量，K-means聚类的目的就是，在给定分类组数k（k ≤ n）值的条件下，将原始数据分成k类
S = {S1, S2, …, Sk}，在数值模型上，即对以下表达式求最小值：
\underset{\mathbf{S}} {\operatorname{arg,min}} \sum_{i=1}^{k} \sum_{\mathbf x_j \in S_i} \left| \mathbf x_j - \boldsymbol\mu_i \right|^2
这里μi 表示分类Si 的平均值。

　　那么在计算机编程中，其又是如何实现的呢？其算法步骤一般如下：

1、从D中随机取k个元素，作为k个簇的各自的中心。

2、分别计算剩下的元素到k个簇中心的相异度，将这些元素分别划归到相异度最低的簇。

3、根据聚类结果，重新计算k个簇各自的中心，计算方法是取簇中所有元素各自维度的算术平均数。

4、将D中全部元素按照新的中心重新聚类。

5、重复第4步，直到聚类结果不再变化。

6、将结果输出。

（2）matlab源代码实现

% 【函数描述】C=bsxfun(fun,A,B)：两个数组间元素逐个计算，fun是函数句柄或者m文件，也可以为如下内置函数 %          @plus 加 %          @minus 减 %          @times 数组乘 %          @rdivide 左除 %          @ldivide 右除 % For example： 如何将一个矩阵的每行或每列元素分别扩大不同的倍数？%               如[1 2 3;4 5 6 ;7 8 9],第一列元素乘以1，第二列元素以2，第三列元素乘以4。% 利用bsxfun函数，可以给出下列代码：% % a = [1,2,3;4,5,6;7,8,9];% acol = bsxfun(@times,a,[1 2 4])function [centroids, labels] = run_kmeans(X, k, max_iter)% 该函数实现Kmeans聚类% 输入参数：%                   X为输入样本集，dxN%                   k为聚类中心个数%                   max_iter为kemans聚类的最大迭代的次数% 输出参数：%                   centroids为聚类中心 dxk%                   labels为样本的类别标记%% 采用K-means++算法初始化聚类中心% k-means++算法选择初始seeds的基本思想就是：初始的聚类中心之间的相互距离要尽可能的远。%% 1. 从输入的数据点集合中随机选择一个点作为第一个聚类中心% 2. 对于数据集中的每一个点x，计算它与最近聚类中心(指已选择的聚类中心)的距离D(x)% 3. 选择一个新的数据点作为新的聚类中心，选择的原则是：D(x)较大的点，被选取作为聚类中心的概率较大% 4. 重复2和3直到k个聚类中心被选出来% 5. 利用这k个初始的聚类中心来运行标准的k-means算法% % % 1. 先从我们的数据库随机挑个随机点当“种子点”% 2. 对于每个点，我们都计算其和最近的一个“种子点”的距离D(x)并保存在一个数组里，然后把这些距离加起来得到Sum(D(x))。% 3. 然后，再取一个随机值，用权重的方式来取计算下一个“种子点”。这个算法的实现是，先取一个能落在Sum(D(x))中的随机值%    Random，然后用Random -= D(x)，直到其<=0，此时的点就是下一个“种子点”。% 4. 重复2和3直到k个聚类中心被选出来% 5. 利用这k个初始的聚类中心来运行标准的k-means算法%% %第3步的核心思想：Sum(D(x))*random时，该值会以较大的概率落入D(x)较大的区间内，所以对应的点会以较大的概率被选中作为% 新的聚类中心。|-----5-----|----------10----------|----4----|  centroids = X(:, 1 + round(rand * (size(X, 2)-1))); % 从输入的数据点集合中随机选择一个点作为第一个聚类中心  labels = ones(1, size(X, 2)); % 初始化每个数据点的标号为1  for i = 2:k%         pp = centroids(:, labels);        D = X - centroids(:, labels); % 计算每一个点与聚类中心的差值        D = sqrt(dot(D, D, 1));% 每一点与聚类中心的差值的二范数        D = cumsum(D); % 累积求和，便于查找下一个种子点        if D(end) == 0,             centroids(:, i:k) = X(:, ones(1, k-i+1));             return;         end        nextCentroid = find(rand < D / D(end), 1); % 查找下一聚类中心点的坐标，Random -= D(x)，直到其<=0        centroids(:, i) = X(:, nextCentroid); % 添加新的聚类中心%         a = 2*real(centroids' * X);%         b = dot(centroids, centroids, 1).';        [~, labels] = max(bsxfun(@minus, 2*real(centroids' * X), dot(centroids, centroids, 1).'));%举例当前中心点为  end%% 标准Kmeans算法  for iter = 1:max_iter        for i = 1:k,             l = labels==i;             centroids(:,i) = sum(X(:, l), 2) / sum(l);         end        [~, labels] = max(bsxfun(@minus, 2*real(centroids' * X), dot(centroids,centroids,1).'), [], 1);  endend

(3)源代码调用示例

该段调用代码是将一副RGB图像转化在HSV颜色模型空间下进行聚类，输出为染色聚类结果

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%功能：利用kmeans算法对rgb图像进行聚类，并显示聚类染色结果%kmeans算法中当聚类的中心距离过小，或者某一类数量过低，可以将这些类合并到旁边类%调用run_kmeans算法%输入参数：%         X：输入样本数据集%         k: 聚类中心个数%         max_iter:kmeans最大迭代次数%输出参数：%        centroids:类中心的坐标 维数与样本集一致%        labels: （1：k）被标记好的分类输出%rgb2ycbcrclear all;close all;clc;%%读图% image=double(imread('lena.bmp'));% % figure(1);% % imshow(uint8(image)),title('Input image');% %%存储图像尺寸% image_rows=size(image,1);% image_cols=size(image,2);%%% %把图片存成样本X% X=(reshape(image,[image_rows*image_cols,3])');% % scatter(X(1,:),'r','filled');% % figure();% %scatter3(X(1,:)', X(2,:)', X(3,:)', c,'filled'),view(0,255);  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%rgb2hsv%读一张rgb图像，在hsv模型里面处理图像%功能：rgb图像转换成HSV模型图像，弄明白HSV图像在matlab中的存储方式，并reshape成图像进行输入clear all;close all;clc;% H=double(I);%H为RGB图像% [hue,s,v]=rgb2hsv(H);image=double(imread('lena.bmp'));image_rows=size(image,1);image_cols=size(image,2);[hue,s,v]=rgb2hsv(image);figure;imshow(hue),title('hue');figure;imshow(s),title('s');figure;imshow(uint8(v)),title('v');% image_hsv=rgb2hsv(image);hue=reshape(hue,1,image_rows*image_cols);s=reshape(s,1,image_rows*image_cols);v=reshape(v,1,image_rows*image_cols);X=[hue;s;v];%%%对样本空间进行聚类k=3;%聚类个数max_iter=100;%最大迭代次数[centroids, labels] = run_kmeans(X, k, max_iter); %kmeans把图像聚成3类% labeled_X=[X;labels];%%显示染色结果% if labeled_X(4,:)==1%     out_image=image(:,:,1)% elseif labeled_X(4,:)==2%         out_image=image(:,:,2)% elseif labeled_X(4,:)==3%     out_image=image(:,:,3)% end;pixels_labeled=reshape(labels,image_rows,image_cols);image_r=zeros(image_rows,image_cols);image_r(pixels_labeled==1)=255;% reshape(image_r,image_rows,image_cols);image_g=zeros(image_rows,image_cols);image_g(pixels_labeled==2)=255;% reshape(image_g,image_rows,image_cols);image_b=zeros(image_rows,image_cols);image_b(pixels_labeled==3)=255;% reshape(image_b,image_rows,image_cols);out_image=zeros(size(image));out_image(:,:,1)=image_r;out_image(:,:,2)=image_g;out_image(:,:,3)=image_b;figure(2)imshow(uint8(out_image));