粒子群算法--背包问题

粒子群算法，即PSO，关于PSO的入门介绍请参看博文：http://www.cnblogs.com/maybe2030/p/5043356.html

PSO的核心是如下两个公式：

V表示粒子的速度，X表示粒子的位置，Pi为当前粒子找到的最优位置，Pj为整个粒子群找到的最优位置。c1、c2为正常数，称为加速因子，一般取0-1之间的某个值。w称为惯性因子，其取值直接关系到PSO的搜索能力，若w取值较大，则PSO搜索范围增加，有比较强的全局搜索能力，若w取值较小，则PSO的搜索范围变小，局部搜索增强。w的取值可以固定在0-1之间的某个值，也可以根据要解决的问题模型自行设计，进行线性变换或正态分布变换等。rand()为[0, 1]之间的随机数。PSO算法就是根据这两个公式，不断迭代、更新整个粒子群的速度、位置信息，直到迭代完成为止。

以下是PSO算法的总流程：

这里还是以01背包问题为例，介绍PSO算法在求解背包问题上的应用。关于背包问题的描述，请参看我的另一篇博文遗传算法--背包问题。

代码如下：（scala语言实现）

package main.scalaimport scala.collection.mutable.ArrayBufferimport scala.io.Sourceimport scala.util.Random// 种群中个体定义class Individual(indivi_bit:Int){  var individual_bit = indivi_bit // 个体中的位数，即背包中物品的数量  var knapsack_data:ArrayBuffer[(Int, Int)] = null // 背包数据的引用  var p_vec = new Array[Int](individual_bit) // 当前位置向量  var v_vec = new Array[Double](individual_bit) // 当前速度向量  var weight_all = 0 // 当前总重量  var fitness = 0 // 当前适应度值  var p_best:Array[Int] = new Array[Int](individual_bit) // 当前个体搜索到的最优位置信息  var best_fitness = 0 // 当前个体搜索到的最优位置的适应度值  var best_fitness_weight = 0 // 当前个体搜索到的最优位置的适应度值对应的总重量  var weight_limit = 1000 // 限重  def this(data:ArrayBuffer[(Int, Int)])  {    this(data.length)    this.knapsack_data = data    this.init  }  // 计算适应度  def calcFitness() = {    weight_all = 0    fitness = 0    for(i <- 0 until individual_bit)    {      p_vec(i) match {        case 0 =>        case _ => {          weight_all = weight_all + knapsack_data(i)._1          fitness = fitness + knapsack_data(i)._2        }      }    }    if(weight_all > weight_limit)      fitness = 0  }  // 初始化  def init = {    for(i <- 0 until individual_bit)    {      p_vec(i) = new Random().nextInt(2)      v_vec(i) = new Random().nextDouble()    }    p_vec.copyToArray(p_best)    calcFitness()    best_fitness = fitness    best_fitness_weight = weight_all  }  def updateCurBestInfo() = {    if(fitness > best_fitness)    {        best_fitness = fitness        best_fitness_weight = weight_all        p_vec.copyToArray(p_best)    }  }  // 复制个体  def copyTo(to:Individual) = {    to.individual_bit = this.individual_bit    to.knapsack_data = this.knapsack_data    to.weight_all = this.weight_all    to.fitness = this.fitness    to.best_fitness = this.best_fitness    to.best_fitness_weight = this.best_fitness_weight    to.weight_limit = this.weight_limit    if(to.p_vec == null)  to.p_vec = this.p_vec.clone()    else  this.p_vec.copyToArray(to.p_vec)    if(to.v_vec == null)  to.v_vec = this.v_vec.clone()    else  this.v_vec.copyToArray(to.v_vec)    if(to.p_best == null)  to.p_best = this.p_best.clone()    else  this.p_best.copyToArray(to.p_best)  }}object PSODemo{  val KNAPSACK_DATA_FILE = "knapsack.data"  val PSO_NUM = 20  // 种群数  val ITERATE_NUM = 30 // 迭代计算次数  val w = 0.8       // 惯性系数  val c1 = 0.7      // 加速因子1  val c2 = 0.7      // 加速因子2  var gbest:Individual = null // 最优个体  // 读取背包问题的初始数据  def readKnapsackData(path:String) = {    val list = new ArrayBuffer[(Int, Int)]()    for(line <- Source.fromFile(path).getLines())    {      val strs = line.split(" ")      list += ((strs(0).toInt, strs(1).toInt))    }    list  }  // 初始化粒子群  def InitPsoPopulation(data:ArrayBuffer[(Int, Int)]) = {    val population = new Array[Individual](PSO_NUM)    for(i <- 0 until PSO_NUM)    {      population(i) = new Individual(data)    }    population  }  // 更新粒子群中的最优个体  def updateBestIndividual(pso:Array[Individual]) = {    var population_best_fitness = 0    var individual_id = -1    for(i <- 0 until PSO_NUM)    {      if(population_best_fitness < pso(i).best_fitness)      {        population_best_fitness = pso(i).best_fitness        individual_id = i      }    }    if(individual_id >= 0)      pso(individual_id).copyTo(gbest)  }  // 更新粒子群中的每个个体  def updateEachIndividual(pso:Array[Individual]) = {    for(i <- 0 until PSO_NUM)    {      val v_new = pso(i).v_vec.clone() // 新的速度信息      val p_new = pso(i).p_vec.clone() // 新的位置信息      val r1 = new Random().nextDouble()      val r2 = new Random().nextDouble()      for(j <- 0 until pso(i).individual_bit)      {        v_new(j) = w * pso(i).v_vec(j) + c1 * r1 * (pso(i).p_best(j) - pso(i).p_vec(j)) +          c2 * r2 * (gbest.p_best(j) - pso(i).p_vec(j))        val p_tmp = pso(i).p_vec(j) + v_new(j)        if(p_tmp < 0.5)          p_new(j) = 0        else          p_new(j) = 1      }      v_new.copyToArray(pso(i).v_vec)      p_new.copyToArray(pso(i).p_vec)      // 计算适应度      pso(i).calcFitness()      // 更新该个体当前最优位置信息      pso(i).updateCurBestInfo()    }  }  def main(args: Array[String]):Unit = {    // 每个元素为 (weight, value)    val knapsack_data = readKnapsackData(KNAPSACK_DATA_FILE)    gbest = new Individual(knapsack_data.length)    // 初始化粒子群，每个个体为 (位置向量, 速度向量)    val pso = InitPsoPopulation(knapsack_data)    // 开始迭代计算    for(i <- 0 to ITERATE_NUM)    {      // 更新最优个体      updateBestIndividual(pso)      println("第" + i + "代，最优个体：" + gbest.p_best.mkString(",") +        "，总重量：" + gbest.best_fitness_weight + "，总价值：" + gbest.best_fitness)      // 更新每个个体的位置信息、速度信息      updateEachIndividual(pso)    }  }}

背包物品数据定义在"knapsack.data"文件中，如下：

每一行数据代表一个物品，其中，第1个值表示物品的重量weight，第2个值表示物品的价值value，背包最大承重量为1000。

源码就是整个PSO流程的最好说明，源码中已给出了比较详尽的注释。在PSODemo.updateEachIndividual()接口中，更新种群中每个粒子的速度、位置信息，其计算方式即按照上述的两个公式来计算。

这里有一个问题，PSO在解决背包问题时，用粒子的位置向量表示背包问题的一个解，即每个粒子的位置向量是一个01数组，其中，数组的每个元素取值只有取0或取1两种情况，但根据公式更新时，位置+速度后的结果可能是0.x或1.x，甚至有可能是-0.x，出现这种情况时，需要把连续的取值给离散化为取0、或取1。代码中我采用了一种比较粗暴的方式，若计算结果<0.5，取0，否则取1。这种离散化方式其实效果很差，求解时收敛很快，很容易陷入局部最优。当然，也可以采用其他的离散化方式来改善PSO算法的性能，比如进行异或、同或等，这就要具体问题具体分析了。从这也可以看出，PSO算法比较适合求解连续性优化问题，如函数优化等，这种离散化问题用PSO求解不是很合适。