quick_pso_svm.m

利用PSO优化SVM

clc
clear
close all
tic
%---------------------------------------------------
load cancer_input;
load cancer_output;

%num = 0;  %%控制迭代次数
num_svm = 0;  %%计算支持向量总数

%%---------------------------------------------------
%%将训练样本分成2的整数倍份。
%%前一份作为训练样本，后一分作为训练样本，分别进行PSO优化
%%找出支持向量

%for n = 1:13
    %while num < 12 & Percent < 0.9
        %num = num+1;
max_iterations  = 20;
no_of_particles = 20;
dimensions      = 2;

%delta_min       = -0.003;
%delta_max       = 0.003;

c1 = 1.3;
c2 = 1.3;
%---------------------------------------------------
% 产生总体测试样本
TestX = cancer_input(301:367,:)
for i = 1:67
    TestY(i) = cancer_output(300+i) - 3;
end

%---------------------------------------------------
%initialise the particles and teir velocity components

for count_x = 1:no_of_particles
    for count_y = 1:dimensions
        particle_position(count_x,count_y) = rand*10;
        particle_velocity(count_x,count_y) = rand*1000;
        p_best(count_x,count_y) = particle_position(count_x,count_y);
    end
end

%initialize the p_best_fitness array
for count = 1:no_of_particles
    p_best_fitness(count) = -1000;
end

for n = 1:30
%particle_position
%particle_velocity
count = 1;
Percent = 0;
%main particle swrm routine
while count < max_iterations & Percent < 0.99
%for count = 1:max_iterations
    
    %find the fitness of each particle
    %change fitness function as per equation requiresd and dimensions
    for count_x = 1:no_of_particles
        x = particle_position(count_x,1);
        y = particle_position(count_x,2);
%---------------------------------------------------
% 产生训练样本与测试样本
% 特别注意：此工具箱用于分类时，只能处理2类分类，且目标值必须为 1 或 -1。

        %n1 = [rand(3,5),rand(3,5)+1];
        %x1 = [1*ones(1,5),-1*ones(1,5)];
%%产生接龙式样本
    if n <= 30    
        for i = 1:10
            for j = 1:9
                n1(j,i) = cancer_input(n*20+i,j);        
            end
        end
        for i = 1:10
            x1(i) = cancer_output(n*20+i)-3;
        end
    
        %n2 = [rand(3,5),rand(3,5)+1];
        %x2 = [1*ones(1,5),-1*ones(1,5)];
        for i = 1:67
            for j = 1:9
                %n2(j,i) = cancer_input((n+1)*20+i,j);   
                n2(j,i) = cancer_input(300+i,j);
            end
        end
        for i = 1:67
            %x2(i) = cancer_output((n+1)*20+i)-3;
            x2(i) = cancer_output(300+i)-3;
        end
    else
        for i = 1:10
            for j = 1:9
                n1(j,i) = cancer_input(n*10+i,j);        
            end
        end
        for i = 1:10
            x1(i) = cancer_output(n*10+i)-3;
        end
    
        %n2 = [rand(3,5),rand(3,5)+1];
        %x2 = [1*ones(1,5),-1*ones(1,5)];
        for i = 1:10
            for j = 1:9
                %n2(j,i) = cancer_input(i,j); 
                n2(j,i) = cancer_input(300+i,j);
            end
        end
        for i = 1:10
            %x2(i) = cancer_output(i)-3;
            x2(i) = cancer_output(300+i)-3;
        end
    end
    
        xn_train = n1;          % 训练样本，每一列为一个样本
        dn_train = x1;          % 训练目标，行向量

        xn_test = n2;           % 测试样本，每一列为一个样本
        dn_test = x2;           % 测试目标，行向量

