K-Means聚类算法原理及其python和matlab实现

4 篇文章 2 订阅
1 篇文章 0 订阅

(一)何谓聚类

  还是那句“物以类聚、人以群分”,如果预先知道人群的标签(如文艺、普通、2B),那么根据监督学习的分类算法可将一个人明确的划分到某一类;如果预先不知道人群的标签,那就只有根据人的特征(如爱好、学历、职业等)划堆了,这就是聚类算法。
  聚类是一种无监督的学习(无监督学习不依赖预先定义的类或带类标记的训练实例),它将相似的对象归到同一个簇中,它是观察式学习,而非示例式的学习,有点像全自动分类。所谓簇就是该集合中的对象有很大的相似性,而不同集合间的对象有很大的相异性。簇识别(cluster identification)给出了聚类结果的含义,告诉我们这些簇到底都是些什么。通常情况下,簇质心可以代表整个簇的数据来做出决策。聚类方法几乎可以应用于所有对象,簇内的对象越相似,聚类的效果越好。  
  说白了,聚类(clustering)是完全可以按字面意思来理解的——将相同、相似、相近、相关的对象实例聚成一类的过程。简单理解,如果一个数据集合包含N个实例,根据某种准则可以将这N个实例划分为m个类别,每个类别中的实例都是相关的,而不同类别之间是区别的也就是不相关的,这就得到了一个聚类模型了。判别新样本点的所属类时,就通过计算该点与这m个类别的相似度,选择最相似的那个类作为该点的归类。
  机器学习中常见的聚类算法包括 k-Means算法、期望最大化算法(Expectation Maximization,EM,参考“EM算法原理”)、谱聚类算法(参考机器学习算法复习-谱聚类)以及人工神经网络算法,本文介绍K-均值(K-means)聚类算法。

(二)K-均值聚类算法

1. 认识K-均值聚类算法

  在聚类分析中,K-均值聚类算法(k-means algorithm)是无监督分类中的一种基本方法,其也称为C-均值算法,其基本思想是:通过迭代的方法,逐次更新各聚类中心的值,直至得到最好的聚类结果。
  k-means算法的基础是最小误差平方和准则。其代价函数是:


这里写图片描述

  式中,μc(i)表示第i个簇的质心,我们希望得到的聚类模型代价函数最小,直观的来说,各簇内的样本越相似,其与该簇质心的误差平方越小。计算所有簇的误差平方之和,即可验证分为k个簇时时的聚类是否是最优的。SSE值越小表示数据点越接近于它们的质心,聚类效果也越好。因为对误差取了平方,因此更加重视那些远离中心的点。一种肯定可以降低SSE值的方法是增加簇的个数,但这违背了聚类的目标,聚类的目标是在保持族数目不变的情况下提高簇的质量。
  k-均值(k-means)聚类算法之所以称之为k-均值是因为它可以发现k个不同的簇,且每个簇的中心采用簇中所含子数据集样本特征的均值计算而成。k-均值是发现给定数据集的k个簇的算法,簇个数k由用户给定,每一个簇通过其质心( centroid) — 即簇中所有点的中心来描述。K-均值聚类算法需要数值型数据来进行相似性度量,也可以将标称型数据映射为二值型数据再用于度量相似性,其优点是容易实现,缺点是可能收敛到局部最小值,在大规模数据集上收敛较慢。
  假设训练样本数据集X为(m, n)维矩阵,m表示样本数、n表示每个样本点的特征数,那么k-均值聚类算法的结果就是得到一个kxn维矩阵,k表示用户指定的簇个数,每一行都是一个长度为n的行向量–第i个元素就是该簇中所有样本第j(j=0,1,…,n-1)个特征的均值。下图是K-均值聚类算法聚类的过程:

这里写图片描述

2. 算法过程

假设要把样本集分为c个类别,算法如下:
(1)适当选择c个类的初始中心;
(2)在第k次迭代中,对任意一个样本,求其到c个中心的距离,将该样本归到距离最短的中心所在的类,
(3)利用均值等方法更新该类的中心值;
(4)对于所有的c个聚类中心,如果利用(2)(3)的迭代法更新后,值保持不变(目标函数收敛),则迭代结束,否则继续迭代。

(三) K-均值算法的实现

