mainpos.m

最基本的粒子群算法实现的数据聚类

 %------基本粒子群优化算法（Particle Swarm Optimization）-----------
%------名称：基本粒子群优化算法（PSO）
%------作用：求解优化问题
%------说明：全局性，并行性，高效的群体智能算法
%------作者：孙明杰（[e-mail]）
%------单位：中国矿业大学理学院计算数学硕2005
%------时间：2006年8月17日 
%------------------------------------------------------------------
%------初始格式化--------------------------------------------------
clear all;
%clc;
format long;
%------给定初始化条件----------------------------------------------
c1=1.5;      %   =1.4962    %学习因子1
c2=2;      % 1.4962    %学习因子2
w=0.8;%.7298;              %惯性权重
MaxDT=10;            %最大迭代次数
D=4;                  %搜索空间维数（未知数个数）
N=30;                  %初始化群体个体数目
eps=10^(-6);           %设置精度(在已知最小值时候用)
Sample=150;     %聚类的样本数
K=3;                    %聚类数目
y=zeros(N,K*D);
vbest=zeros(N,K*D);  %记录最好速度的
q=0.14;

s=[5.1,3.5,1.4,0.2;
4.9,3.0,1.4,0.2;
4.7,3.2,1.3,0.2;
4.6,3.1,1.5,0.2;
5.0,3.6,1.4,0.2;
5.4,3.9,1.7,0.4;
4.6,3.4,1.4,0.3;
5.0,3.4,1.5,0.2;
4.4,2.9,1.4,0.2;
4.9,3.1,1.5,0.1;
5.4,3.7,1.5,0.2;
4.8,3.4,1.6,0.2;
4.8,3.0,1.4,0.1;
4.3,3.0,1.1,0.1;
5.8,4.0,1.2,0.2;
5.7,4.4,1.5,0.4;
5.4,3.9,1.3,0.4;
5.1,3.5,1.4,0.3;
5.7,3.8,1.7,0.3;
5.1,3.8,1.5,0.3;
5.4,3.4,1.7,0.2;
5.1,3.7,1.5,0.4;
4.6,3.6,1.0,0.2;
5.1,3.3,1.7,0.5;
4.8,3.4,1.9,0.2;
5.0,3.0,1.6,0.2;
5.0,3.4,1.6,0.4;
5.2,3.5,1.5,0.2;
5.2,3.4,1.4,0.2;
4.7,3.2,1.6,0.2;
4.8,3.1,1.6,0.2;
5.4,3.4,1.5,0.4;
5.2,4.1,1.5,0.1;
5.5,4.2,1.4,0.2;
4.9,3.1,1.5,0.1;
5.0,3.2,1.2,0.2;
5.5,3.5,1.3,0.2;
4.9,3.1,1.5,0.1;
4.4,3.0,1.3,0.2;
5.1,3.4,1.5,0.2;
5.0,3.5,1.3,0.3;
4.5,2.3,1.3,0.3;
4.4,3.2,1.3,0.2;
5.0,3.5,1.6,0.6;
5.1,3.8,1.9,0.4;
4.8,3.0,1.4,0.3;
5.1,3.8,1.6,0.2;
4.6,3.2,1.4,0.2;
5.3,3.7,1.5,0.2;
5.0,3.3,1.4,0.2;
7.0,3.2,4.7,1.4;
6.4,3.2,4.5,1.5;
6.9,3.1,4.9,1.5;
5.5,2.3,4.0,1.3;
6.5,2.8,4.6,1.5;
5.7,2.8,4.5,1.3;
6.3,3.3,4.7,1.6;
4.9,2.4,3.3,1.0;
6.6,2.9,4.6,1.3;
5.2,2.7,3.9,1.4;
5.0,2.0,3.5,1.0;
5.9,3.0,4.2,1.5;
6.0,2.2,4.0,1.0;
6.1,2.9,4.7,1.4;
5.6,2.9,3.6,1.3;
6.7,3.1,4.4,1.4;
5.6,3.0,4.5,1.5;
5.8,2.7,4.1,1.0;
6.2,2.2,4.5,1.5;
5.6,2.5,3.9,1.1;
5.9,3.2,4.8,1.8;
6.1,2.8,4.0,1.3;
6.3,2.5,4.9,1.5;
6.1,2.8,4.7,1.2;
6.4,2.9,4.3,1.3;
6.6,3.0,4.4,1.4;
6.8,2.8,4.8,1.4;
6.7,3.0,5.0,1.7;
6.0,2.9,4.5,1.5;
5.7,2.6,3.5,1.0;
5.5,2.4,3.8,1.1;
5.5,2.4,3.7,1.0;
5.8,2.7,3.9,1.2;
6.0,2.7,5.1,1.6;
5.4,3.0,4.5,1.5;
6.0,3.4,4.5,1.6;
6.7,3.1,4.7,1.5;
6.3,2.3,4.4,1.3;
5.6,3.0,4.1,1.3;
5.5,2.5,4.0,1.3;
5.5,2.6,4.4,1.2;
6.1,3.0,4.6,1.4;
5.8,2.6,4.0,1.2;
5.0,2.3,3.3,1.0;
5.6,2.7,4.2,1.3;
5.7,3.0,4.2,1.2;
5.7,2.9,4.2,1.3;
6.2,2.9,4.3,1.3;
5.1,2.5,3.0,1.1;
5.7,2.8,4.1,1.3;
6.3,3.3,6.0,2.5;
5.8,2.7,5.1,1.9;
7.1,3.0,5.9,2.1;
6.3,2.9,5.6,1.8;
6.5,3.0,5.8,2.2;
7.6,3.0,6.6,2.1;
4.9,2.5,4.5,1.7;
7.3,2.9,6.3,1.8;
6.7,2.5,5.8,1.8;
7.2,3.6,6.1,2.5;
6.5,3.2,5.1,2.0;
6.4,2.7,5.3,1.9;
6.8,3.0,5.5,2.1;
5.7,2.5,5.0,2.0;
5.8,2.8,5.1,2.4;
6.4,3.2,5.3,2.3;
6.5,3.0,5.5,1.8;
7.7,3.8,6.7,2.2;
7.7,2.6,6.9,2.3;
6.0,2.2,5.0,1.5;
6.9,3.2,5.7,2.3;
5.6,2.8,4.9,2.0;
7.7,2.8,6.7,2.0;
6.3,2.7,4.9,1.8;
6.7,3.3,5.7,2.1;
7.2,3.2,6.0,1.8;
6.2,2.8,4.8,1.8;
6.1,3.0,4.9,1.8;
6.4,2.8,5.6,2.1;
7.2,3.0,5.8,1.6;
7.4,2.8,6.1,1.9; 
7.9,3.8,6.4,2.0;
6.4,2.8,5.6,2.2;
6.3,2.8,5.1,1.5;
6.1,2.6,5.6,1.4;
7.7,3.0,6.1,2.3;
6.3,3.4,5.6,2.4;
6.4,3.1,5.5,1.8;
6.0,3.0,4.8,1.8;
6.9,3.1,5.4,2.1;
6.7,3.1,5.6,2.4;
6.9,3.1,5.1,2.3;
5.8,2.7,5.1,1.9;
6.8,3.2,5.9,2.3;
6.7,3.3,5.7,2.5;
6.7,3.0,5.2,2.3;
6.3,2.5,5.0,1.9;
6.5,3.0,5.2,2.0;
6.2,3.4,5.4,2.3;
5.9,3.0,5.1,1.8;];