%---------------------------------------------------
% 参数设置

        trnX = xn_train';
        trnY = dn_train';
        tstX = xn_test';
        tstY = dn_test';

        ker = 'rbf';        % 核函数 k = exp(-(u-v)*(u-v)'/(2*p1^2))
        global p1 ;
        p1 = x;             % p1 is width of rbfs (sigma)
        C = y;             % 折衷系数

%---------------------------------------------------
% 训练与测试

        [nsv,alpha,bias] = svc(trnX,trnY,ker,C);                        % 训练

        actfunc = 0;                                                    % 1 为实际输出，0 为取sign输出 
        predictedY = svcoutput(trnX,trnY,tstX,ker,alpha,bias,actfunc);  % 测试

%---------------------------------------------------
% 结果统计

        Result = ~abs(predictedY-tstY)               % 正确分类显示为1
        Percent = sum(Result)/length(Result)   % 正确分类率
        soln = 1-Percent
        
        %x = particle_position(count_x,1);
        %y = particle_position(count_x,2);
        %z = particle_position(count_x,3);
        %soln = x^2 - 3*y*x + z;
        
        %x = particle_position(count_x); 
        %soln = x^2-2*x+1;
        
       % x = particle_position(count_x);
   % soln = x-7;
        
        if soln~=0 
            current_fitness(count_x) = 1/abs(soln)+0.0001;
        else
            current_fitness(count_x) =1000;
        end
    end
    
    %decide on p_best etc for each particle
    for count_x = 1:no_of_particles
        if current_fitness(count_x) > p_best_fitness(count_x)
            p_best_fitness(count_x) = current_fitness(count_x);
            for count_y = 1:dimensions
                p_best(count_x,count_y) = particle_position(count_x,count_y);
            end
        end
    end
        
    %decide on the global best among all the particles
    [g_best_val,g_best_index] = max(current_fitness);
    
    %g_best contains the position of teh global best
    for count_y = 1:dimensions
        g_best(count_y) = particle_position(g_best_index,count_y);      
    end

    %update the position and velocity compponents
    for count_x = 1:no_of_particles
        for count_y = 1:dimensions
            p_current(count_y) = particle_position(count_x,count_y);
        end

        for count_y = 1:dimensions
            particle_velocity(count_y) = particle_velocity(count_y) +  c1*rand*(p_best(count_y)-p_current(count_y)) + c2*rand*(g_best(count_y)-p_current(count_y));
            particle_positon(count_x,count_y) = p_current(count_y) +particle_velocity(count_y);
        end
    end  
    count = count + 1;
end
%%------------------------------------------------------
%%计算支持向量
        p1 = g_best(1);             % p1 is width of rbfs (sigma)
        C =  g_best(2);             % 折衷系数
        [nsv,alpha,bias] = svc(trnX,trnY,ker,C);                        % 训练
        actfunc = 0;                                                    % 1 为实际输出，0 为取sign输出 
        PredictedY(:,n) = svcoutput(trnX,trnY,TestX,ker,alpha,bias,actfunc);  % 测试
        Result = ~abs(PredictedY(:,n)-tstY)               % 正确分类显示为1
        PercentY(n) = sum(Result)/length(Result)   % 正确分类率
       
        for i = 1:20
         %   afa(n,i) = alpha(i);
          %  b0(n) = bias;
            if alpha(i)<0.0001
                num_svm = num_svm + 1;
             %   afa(n,i) = 0;
                for j = 1:9
                    sv_in(j,num_svm) = cancer_input(n*20+i,j);            
                end
                sv_out(num_svm) = cancer_output(n*20+i)-3;        
            end                 
        end
        for i = 2
            best(n,i) = g_best(i);
        end            
end
%------------------------
%加权输出

p_max = max(PercentY);
p_min = min(PercentY);
for i = 1:67
    predicted_Y(i) = 0;
end
for i = 1:30
    P(i) = (PercentY(i)-p_min)/(p_max-p_min);
    PredictedY(:,i) = PredictedY(:,i)*P(i) ;
end
for i = 1:67
    predicted_Y(i) = predicted_Y(i)+PredictedY(:,i);
end
predicted_Y = sign(predicted_Y);

Result = ~abs(predicted_Y-TestY)               % 正确分类显示为1
Percent = sum(Result)/length(Result)   % 正确分类率



t=toc;