实现步骤:

 1.数据集的赋值
 2.k均值算法的python实现
 3.k均值算法的matlab实现

1.数据集的赋值

  开始之前,我们可以看一下数据,如图:
  

这里写图片描述
  

  测试数据是二维的,共80个样本。有4个类。如下:
testSet.txt

1.658985    4.285136  
-3.453687   3.424321  
4.838138    -1.151539  
-5.379713   -3.362104  
0.972564    2.924086  
-3.567919   1.531611  
0.450614    -3.302219  
-3.487105   -1.724432  
2.668759    1.594842  
-3.156485   3.191137  
3.165506    -3.999838  
-2.786837   -3.099354  
4.208187    2.984927  
-2.123337   2.943366  
0.704199    -0.479481  
-0.392370   -3.963704  
2.831667    1.574018  
-0.790153   3.343144  
2.943496    -3.357075  
-3.195883   -2.283926  
2.336445    2.875106  
-1.786345   2.554248  
2.190101    -1.906020  
-3.403367   -2.778288  
1.778124    3.880832  
-1.688346   2.230267  
2.592976    -2.054368  
-4.007257   -3.207066  
2.257734    3.387564  
-2.679011   0.785119  
0.939512    -4.023563  
-3.674424   -2.261084  
2.046259    2.735279  
-3.189470   1.780269  
4.372646    -0.822248  
-2.579316   -3.497576  
1.889034    5.190400  
-0.798747   2.185588  
2.836520    -2.658556  
-3.837877   -3.253815  
2.096701    3.886007  
-2.709034   2.923887  
3.367037    -3.184789  
-2.121479   -4.232586  
2.329546    3.179764  
-3.284816   3.273099  
3.091414    -3.815232  
-3.762093   -2.432191  
3.542056    2.778832  
-1.736822   4.241041  
2.127073    -2.983680  
-4.323818   -3.938116  
3.792121    5.135768  
-4.786473   3.358547  
2.624081    -3.260715  
-4.009299   -2.978115  
2.493525    1.963710  
-2.513661   2.642162  
1.864375    -3.176309  
-3.171184   -3.572452  
2.894220    2.489128  
-2.562539   2.884438  
3.491078    -3.947487  
-2.565729   -2.012114  
3.332948    3.983102  
-1.616805   3.573188  
2.280615    -2.559444  
-2.651229   -3.103198  
2.321395    3.154987  
-1.685703   2.939697  
3.031012    -3.620252  
-4.599622   -2.185829  
4.196223    1.126677  
-2.133863   3.093686  
4.668892    -2.562705  
-2.793241   -2.149706  
2.884105    3.043438  
-2.967647   2.848696  
4.479332    -1.764772  
-4.905566   -2.911070 

2、k均值算法的python实现

kmeans.py

#################################################  
# kmeans: k-means cluster  
# Author : stu_why
# Date   : 2016-12-06
# HomePage : http://blog.csdn.net/zpp1994
# Email  : 1620009136@qq.com
#################################################  
from sklearn.cluster import KMeans
import numpy
import matplotlib.pyplot as plt  

# step 1: load data
print('step 1: load data...') 

#读取testSet.txt数据并存储到dataSet中
dataSet = []
fileIn = open('testSet.txt')
for line in fileIn.readlines():  
    lineArr = line.strip().split()
    dataSet.append('%0.6f' % float(lineArr[0]))
    dataSet.append('%0.6f' % float(lineArr[1]))

# step 2: clustering...
print('step 2: clustering...')

#调用sklearn.cluster中的KMeans类
dataSet = numpy.array(dataSet).reshape(80,2)
kmeans = KMeans(n_clusters=4, random_state=0).fit(dataSet)

#求出聚类中心
center=kmeans.cluster_centers_
center_x=[]
center_y=[]
for i in range(len(center)):
    center_x.append('%0.6f' % center[i][0])
    center_y.append('%0.6f' % center[i][1])