%初始化Vmax、Vmin
  Vmin1(1,1:4)=s(1,1:4);
  Vmax1(1,1:4)=s(1,1:4);
  Vmin2(1,1:4)=s(51,1:4);
  Vmax2(1,1:4)=s(51,1:4);
  Vmin3(1,1:4)=s(101,1:4);
  Vmax3(1,1:4)=s(101,1:4);

for j=1:4
    for i=1:50
        if(Vmin1(1,j)>s(i,j))
            Vmin1(1,j)=s(i,j);
        end
        if(Vmax1(1,j)<s(i,j))
           Vmax1(1,j)=s(i,j);
        end
    end
end
for j=1:4
     for i=51:100
         if(Vmin2(1,j)>s(i,j))
            Vmin2(1,j)=s(i,j);
         end
         if(Vmax2(1,j)<s(i,j))
            Vmax2(1,j)=s(i,j);
         end
     end
 end
 for j=1:4
     for i=101:150
         if(Vmin3(1,j)>s(i,j))
            Vmin3(1,j)=s(i,j);
         end
         if(Vmax3(1,j)<s(i,j))
            Vmax3(1,j)=s(i,j);
         end
     end
 end
Vmin(1:4)=Vmin1(1:4);
Vmin(5:8)=Vmin2(1:4);
Vmin(9:12)=Vmin3(1:4);
Vmax(1:4)=Vmax1(1:4);
Vmax(5:8)=Vmax2(1:4);
Vmax(9:12)=Vmax3(1:4);
x=zeros(N,K*D);
v=zeros(N,K*D);

for i=1:N
  %  x(i,1:D)=s(i,:);
  %  v(i,1:D)=0.15*s(i,:);
  %  x(i,1+D:2*D)=s(51+i,:);
  %  v(i,1+D:2*D)=0.15*s(50+i,:);
  %  x(i,2*D+1:3*D)=s(101+i,:);
 %   v(i,2*D+1:3*D)=0.15*s(100+i,:);
