C_Meng PSNA

Never wait for the storm to pass, just dance in the rain.

0%

启发式算法学习（一）：粒子群算法PSO

Posted on Edited on In Valine:

算法简介

粒子群算法（Particle swarm optimization，PSO）是模拟群体智能所建立起来的一种优化算法，主要用于解决最优化问题（optimization problems）。1995年由 Eberhart和Kennedy 提出，是基于对鸟群觅食行为的研究和模拟而来的。

假设一群鸟在觅食，在觅食范围内，只在一个地方有食物，所有鸟儿都看不到食物（即不知道食物的具体位置。当然不知道了，知道了就不用觅食了），但是能闻到食物的味道（即能知道食物距离自己是远是近。鸟的嗅觉是很灵敏的）。

假设鸟与鸟之间能共享信息（即互相知道每个鸟离食物多远。这个是人工假定，实际上鸟们肯定不会也不愿意），那么最好的策略就是结合自己离食物最近的位置和鸟群中其他鸟距离食物最近的位置这2个因素综合考虑找到最好的搜索位置。

粒子群算法与《遗传算法》等进化算法有很多相似之处。也需要初始化种群，计算适应度值，通过进化进行迭代等。但是与遗传算法不同，它没有交叉，变异等进化操作。与遗传算法比较，PSO的优势在于很容易编码，需要调整的参数也很少。

核心概念

PSO有几个核心概念：

  1. 粒子（particle）：一只鸟。类似于遗传算法中的个体。
  2. 种群（population）：一群鸟。类似于遗传算法中的种群。
  3. 位置（position）：一个粒子（鸟）当前所在的位置。
  4. 经验（best）：一个粒子（鸟）自身曾经离食物最近的位置。
  5. 速度（velocity ）：一个粒子（鸟）飞行的速度。
  6. 适应度（fitness）：一个粒子（鸟）距离食物的远近。与遗传算法中的适应度类似。

算法过程

算法说明

两个核心公式

加速度更新公式：

$$v[i] = w * v[i] + c1 * rand() *(pbest[i] - present[i]) + c2 * rand() * (gbest - present[i])$$

其中v[i]代表第i个粒子的速度，w代表惯性权值,c1和c2表示学习参数，rand()表示在0-1之间的随机数,pbest[i]代表第i个粒子搜索到的最优值,gbest代表整个集群搜索到的最优值,present[i]代表第i个粒子的当前位置。

位置更新公式：

$$present[i]=present[i]+v[i]$$

解释说明

1.粒子数：粒子数的选取一般在20个到40个之间，但是需要具体问题具体对待，如果对于复杂问题，则需要设置更多的粒子，粒子数量越多，其搜索范围就越大。

2.惯性因子 $w$ ：用来控制继承多少粒子当前的速度的，越大则对于当前速度的继承程度越小，越小则对于当前速度的继承程度越大。有些同学可能会产生疑问，是不是说反了。其实不是，从公式中可以明确看出，其值越大，则速度的改变幅度就越大，则对于粒子的当前速度继承越小；反之，速度的改变幅度越小，则对于粒子当前速度继承越大。因此如果的值越大，则解的搜索范围越大，可以提高算法的全局搜索能力，但也损失了局部搜索能力，有可能错失最优解；反之如果的值越小，则解的搜索范围也就越小，算法的全局搜索能力也就越小，容易陷入局部最优。如果是变量，则其值应该随着迭代次数的增加而减小（类似于梯度下降当中的学习率）。如果为定值，则建议在0.6到0.75之间进行选取。

3.加速常数 $c1,c2$ ：通过公式一可以看出，加速常数控制着飞翔速度的计算是更加看重自身经验还是群体经验。公式一中的第二项就是自身经验的体现，加速常数可以看做是用来调整自身经验在计算粒子飞翔速度上的权重。同理是用来控制群体经验在计算粒子飞翔速度过程中的权重的。如果为0，则自身经验对于速度的计算不起作用，如果为0，则群体经验对于粒子飞翔速度的计算不起作用。的取值在学术界分歧很大主要有如下几种情况：

学者 参数取值
Clerc c1=c2=2.05
Carlisle c1=2.8, c2=1.3
Trelea w=0.6, c1=c2=1.7
Eberhart w=0.729, c1=c2=1.494

python实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
# -*- coding: utf-8 -*-
"""
f(x1,x2) = x1**2 + x2**2, x1,x2 belongs to [-10,10],求Min f
"""

import matplotlib.pyplot as plt
import numpy as np
 
 
class PSO(object):
    def __init__(self, population_size, max_steps):
        self.w = 0.6  # 惯性权重
        self.c1 = self.c2 = 2
        self.population_size = population_size  # 粒子群数量
        self.dim = 2  # 搜索空间的维度
        self.max_steps = max_steps  # 迭代次数
        self.x_bound = [-10, 10]  # 解空间范围
        self.x = np.random.uniform(self.x_bound[0], self.x_bound[1],
                                   (self.population_size, self.dim))  # 初始化粒子群位置
        self.v = np.random.rand(self.population_size, self.dim)  # 初始化粒子群速度
        fitness = self.calculate_fitness(self.x)
        self.p = self.x  # 个体的最佳位置
        self.pg = self.x[np.argmin(fitness)]  # 全局最佳位置
        self.individual_best_fitness = fitness  # 个体的最优适应度
        self.global_best_fitness = np.max(fitness)  # 全局最佳适应度
 
    def calculate_fitness(self, x):
        return np.sum(np.square(x), axis=1)
 
    def evolve(self):
        fig = plt.figure()
        for step in range(self.max_steps):
            r1 = np.random.rand(self.population_size, self.dim)
            r2 = np.random.rand(self.population_size, self.dim)
            # 更新速度和权重
            self.v = self.w*self.v+self.c1*r1*(self.p-self.x)+self.c2*r2*(self.pg-self.x)
            self.x = self.v + self.x
            plt.clf()
            plt.scatter(self.x[:, 0], self.x[:, 1], s=30, color='k')
            plt.xlim(self.x_bound[0], self.x_bound[1])
            plt.ylim(self.x_bound[0], self.x_bound[1])
            plt.pause(0.01)
            # plt.ion()
            # plt.show()
           
            fitness = self.calculate_fitness(self.x)
            # 需要更新的个体
            update_id = np.greater(self.individual_best_fitness, fitness)
            self.p[update_id] = self.x[update_id]
            self.individual_best_fitness[update_id] = fitness[update_id]
            # 新一代出现了更小的fitness，所以更新全局最优fitness和位置
            if np.min(fitness) < self.global_best_fitness:
                self.pg = self.x[np.argmin(fitness)]
                self.global_best_fitness = np.min(fitness)
            print('best fitness: %.5f, mean fitness: %.5f' % (self.global_best_fitness, np.mean(fitness)))
            
pso = PSO(10, 100)
pso.evolve()

Reference:
https://blog.csdn.net/zhaozx19950803/article/details/79854466
https://blog.csdn.net/yy2050645/article/details/80740641
https://blog.csdn.net/zj15527620802/article/details/81366105

点击这个按钮，会有好事发生哦~(>w<)~