#标注每个点的聚类结果
labels=kmeans.labels_
type1_x = []
type1_y = []
type2_x = []
type2_y = []
type3_x = []
type3_y = []
type4_x = []
type4_y = []
for i in range(len(labels)):
    if labels[i] == 0:
        type1_x.append(dataSet[i][0])
        type1_y.append(dataSet[i][1])
    if labels[i] == 1:
        type2_x.append(dataSet[i][0])
        type2_y.append(dataSet[i][1])
    if labels[i] == 2:
        type3_x.append(dataSet[i][0])
        type3_y.append(dataSet[i][1])
    if labels[i] == 3:
        type4_x.append(dataSet[i][0])
        type4_y.append(dataSet[i][1])

#画出四类数据点及聚类中心
plt.figure(figsize=(10,8), dpi=80)
axes = plt.subplot(111)
type1 = axes.scatter(type1_x, type1_y, s=40, c='red')
type2 = axes.scatter(type2_x, type2_y, s=40, c='green')
type3 = axes.scatter(type3_x, type3_y,s=40, c='pink' )
type4 = axes.scatter(type4_x, type4_y, s=40, c='yellow')
type_center = axes.scatter(center_x, center_y, s=40, c='blue')
plt.xlabel('x')
plt.ylabel('y')

axes.legend((type1, type2, type3, type4,type_center), ('0','1','2','3','center'),loc=1)
plt.show()  

代码运行结果:

其中,不同的类用不同的颜色来表示,其中的蓝色圆圈是对应类的均值质心点。

3、k均值算法的matlab实现

用matlab编写k-均值聚类程序:
kmean.m

% k-均聚类算法
clc
clear;

% main variables
dim = 2; % 模式样本维数
k = 4;   % 设有k个聚类中心

load('testSet.txt');
PM=testSet;% 模式样本矩阵 
N = size(PM,1);
figure();
subplot(1,2,1);
for(i=1:N)
   plot(PM(i,1),PM(i,2), '*r'); % 绘出原始的数据点
   hold on
end
xlabel('X');
ylabel('Y');
title('聚类之前的数据点');

CC = zeros(k,dim); % 聚类中心矩阵,CC(i,:)初始值为i号样本向量
D = zeros(N,k); % D(i,j)是样本i和聚类中心j的距离

C = cell(1,k); %% 聚类矩阵,对应聚类包含的样本。初始状况下,聚类i(i<k)的样本集合为[i],聚类k的样本集合为[k,k+1,...N]
for i = 1:k-1
    C{i} = [i];
end
C{k} = k:N;

B = 1:N; % 上次迭代中,样本属于哪一聚类,设初值为1
B(k:N) = k;

for i = 1:k
    CC(i,:) = PM(i,:);
end

while 1   
    change = 0;%用来标记分类结果是否变化
    % 对每一个样本i,计算到k个聚类中心的距离
    for i = 1:N
        for j = 1:k
%             D(i,j) = eulerDis( PM(i,:), CC(j,:) );
              D(i,j) = sqrt((PM(i,1) - CC(j,1))^2 + (PM(i,2) - CC(j,2))^2);
        end
        t = find( D(i,:) == min(D(i,:)) ); % i属于第t类
        if B(i) ~= t % 上次迭代i不属于第t类
            change = 1;
            % 将i从第B(i)类中去掉
            t1 = C{B(i)};
            t2 = find( t1==i );            
            t1(t2) = t1(1);
            t1 = t1(2:length(t1)); 
            C{B(i)} = t1;

            C{t} = [C{t},i]; % 将i加入第t类

            B(i) = t;
        end        
    end

    if change == 0 %分类结果无变化,则迭代停止
        break;
    end

    % 重新计算聚类中心矩阵CC
    for i = 1:k
        CC(i,:) = 0;
        iclu = C{i};
        for j = 1:length(iclu)
            CC(i,:) = PM( iclu(j),: )+CC(i,:);
        end
        CC(i,:) = CC(i,:)/length(iclu);
    end
end

subplot(1,2,2);
plot(CC(:,1),CC(:,2),'o')
hold on
for(i=1:N)
    if(B(1,i)==1)
        plot(PM(i,1),PM(i,2),'*b'); %作出第一类点的图形
        hold on
    elseif(B(1,i)==2)
        plot(PM(i,1),PM(i,2), '*r'); %作出第二类点的图形
     hold on
     elseif(B(1,i)==3)
        plot(PM(i,1),PM(i,2),'*g'); %作出第三类点的图形
        hold on
    else
        plot(PM(i,1),PM(i,2), '*m'); %作出第四类点的图形
        hold on
    end
