基于粒子群优化算法的阈值分割方法

一、算法原理与数学模型

1. PSO阈值分割框架

\(\text{分割目标}：\min_{\tau} f(\tau) = \text{类间方差} \ \text{或} \ \text{信息熵}\)

• 粒子表示：每个粒子位置\(x_i\)对应候选阈值\(\tau\)
• 速度更新：
  \(v_{i}^{(k+1)} = w \cdot v_{i}^{(k)} + c_1 r_1 (pbest_i - x_i^{(k)}) + c_2 r_2 (gbest - x_i^{(k)})\)
• 位置更新：
  \(x_{i}^{(k+1)} = x_i^{(k)} + v_{i}^{(k+1)}\)
  （\(w\)：惯性权重，\(c_1,c_2\)：加速因子，\(r_1,r_2\)：[0,1]随机数）

2. 适应度函数设计

(1) Otsu准则（类间方差最大化）
\(f(\tau) = \omega_0(\tau)\omega_1(\tau)(\mu_0(\tau)-\mu_1(\tau))^2\)
其中\(\omega_0,\omega_1\)为两类的概率，\(\mu_0,\mu_1\)为均值

(2) 熵准则（信息最大化）
\(f(\tau) = H_L + H_R = -\sum_{i=0}^\tau p_i \log p_i - \sum_{i=\tau+1}^{L-1} p_i \log p_i\)
（\(H_L,H_R\)为左右区域熵，\(p_i\)为灰度级概率）

---

二、算法实现流程

1. 初始化阶段

% 参数设置
n_particles = 30;    % 粒子数量
max_iter = 100;      % 最大迭代次数
w = 0.729;           % 惯性权重
c1 = 1.494; c2 = 1.494;% 加速因子

% 粒子初始化（阈值范围[0,255]）
particles = randi([0,255],n_particles,1);
velocities = 0.1*randn(n_particles,1);
pbest = particles;   % 个体最优
gbest = particles(1);// 全局最优

2. 迭代优化过程

for iter = 1:max_iter
    % 计算适应度
    for i = 1:n_particles
        tau = particles(i);
        fitness(i) = otsu_fitness(I,tau); % 调用适应度函数
    end
    
    % 更新个体最优
    update_idx = fitness > prev_fitness;
    pbest(update_idx) = particles(update_idx);
    
    % 更新全局最优
    [best_fitness, best_idx] = max(fitness);
    if best_fitness > gbest_fitness
        gbest = particles(best_idx);
        gbest_fitness = best_fitness;
    end
    
    % 速度与位置更新
    velocities = w*velocities + ...
        c1*rand(n_particles,1).*(pbest - particles) + ...
        c2*rand(n_particles,1).*(gbest - particles);
    particles = particles + velocities;
    
    % 边界处理
    particles = max(0,min(255,particles));
end

3. 关键函数实现

function fitness = otsu_fitness(I,tau)
    % 计算类间方差
    [counts,~] = imhist(I);
    total = sum(counts);
    w0 = sum(counts(1:tau+1))/total;
    w1 = sum(counts(tau+2:end))/total;
    
    mu0 = sum((0:tau)' .* counts(1:tau+1)) / (w0*total);
    mu1 = sum((tau+1:255)' .* counts(tau+2:end)) / (w1*total);
    
    fitness = w0*w1*(mu0 - mu1)^2;
end

---

三、算法优化策略

1. 自适应参数调整

• 动态惯性权重：
  \(w = w_{max} - \frac{iter}{max\_iter}(w_{max}-w_{min})\)
  （初始\(w_{max}=0.9\), \(w_{min}=0.4\)）
• 加速因子调整：
  \(c_1 = c_{1\_init} \cdot e^{-\lambda iter}\)
  （\(\lambda=0.05\)控制收敛速度）

2. 混合优化策略

• K-Means+PSO初始化：

  % K-Means预分割
  [idx,centers] = kmeans(im2double(I),3);
  % PSO初始化粒子
  particles = centers(:,1)*255; 
  

• 遗传算法交叉：引入SBX交叉算子增强多样性

  function child = sbx_crossover(parent1,parent2,pc)
      if rand < pc
          eta = 2; % 分布指数
          u = rand(size(parent1));
          beta = (2*u).^(1/(eta+1));
          child = 0.5*((1+beta).*parent1 + (1-beta).*parent2);
      else
          child = parent1;
      end
  end
  

3. 并行计算加速

% 使用parfor并行计算适应度
parfor i = 1:n_particles
    fitness(i) = otsu_fitness(I,particles(i));
end

---

六、Matlab代码

1. 完整代码框架

function [threshold,seg] = pso_threshold_segmentation(I)
    % 参数设置
    n_particles = 30;
    max_iter = 100;
    w = 0.7; c1=1.5; c2=1.5;
    
    % 初始化粒子群
    particles = 255*rand(n_particles,1);
    velocities = 0.1*randn(n_particles,1);
    pbest = particles;
    [gbest,gbest_fitness] = update_global_best(I,particles);
    
    % 迭代优化
    for iter = 1:max_iter
        % 更新速度位置
        velocities = w*velocities + ...
            c1*rand(n_particles,1).*(pbest - particles) + ...
            c2*rand(n_particles,1).*(gbest - particles);
        particles = particles + velocities;
        particles = max(0,min(255,particles));
        
        % 更新最优解
        [current_best,fitness] = update_global_best(I,particles);
        if fitness > gbest_fitness
            gbest = current_best;
            gbest_fitness = fitness;
        end
        
        % 显示进度
        fprintf('Iter %d: Best Fitness=%.4f\n',iter,gbest_fitness);
    end
    
    % 输出结果
    threshold = gbest;
    seg = I > threshold;
end

function [best,fitness] = update_global_best(I,particles)
    n = numel(particles);
    fitness = zeros(n,1);
    for i=1:n
        fitness(i) = otsu_fitness(I,particles(i));
    end
    [best_fitness, idx] = max(fitness);
    best = particles(idx);
end

---

八、参考文献

刘桂红等. 改进PSO的Otsu图像分割方法. 计算机科学,2016
徐小慧等. 基于PSO的最佳熵阈值分割. 计算机工程与应用,2006
代码 基于粒子群算法的图像分割算法代码源程序 youwenfan.com/contentcna/52366.html

通过上述方法，可有效解决传统阈值分割方法依赖人工参数、对复杂图像适应性差等问题。实际应用中需根据图像特性调整PSO参数，并配合形态学后处理提升分割质量。