基于花授粉的模糊聚类的图像分割（MATLAB）

基于花授粉的模糊聚类

（自学整理代码，如有错误之处，还望各位大神留言指导）

function [best,fmin,U,N_iter]=C_FPA(para, data, k, expo, min_impro)

% min_impro:隶属度最小变化量，迭代终止条件
% expo: 隶属度矩阵U的加权指数

n=para(1);           % Population size, typically 10 to 25
p=para(2);           % 转换概率

% Iteration parameters
N_iter=1000;            % Total number of iterations

% Dimension of the search variables
d=2;
Lb=min(data)*ones(1,d);
Ub=max(data)*ones(1,d);
% Initialize the population/solutions
for i=1:n
%   Sol(i,:)=Lb+(Ub-Lb).*rand(1,d);
  Fitness1 = [];
  U = rand(k,size(data,1));
  col_sum = sum(U);
  U = U./col_sum(ones(k,1),:);
  mf = U.^expo;
  center(:,i) = mf*data./((ones(size(data,2),1)*sum(mf'))');
  Fitness1 = distfcm(center,data);
  Fitness(i) = sum(sum((Fitness1.^2)*mf'));
end

% Find the current best
[fmin,I]=min(Fitness);
best=center(:,I);
S=center; 

% Start the iterations -- Flower Algorithm 
for t=1:N_iter
        % Loop over all bats/solutions
        for i=1:n
          if rand>p %局部授粉
          L=Levy(d);
          dS=L'.*(center(:,i)-best);
          S(:,i)=center(:,i)+dS;
          S(:,i)=simplebounds(S(:,i),Lb,Ub);
          
          else  %全局授粉
              epsilon=rand;
              JK=randperm(n);
              S(:,i)=S(:,i)+epsilon*(center(:,JK(1))-center(:,JK(2)));
              S(:,i)=simplebounds(S(:,i),Lb,Ub);
          end
          
          % 更新新一代的花朵与适应度（函数变量和函数值）
            Fitness1 = [];
            Fitness1 = distfcm(S(:,i),data);
            Fnew(i) = sum(sum((Fitness1.^2)*mf'));
         
          if (Fnew(i)<=Fitness(i))
               center(:,i)=S(:,i);
               Fitness(i)=Fnew(i);
          end
           
          if Fnew(i)<=fmin
                best=S(:,i)   ;
                fmin=Fnew(i)   ;
          end
        end

        if round(t/10)==t/10
        best
        fmin
        end
        %终止条件判别
       if i>1
           if abs(Fnew(i)-Fnew(i-1))<min_impro
               break;
           end 
       end 
       tmp = Fitness1.^(-2/(expo-1));
       U = tmp./sum(tmp);   %隶属度矩阵
end


disp(['Total number of evaluations: ',num2str(N_iter*n)]);
disp(['fmin=', num2str(fmin)]);

% Application of simple constraints
function s=simplebounds(s,Lb,Ub)
  % Apply the lower bound
  ns_tmp=s;
  I=ns_tmp<Lb;
  ns_tmp(I)=Lb(I);
  
  % Apply the upper bounds 
  J=ns_tmp>Ub;
  ns_tmp(J)=Ub(J);
  % Update this new move 
  s=ns_tmp;


% 飞行参数
function L=Levy(d)
beta=3/2;
sigma=(gamma(1+beta)*sin(pi*beta/2)/(gamma((1+beta)/2)*beta*2^((beta-1)/2)))^(1/beta);
    u=randn(1,d)*sigma;
    v=randn(1,d);
    step=u./abs(v).^(1/beta);
L=0.01*step; 

%%子函数 计算距离
function out=distfcm(center,data)
out=zeros(size(center,1),size(data,1));
   for k=1:size(center,1)
       out(k,:)=sqrt(sum(((data-ones(size(data,1),1)*center(k,:)).^2)',1));
   end 

调用

close all
clc;
clear;
I=imread('3096.jpg');
%  figure; imshow(I,[]);title('yuan图像'); 
 if size(I,3) == 3
  I=rgb2gray(I);
   else
end
% figure;imshow(I);title('图像'); %imwrite(I,'1.tiff','tiff','Resolution',300);%保存为tif
g = double(I);
[p q] = size(g);
data = g(:);
data_n = size(data,1);
in_n = size(data,2);
k = 2;%类别数
expo = 2;
min_impro = 1e-5;
para = [25 0.8];  

[best,fmin,U,N_iter] = C_FPA(para, data, k, expo, min_impro);

maxU = max(U);
for j=1:k
    index = find(U(j, :) == maxU);
    fcm_label(index) = j;     
end
fcm_result=reshape(fcm_label,[p q]);
figure, imshow(fcm_result,[]);
title('fpa_fcm分割结果');