end
xlabel('X');
ylabel('Y');
title('聚类之后的数据点');
      % 打印C,CC     
for i = 1:k    %输出每一类的样本点标号
    str=['第' num2str(i) '类包含点:  ' num2str(C{i})];
    disp(str);
end;              

利用testSet.tex的数据,代码运行结果如下:

并在命令行输出每个点的类标记:

第1类包含点: 5 9 13 17 1 21 25 29 33 37 41 45 49 53 57 61 65 69 77 73
第2类包含点: 2 6 10 14 22 26 30 34 38 42 46 50 54 58 62 66 70 74 78 18
第3类包含点: 7 11 15 19 23 27 31 35 39 43 47 51 55 59 63 67 71 3 75 79
第4类包含点: 16 36 48 64 28 8 4 68 52 40 24 72 32 56 12 76 44 20 60 80

  • 48
    点赞
  • 339
    收藏
    觉得还不错? 一键收藏
  • 23
    评论
function [idx, C, sumD, D] = kmeans(X, k, varargin) % varargin:实际输入参量 if nargin 1 % 大于1刚至少有一种距离 error(sprintf('Ambiguous ''distance'' parameter value: %s.', distance)); elseif isempty(i) % 如果是空的,则表明没有合适的距离 error(sprintf('Unknown ''distance'' parameter value: %s.', distance)); end % 针对不同的距离,处理不同 distance = distNames{i}; switch distance case 'cityblock' % sort 列元素按升序排列,Xord中存的是元素在原始矩阵中的列中对应的大小位置 [Xsort,Xord] = sort(X,1); case 'cosine' % 余弦 % 计算每一行的和的平方根 Xnorm = sqrt(sum(X.^2, 2)); if any(min(Xnorm) <= eps * max(Xnorm)) error(['Some points have small relative magnitudes, making them ', ... 'effectively zero.\nEither remove those points, or choose a ', ... 'distance other than ''cosine''.'], []); end % 标量化 Xnorm(:,ones(1,p))得到n*p的矩阵 X = X ./ Xnorm(:,ones(1,p)); case 'correlation' % 线性化 X = X - repmat(mean(X,2),1,p); % 计算每一行的和的平方根 Xnorm = sqrt(sum(X.^2, 2)); if any(min(Xnorm) 1 error(sprintf('Ambiguous ''start'' parameter value: %s.', start)); elseif isempty(i) error(sprintf('Unknown ''start'' parameter value: %s.', start)); elseif isempty(k) error('You must specify the number of clusters, K.'); end start = startNames{i}; % strcmp比较两个字符串 uniform是在X中随机选择K个点 if strcmp(start, 'uniform') if strcmp(distance, 'hamming') error('Hamming distance cannot be initialized with uniform random values.'); end % 标量化后的X Xmins = min(X,1); Xmaxs = max(X,1); end elseif isnumeric(start) % 判断输入是否是一个数 这里的start是一个K*P的矩阵,表示初始聚类中心 CC = start; % CC表示初始聚类中心 start = 'numeric'; if isempty(k) k = size(CC,1); elseif k ~= size(CC,1); error('The ''start'' matrix must have K rows.'); elseif size(CC,2) ~= p error('The ''start'' matrix must have the same number of columns as X.'); end if isempty(reps) reps = size(CC,3); elseif reps ~= size(CC,3); error('The third dimension of the ''start'' array must match the ''replicates'' parameter value.'); end % Need to center explicit starting points for 'correlation'. % 线性距离需要指定具体的开始点 % (Re)normalization for 'cosine'/'correlation' is done at each % iteration.每一次迭代进行“余弦和线性”距离正规化 % 判断是否相等 if isequal(distance, 'correlation') % repmat复制矩阵1*P*1 线性化 CC = CC - repmat(mean(CC,2),[1,p,1]); end else error('The ''start'' parameter value must be a string or a numeric matrix or array.'); end % ------------------------------------------------------------------ % 如果一个聚类丢失了所有成员,这个进程将被调用 if ischar(emptyact) emptyactNames = {'error','drop','singleton'}; i = strmatch(lower(emptyact), emptyactNames); if length(i) > 1 error(sprintf('Ambiguous ''emptyaction'' parameter value: %s.', emptyact)); elseif isempty(i) error(sprintf('Unknown ''emptyaction'' parameter value: %s.', emptyact)); end emptyact = emptyactNames{i}; else error('The ''emptyaction'' parameter value must be a string.'); end % ------------------------------------------------------------------ % 控制输出的显示示信息 if ischar(display) % strvcat 垂直连接字符串 i = strmatch(lower(display), strvcat('off','notify','final','iter')); if length(i) > 1 error(sprintf('Ambiguous ''display'' parameter value: %s.', display)); elseif isempty(i) error(sprintf('Unknown ''display'' parameter value: %s.', display)); end display = i-1; else error('The ''display'' parameter value must be a string.'); end % ------------------------------------------------------------------ % 输入参数K的控制 if k == 1 error('The number of clusters must be greater than 1.'); elseif n 2 % 'iter',\t 表示向后空出8个字符 disp(sprintf(' iter\t phase\t num\t sum')); end % ------------------------------------------------------------------ % Begin phase one: batch reassignments 第一队段:分批赋值 converged = false; iter = 0; while true % Compute the distance from every point to each cluster centroid % 计算每一个点到每一个聚类中心的距离,D中存放的是N*K个数值 D(:,changed) = distfun(X, C(changed,:), distance, iter); % Compute the total sum of distances for the current configuration. % Can't do it first time through, there's no configuration yet. % 计算当前配置的总距离,但第一次不能执行,因为还没有配置 if iter > 0 totsumD = sum(D((idx-1)*n + (1:n)')); % Test for a cycle: if objective is not decreased, back out % the last step and move on to the single update phase % 循环测试:如果目标没有减少,退出到最后一步,移动到第二阶段 % prevtotsumD表示前一次的总距离,如果前一次的距离比当前的小,则聚类为上一次的,即不发生变化 if prevtotsumD 2 % 'iter' disp(sprintf(dispfmt,iter,1,length(moved),totsumD)); end if iter >= maxit, % break(2) break; % 跳出while end end % Determine closest cluster for each point and reassign points to clusters % 决定每一个点的最近分簇,重新分配点到各个簇 previdx = idx; % 大小为n*1 % totsumD 被初始化为无穷大,这里表示总距离 prevtotsumD = totsumD; % 返回每一行中最小的元素,d的大小为n*1,nidx为最小元素在行中的位置,其大小为n*1,D为n*p [d, nidx] = min(D, [], 2); if iter == 0 % iter==0,表示第一次迭代 % Every point moved, every cluster will need an update % 每一个点需要移动,每一个簇更新 moved = 1:n; idx = nidx; changed = 1:k; else % Determine which points moved 决定哪一个点移动 % 找到上一次和当前最小元素不同的位置 moved = find(nidx ~= previdx); if length(moved) > 0 % Resolve ties in favor of not moving % 重新分配而不是移动 括号中是一个逻辑运算 确定需要移动点的位置 moved = moved(D((previdx(moved)-1)*n + moved) > d(moved)); end % 如果没有不同的,即当前的是最小元素,跳出循环,得到的元素已经是各行的最小值 if length(moved) == 0 % break(3) break; end idx(moved) = nidx(moved); % Find clusters that gained or lost members 找到获得的或者丢失的成员的分簇 % 得到idx(moved)和previdx(moved)中不重复出现的所有元素,并按升序排列 changed = unique([idx(moved); previdx(moved)])'; end % Calculate the new cluster centroids and counts. 计算新的分簇中心和计数 % C(changed,:)表示新的聚类中心,m(changed)表示聚类标号在idx中出现的次数 % sort 列元素按升序排列,Xsort存放对的元素,Xord中存的是元素在原始矩阵中的列中对应的大小位置 [C(changed,:), m(changed)] = gcentroids(X, idx, changed, distance, Xsort, Xord); iter = iter + 1; % Deal with clusters that have just lost all their members % 处理丢失所有成员的分簇,empties表示丢失元素的聚类标号 不用考虑 empties = changed(m(changed) == 0); if ~isempty(empties) switch emptyact case 'error' % 默认值,一般情况下不会出现 error(sprintf('Empty cluster created at iteration %d.',iter)); case 'drop' % Remove the empty cluster from any further processing % 移走空的聚类 D(:,empties) = NaN; changed = changed(m(changed) > 0); if display > 0 warning(sprintf('Empty cluster created at iteration %d.',iter)); end case 'singleton' if display > 0 warning(sprintf('Empty cluster created at iteration %d.',iter)); end for i = empties % Find the point furthest away from its current cluster. % Take that point out of its cluster and use it to create % a new singleton cluster to replace the empty one. % 找到离聚类中心最远距离的点,把该点从它的聚类中移走,用它来创建一个新的聚类,来代替空的聚类 % 得到列的最大元素(dlarge),以及该元素在列中的标号(lonely) [dlarge, lonely] = max(d); from = idx(lonely); % taking from this cluster 从当前聚类移走 C(i,:) = X(lonely,:); m(i) = 1; idx(lonely) = i; d(lonely) = 0; % Update clusters from which points are taken % 更新那些点被移走的聚类 [C(from,:), m(from)] = gcentroids(X, idx, from, distance, Xsort, Xord); changed = unique([changed from]); end end end end % phase one % ------------------------------------------------------------------ % Initialize some cluster information prior to phase two % 为第二阶段初始化一些先验聚类信息 针对特定的距离,默认的是欧氏距离 switch distance case 'cityblock' Xmid = zeros([k,p,2]); for i = 1:k if m(i) > 0 % Separate out sorted coords for points in i'th cluster, % and save values above and below median, component-wise % 分解出第i个聚类中挑选的点的坐标,保存它的上,下中位数 % reshape把矩阵分解为要求的行列数m*p % sort 列元素按升序排列,Xord中存的是元素在原始矩阵中的列中对应的大小位置 Xsorted = reshape(Xsort(idx(Xord)==i), m(i), p); % floor取比值小或者等于的最近的值 nn = floor(.5*m(i)); if mod(m(i),2) == 0 Xmid(i,:,1:2) = Xsorted([nn, nn+1],:)'; elseif m(i) > 1 Xmid(i,:,1:2) = Xsorted([nn, nn+2],:)'; else Xmid(i,:,1:2) = Xsorted([1, 1],:)'; end end end case 'hamming' Xsum = zeros(k,p); for i = 1:k if m(i) > 0 % Sum coords for points in i'th cluster, component-wise % 对基于分量对第i个聚类的坐标点求和 Xsum(i,:) = sum(X(idx==i,:), 1); end end end % ------------------------------------------------------------------ % Begin phase two: single reassignments 第二阶段:单独赋值 % m中保存的是每一个聚类的个数,元素和为n % find(m' > 0)得到m'中大于0的元素的位置(索引) % 实际情况(默认情况下)changed=1:k changed = find(m' > 0); lastmoved = 0; nummoved = 0; iter1 = iter; while iter < maxit % Calculate distances to each cluster from each point, and the % potential change in total sum of errors for adding or removing % each point from each cluster. Clusters that have not changed % membership need not be updated. % 计算从每一个点到每一个聚类的距离,潜在由于总距离发生变化移除或添加引起的错误。 % 那些成员没有改变的聚类不需要更新。 % % Singleton clusters are a special case for the sum of dists % calculation. Removing their only point is never best, so the % reassignment criterion had better guarantee that a singleton % point will stay in its own cluster. Happily, we get % Del(i,idx(i)) == 0 automatically for them. % 单独聚类在计算距离时是一个特殊情况,仅仅移除点不是最好的方法,因此重新赋值准则能够保证一个 % 单独的点能够留在它的聚类中,可喜的是,自动有Del(i,idx(i)) == 0。 switch distance case 'sqeuclidean' for i = changed % idx存放的距离矩阵D中每一行的最小元素所处的列,也即位置 mbrs = (idx == i); sgn = 1 - 2*mbrs; % -1 for members, 1 for nonmembers % 表示只有一个聚类 if m(i) == 1 % prevent divide-by-zero for singleton mbrs 防止单个成员做0处理 sgn(mbrs) = 0; end Del(:,i) = (m(i) ./ (m(i) + sgn)) .* sum((X - C(repmat(i,n,1),:)).^2, 2); end case 'cityblock' for i = changed if mod(m(i),2) == 0 % this will never catch singleton clusters ldist = Xmid(repmat(i,n,1),:,1) - X; rdist = X - Xmid(repmat(i,n,1),:,2); mbrs = (idx == i); sgn = repmat(1-2*mbrs, 1, p); % -1 for members, 1 for nonmembers Del(:,i) = sum(max(0, max(sgn.*rdist, sgn.*ldist)), 2); else Del(:,i) = sum(abs(X - C(repmat(i,n,1),:)), 2); end end case {'cosine','correlation'} % The points are normalized, centroids are not, so normalize them normC(changed) = sqrt(sum(C(changed,:).^2, 2)); if any(normC 0 % Resolve ties in favor of not moving % 重新分配而不是移动 确定移动的位置 moved = moved(Del((previdx(moved)-1)*n + moved) > minDel(moved)); end if length(moved) == 0 % Count an iteration if phase 2 did nothing at all, or if we're % in the middle of a pass through all the points if (iter - iter1) == 0 | nummoved > 0 iter = iter + 1; if display > 2 % 'iter' disp(sprintf(dispfmt,iter,2,nummoved,totsumD)); end end converged = true; break; end % Pick the next move in cyclic order moved = mod(min(mod(moved - lastmoved - 1, n) + lastmoved), n) + 1; % If we've gone once through all the points, that's an iteration if moved 2 % 'iter' disp(sprintf(dispfmt,iter,2,nummoved,totsumD)); end if iter >= maxit, break; end nummoved = 0; end nummoved = nummoved + 1; lastmoved = moved; oidx = idx(moved); nidx = nidx(moved); totsumD = totsumD + Del(moved,nidx) - Del(moved,oidx); % Update the cluster index vector, and rhe old and new cluster % counts and centroids idx(moved) = nidx; m(nidx) = m(nidx) + 1; m(oidx) = m(oidx) - 1; switch distance case 'sqeuclidean' C(nidx,:) = C(nidx,:) + (X(moved,:) - C(nidx,:)) / m(nidx); C(oidx,:) = C(oidx,:) - (X(moved,:) - C(oidx,:)) / m(oidx); case 'cityblock' for i = [oidx nidx] % Separate out sorted coords for points in each cluster. % New centroid is the coord median, save values above and % below median. All done component-wise. Xsorted = reshape(Xsort(idx(Xord)==i), m(i), p); nn = floor(.5*m(i)); if mod(m(i),2) == 0 C(i,:) = .5 * (Xsorted(nn,:) + Xsorted(nn+1,:)); Xmid(i,:,1:2) = Xsorted([nn, nn+1],:)'; else C(i,:) = Xsorted(nn+1,:); if m(i) > 1 Xmid(i,:,1:2) = Xsorted([nn, nn+2],:)'; else Xmid(i,:,1:2) = Xsorted([1, 1],:)'; end end end case {'cosine','correlation'} C(nidx,:) = C(nidx,:) + (X(moved,:) - C(nidx,:)) / m(nidx); C(oidx,:) = C(oidx,:) - (X(moved,:) - C(oidx,:)) / m(oidx); case 'hamming' % Update summed coords for points in each cluster. New % centroid is the coord median. All done component-wise. Xsum(nidx,:) = Xsum(nidx,:) + X(moved,:); Xsum(oidx,:) = Xsum(oidx,:) - X(moved,:); C(nidx,:) = .5*sign(2*Xsum(nidx,:) - m(nidx)) + .5; C(oidx,:) = .5*sign(2*Xsum(oidx,:) - m(oidx)) + .5; end changed = sort([oidx nidx]); end % phase two % ------------------------------------------------------------------ if (~converged) & (display > 0) warning(sprintf('Failed to converge in %d iterations.', maxit)); end % Calculate cluster-wise sums of distances nonempties = find(m(:)'>0); D(:,nonempties) = distfun(X, C(nonempties,:), distance, iter); d = D((idx-1)*n + (1:n)'); sumD = zeros(k,1); for i = 1:k sumD(i) = sum(d(idx == i)); end if display > 1 % 'final' or 'iter' disp(sprintf('%d iterations, total sum of distances = %g',iter,totsumD)); end % Save the best solution so far if totsumD 3 Dbest = D; end end end % Return the best solution
评论 23
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值