clear;
% Sample = [10 + 10*rand(20,1) 10 + 10*rand(20,1)
% 10 + 10*rand(20,1) 10 - 10*rand(20,1)
% 10 - 10*rand(20,1) 10 + 10*rand(20,1)
% 10 - 10*rand(20,1) 10 - 10*rand(20,1)%
% 50 + 10*rand(5,1) 30 + 10*rand(5,1)
% 50 + 10*rand(5,1) 30 - 10*rand(5,1)
% 50 - 10*rand(5,1) 30 + 10*rand(5,1)
% 50 - 10*rand(5,1) 30 - 10*rand(5,1)%
% 10 + 10*rand(1,1) 50 + 10*rand(1,1)
% 10 + 10*rand(1,1) 50 - 10*rand(1,1)
% 10 - 10*rand(1,1) 50 + 10*rand(1,1)
% 10 - 10*rand(1,1) 50 - 10*rand(1,1)%
% ];
Sample = [12.807 16.872
10.219 16.307
14.33 16.837
16.551 16.3
14.869 19.618
13.31 10.712
11.684 16.057
17.708 15.642
12.424 18.115
10.505 17.92
10.663 16.198
16.665 18.592
12.116 13.647
12.978 16.136
18.958 10.925
10.958 19.101
19.464 14.904
12.083 10.118
17.566 18.728
14.24 15.715
14.119 1.2928
12.787 1.2363
15.636 7.5374
18.502 2.397
11.74 5.0391
10.744 8.2153
18.544 7.2736
15.147 4.3867
18.577 2.7317
10.358 5.0572
15.841 3.0786
10.812 7.1301
12.725 4.5619
17.584 8.0403
19.557 5.564
10.255 9.5455
11.751 7.2525
12.431 5.4754
16.614 8.3663
14.566 7.8058
1.0642 11.901
7.384 18.152
1.4097 17.176
7.6972 16.134
1.4571 14.624
6.9975 18.817
4.7723 10.235
9.2443 14.239
7.1532 10.142
8.443 16.591
2.9151 13.132
2.5722 16.83
2.8663 15.764
4.4414 15.37
8.3199 12.694
6.2803 10.935
6.3535 14.389
2.195 11.53
2.5685 12.091
8.3686 16.813
9.583 4.4095
2.01 7.9818
9.6902 2.0053
8.1312 3.7911
4.5517 7.761
6.0212 0.3755
7.9374 9.1311
2.1714 7.9006
1.5215 7.8467
3.9348 1.4959
3.0766 7.3091
8.0265 6.2431
7.0844 6.4797
0.9018 3.7683
3.7333 3.3184
7.7486 7.0959
1.0079 2.8365
3.051 2.8744
3.0352 3.0705
7.2478 3.257
53.7 34.891
51.355 33.999
53.644 36.623
54.035 33.214
53.549 31.827
50.683 26.372
56.123 22.695
56.357 21.741
51.085 22.445
51.643 28.57
46.435 31.286
41.048 32.852
41.445 36.567
48.012 33.297
43.445 39.583
42.471 21.859
40.765 21.814
48.837 29.574
47.299 20.363
48.535 27.975
11.708 53.722
12.146 48.595
6.8044 53.062
9.493 46.571
];
% Sample = [12.358 14.974
% 14.945 10.986
% 16.853 15.544
% 19.446 14.357
% 14.393 11.813
% 11.538 0.38958
% 19.337 2.245
% 11.363 9.9928
% 19.985 6.7453
% 16.906 4.559
% 6.9349 10.195
% 3.273 18.128
% 6.5429 19.828
% 0.19401 13.118
% 6.4603 13.789
% 2.5197 8.6545
% 5.5717 8.687
% 9.0011 1.3601
% 0.073775 2.9173
% 4.6374 0.017593
% 58.158 32.117
% 57.822 38.511
% 57.636 37.43
% 59.05 34.051
% 51.031 33.05
% 54.742 27.725
% 51.447 20.65
% 54.943 24.823
% 55.301 25.407
% 55.102 22.233
% 42.662 38.763
% 43.976 34.275
% 49.431 37.084
% 49.28 30.702
% 43.643 30.191
% 44.917 22.856
% 41.306 27.85
% 42.015 22.135
% 42.008 28.279
% 48.039 27.939
% 14.645 53.075
% 10.745 54.619
% 15.966 58.298
% 11.91 53.001
% 19.348 57.096
% 10.746 49.458
% 17.502 42.616
% 14.689 43.644
% 19.449 40.69
% 10.14 45.39
% 4.8711 58.23
% 5.7552 59.139
% 5.8408 58.36
% 6.1177 58.681
% 1.4446 57.756
% 5.6807 41.614
% 5.639 47.138
% 7.757 47.159
% 9.8682 46.729
% 6.2403 45.238
% ];
n = length(Sample);%样本总数
c = 3; %分类数
N = zeros(c);
centerx = zeros(c);
centery = zeros(c);
while(N(1)*N(2)*N(3) == 0) %若后面聚类的类数少于c，则重新选择代表点
N = zeros(c);
class1 = [1 1];
class2 = [1 1];
class3 = [1 1];
sumx = zeros(c);
sumy = zeros(c);
%将全部样本随机分成k类，计算每类重心，把这些重心作为每类的代表点
for i = 1:n
randnum = ceil(c*rand(1,1));
if(randnum == 1)
class1(N(1)+1,:) = Sample(i,:);
sumx(1)= sumx(1) + Sample(i,1);
sumy(1) = sumy(1) + Sample(i,2);
N(1) = N(1) + 1;
end
if(randnum == 2)
class2(N(2)+1,:) = Sample(i,:);
sumx(2) = sumx(2) + Sample(i,1);
sumy(2) = sumy(2) + Sample(i,2);
N(2) = N(2) + 1;
end
if(randnum == 3)
class3(N(3)+1,:) = Sample(i,:);
sumx(3) = sumx(3) + Sample(i,1);
sumy(3) = sumy(3) + Sample(i,2);
N(3) = N(3) + 1;
end
end
for i = 1:c
centerx(i) = sumx(i)/N(i);
centery(i) = sumy(i)/N(i);
end
% Sx = Sample(:,1);
% Sy = Sample(:,2);
% plot(Sx,Sy,'r.');
class = [class1
class2
class3];
class1 = [1 1];
class2 = [1 1];
class3 = [1 1];
N = zeros(c);
%开始以上面找到的代表点聚类
%计算每个样本与聚合中心的距离
for i = 1:n
d1 = (Sample(i,1) - centerx(1))^2 + (Sample(i,2) - centery(1))^2;
d2 = (Sample(i,1) - centerx(2))^2 + (Sample(i,2) - centery(2))^2;
d3 = (Sample(i,1) - centerx(3))^2 + (Sample(i,2) - centery(3))^2;
if(d1 <= d2 && d1 <= d3)
class1(N(1)+1,:) = Sample(i,:);
N(1) = N(1) + 1;
end
if(d2 <= d1 && d2 <= d3)
class2(N(2)+1,:) = Sample(i,:);
N(2) = N(2) + 1;
end
if(d3 <= d1 && d3 <= d2)
class3(N(3)+1,:) = Sample(i,:);
N(3) = N(3) + 1;
end
end
end %while(n1*n2*n3 == 0)
C1x = class1(:,1);
C1y = class1(:,2);
plot(C1x,C1y,'rs');
hold on;
plot(centerx(1),centery(1),'bs');
C2x = class2(:,1);
C2y = class2(:,2);
plot(C2x,C2y,'ro');
plot(centerx(2),centery(2),'bo');
C3x = class3(:,1);
C3y = class3(:,2);
plot(C3x,C3y,'r^');
plot(centerx(3),centery(3),'b^');
hold off;
I = 1;%迭代次数
sumx = zeros(c);
sumy = zeros(c);
%计算新的聚合中心
for i = 1:N(1)
sumx(1) = sumx(1) + class1(i,1);
sumy(1) = sumy(1) + class1(i,2);
end
for i = 1:N(2)
sumx(2) = sumx(2) + class2(i,1);
sumy(2) = sumy(2) + class2(i,2);
end
for i = 1:N(3)
sumx(3) = sumx(3) + class3(i,1);
sumy(3) = sumy(3) + class3(i,2);
end
for i = 1:c
centerx(i) = sumx(i)/N(i);
centery(i) = sumy(i)/N(i);
end
C1x = class1(:,1);
C1y = class1(:,2);
plot(C1x,C1y,'rs');
hold on;
plot(centerx(1),centery(1),'bs');
C2x = class2(:,1);
C2y = class2(:,2);
plot(C2x,C2y,'ro');
plot(centerx(2),centery(2),'bo');
C3x = class3(:,1);
C3y = class3(:,2);
plot(C3x,C3y,'r^');
plot(centerx(3),centery(3),'b^');
hold off;
Jc = Err(class1,centerx(1),centery(1),N(1)) + Err(class2,centerx(2),centery(2),N(2)) + Err(class3,centerx(3),centery(3),N(3));
flag = 1;%标记有无样本在类间移动标志
count = 0;
while(flag)
flag = 0;
%class 1
for i = 1:N(1)
[N(1) a] = size(class1);
if(N(1) > 1 && i <= N(1))
rou11 = N(1)*Distance(class1(i,:),[centerx(1),centery(1)])/(N(1)-1);
rou12 = N(2)*Distance(class1(i,:),[centerx(2),centery(2)])/(N(2)+1);
rou13 = N(3)*Distance(class1(i,:),[centerx(3),centery(3)])/(N(3)+1);
if(rou12 < rou13)
if(rou12 < rou11)
ce = [centerx(1) centery(1)] + ([centerx(1),centery(1)] - class1(i,:))/(N(1)-1);
centerx(1) = ce(1);
centery(1) = ce(2);
ce = [centerx(2),centery(2)] - ([centerx(2),centery(2)] - class1(i,:))/(N(2)+1);
centerx(2) = ce(1);
centery(2) = ce(2);
class2(N(2)+1,:) = class1(i,:);
%从class1中删除
class1 = delet(class1,i);
N(2) = N(2) + 1;
[N(1) a] = size(class1);
flag = 1;
count = count + 1;
end
else
if(rou13 < rou11)
ce = [centerx(1),centery(1)] + ([centerx(1),centery(1)] - class1(i,:))/(N(1)-1);
centerx(1) = ce(1);
centery(1) = ce(2);
ce = [centerx(3),centery(3)] - ([centerx(3),centery(3)] - class1(i,:))/(N(3)+1);