end
pos_temper=zeros(1,12);
%pos_temper1(1,1:12)=[4.8987 3.4133 1.5003 0.2664 5.8763 2.8459 4.2884 1.4052 6.4556 2.9685 5.6210 2.0203]; % 99.3689——99.9866
%pos_temper(1,1:12)=[4.9237 3.3702 1.4871 0.2633 5.8992 2.8553 4.3065 1.4326 6.4453 2.9356 5.6084 2.0124];%s(129,:);
%pos_temper1(1,1:12)=[4.8779 3.4490 1.5112 0.2633 5.8573 2.8381 4.2734 1.3653 6.4641 2.9959 5.4974 2.0268];%99.2019_99.5274
%pos_temper(1,1:12)=[4.8607 3.4199 1.4692 0.2716 5.8416 2.7620 4.2609 1.3561 6.4712 3.0186 5.5288 2.0069];%99.8-99.6

%pos_temper(1,1:12)=[4.8607 3.3862 1.4817 0.2470 5.8416 2.8094 4.2302 1.3618 6.4712 2.9982 5.5288 2.0110];%99.1905_99.1345
pos_temper(1,1:12)=[4.9436 3.4024 1.4564 0.2498 5.8171 2.8550 4.3097 1.4095 6.6065 3.0805 5.5740 2.0594];%97.8979   50代30个时为97.4308——97.5200
%pos_temper(1,1:12)=[5.0141 3.4608 1.4841 0.2602 5.9245 2.8619 4.3515 1.4359 6.6890 3.1282 5.7039 2.0625];%97.6001
%pos_temper(1,1:12)=(vtemper1(1,:)+vtemper2(1,:))/2;
%97.5200
%pos_temper(1,1:12)=[5.0141 3.4076 1.4629 0.2491 5.9245 2.8667 4.3232 1.4196 6.6801 3.1039 5.6998 2.0602];
%97.4308
%pos_temper(1,1:12)=[5.0141 3.4108 1.4702 0.2449 5.9245 2.8768 4.3384 1.4267 6.7027 3.1036 5.7039 2.0699];
for i=1:N
    for j=1:K*D   %粒子长度为：（聚类数目*每类的维数）
        x(i,j)=rand*(Vmax(1,j)+Vmin(1,j))/2; %随机初始化位置
        v(i,j)=rand*q*(Vmax(1,j)+Vmin(1,j))/2; %随机初始化速度
    end
end
%Vmin(1:12)=Vmin(1:12)MaxDT;
%Vmax(1:12)=Vmax(1:12)/MaxDT;

%------初始化种群的个体(可以在这里限定位置和速度的范围)------------
%x=crtrp(N,FieldDR);

%x

%------先计算各个粒子的适应度，并初始化Pi和Pg----------------------
format long;
for i=1:N
    p(i)=fitness(x(i,:),s,D,K,Sample);
    y(i,:)=x(i,:);
end

pg=x(1,:);             %Pg为全局最优
for i=2:N
    if(fitness(x(i,:),s,D,K,Sample)>fitness(pg,s,D,K,Sample))
        pg=x(i,:);
    end
end

%------进入主要循环，按照公式依次迭代，直到满足精度要求------------
for t=1:MaxDT
    for i=1:N
         %更新速度,并限制在[Vmin,Vmax]内
        w=0.8-(1+t)*0.2/MaxDT;
       %w=0.9-(0.5/MaxDT)*t;
      % w=-0.8;
        v(i,:)=w*v(i,:)+c1*rand*(y(i,:)-x(i,:))+c2*rand*(pg-x(i,:));
        for j=1:D*K
           if(v(i,j)>Vmax(1,j)*q)
               v(i,j)=Vmax(1,j)*q;
            elseif(v(i,j)<Vmin(1,j)*q)
               v(i,j)=Vmin(1,j)*q;
           end
        end
        pos(i,:)=x(i,:)+v(i,:);   %更新位置
        for j=1:D*K
            if(pos(i,j)>Vmax(1,j))
                pos(i,:)=pos_temper*0.85;
            elseif(pos(i,j)<Vmin(1,j))
                pos(i,:)=pos_temper*0.85;
            end
        end
        x(i,:)=pos(i,:);      
        if(fitness(x(i,:),s,D,K,Sample)>p(i))
            p(i)=fitness(x(i,:),s,D,K,Sample);
            y(i,:)=x(i,:);
        end        
        if(fitness(x(i,:),s,D,K,Sample)>fitness(pg,s,D,K,Sample))
            pg=x(i,:);
            vbest=v(i,:);
        end
       % disp('个体的适应度：');
       % 1/fitness(x(i,:),s,D,K,Sample) 
    end
%    Pbest(t)=fitness(pg,s,D,K);
    %disp('个体最好解：');
    % disp('全局最好解：');
    %   1/fitness(pg,s,D,K,Sample)
    
end
%------最后给出计算结果
format short;%long;
disp('*************************************************************');
disp('pos函数的全局最优位置为：')
Solution=pg'
disp('pos函数的全局最优速度为：')
vbest
disp('最后得到的优化极值为：')
Result=fitness(pg,s,D,K,Sample)
disp('最后结果是：')
dis=(1/Result)-1
disp('*************************************************************');