读书笔记6

PSO粒子群算法

算法简介：

粒子群算法（Particle Swarm Optimization，简称PSO），或称粒子群优化，是属于人工智能算法，公元1995年由肯尼迪（Kennedy）与埃伯哈特（Eberhart）（1995）两位学者所提出，这两位学者借由观察鸟类族群觅食的讯息传递所得到的一个启发，粒子群算法的理论基础是以单一粒子来做为鸟类族群之中的单一个体，于算法中赋予该粒子（个体）拥有记忆性，并能够透过与粒子群体中的其他粒子之间的互动而寻求到最适解。因此在粒子群算法的基础理论可以理解，任一个体（粒子）皆可用有自身移动过程中所产生的记忆与经验，当个体移动的同时，能依造自身的经验与记忆来学习调整自身的移动方向，由于在粒子群算法中，多个粒子是同时移动的，且同时以自身经验与其他粒子所提供的经验进行比对找寻最适当的解，并使自己处于最适解中，该粒子群算法的特性使得粒子不单单受自身演化的影响，同时会对群体间的演化拥有学习性、记忆性，并使粒子本身达到最佳调整。

基本思想：

粒子群算法通过设计一种无质量的粒子来模拟鸟群中的鸟，粒子仅具有两个属性：速度和位置，速度代表移动的快慢，位置代表移动的方向。每个粒子在搜索空间中单独的搜寻最优解，并将其记为当前个体极值，并将个体极值与整个粒子群里的其他粒子共享，找到最优的那个个体极值作为整个粒子群的当前全局最优解，粒子群中的所有粒子根据自己找到的当前个体极值和整个粒子群共享的当前全局最优解来调整自己的速度和位置。

算法的实现：

假设目标函数有n个输入参数：f(x1, x2, ..., xn)。

对于每个粒子，它包含的元素：

(1) 当前时刻位置：(x1, x2, ..., xn)

(2) 当前时刻的目标函数值（也称为适应度）：f(x1, x2, ..., xn)

(3) 该粒子的历史最优位置：(x1_pbest, x2_pbest, ..., xn_pbest)

(4) 该粒子的历史最优目标函数值：f_pbest = f(x1_pbest, x2_pbest, ..., xn_pbest)

对于全部粒子，它们共有的元素：

(1) 粒子总数：num

(2) 迭代总次数：cnt。粒子群算法也是一个迭代的过程，需要多次迭代才能获取到理想的最优解。

(3) 全局最优位置：(x1_gbest, x2_gbest, ..., xn_gbest)

(4) 全局最优目标函数值：f_gbest = f(x1_gbest, x2_gbest, ..., xn_gbest)

(5) 位置随机化的上下限：xmin，xmax。迭代开始的时候以及迭代的过程中，均需要对粒子的位置进行随机分布，需要设置随机分布的上下限，不然随机分布偏离得太远，会严重影响优化结果。

(6) 速度的上下限：Vmin，Vmax。迭代过程中，速度也具有一定的随机性，需要限制速度的大小在一定范围内，不然如果速度值太大也会严重影响优化结果。

(7) 速度计算参数：c1、c2。通常取1.0~1.8的值。

粒子位置的更新公式：

参考上述提到的博客，对于每一个粒子的每一个位置参数xi，其位置更新公式如下，其中randf为0~1的随机数，c1与c2通常取1.0~1.8之间的固定值。

为了加快收敛速度，根据速度具有惯性的原理，后来人们提出了在速度计算中增加惯性，也即加上上一轮迭代时该位置参数的速度。如下式，其中w为0~1的参数，通常随着迭代次数的增加，逐渐减小w，因为越到后面可能解就越接近真实解，迭代收敛就越慢，所以需要减小w来放慢速度，否则容易错过最优解。 

 算法流程图如下：

#include <iostream>
#include <vector>
#include <cmath>
#include <map>
#include <algorithm>
#include <random>
#include <ctime>
#include <Eigen/Dense>
using namespace Eigen;
using namespace std;

const int dim = 1;//维数
const int p_num = 10;//粒子数量
const int iters = 100;//迭代次数
const int inf = 999999;//极大值
const double pi = 3.1415;
//定义粒子的位置和速度的范围
const double v_max = 4;
const double v_min = -2;
const double pos_max = 2;
const double pos_min = -1;
//定义位置向量和速度向量
vector<double> pos;
vector<double> spd;
//定义粒子的历史最优位置和全局最优位置
vector<double> p_best;
double g_best;
//使用eigen库定义函数值矩阵和位置矩阵
Matrix<double, iters, p_num> f_test;
Matrix<double, iters, p_num> pos_mat;

//定义适应度函数
double fun_test(double x)
{
    double res = x * x + 1;
    return res;
}

//初始化粒子群的位置和速度
void init()
{
    //矩阵中所有元素初始化为极大值
    f_test.fill(inf);
    pos_mat.fill(inf);
    //生成范围随机数
    static std::mt19937 rng;
    static std::uniform_real_distribution<double> distribution1(-1, 2);
    static std::uniform_real_distribution<double> distribution2(-2, 4);
    for (int i = 0; i < p_num; ++i)
    {
        pos.push_back(distribution1(rng));
        spd.push_back(distribution2(rng));
    }
    vector<double> vec;
    for (int i = 0; i < p_num; ++i)
    {
        auto temp = fun_test(pos[i]);//计算函数值
        //初始化函数值矩阵和位置矩阵
        f_test(0, i) = temp;
        pos_mat(0, i) = pos[i];
        p_best.push_back(pos[i]);//初始化粒子的历史最优位置
    }
    std::ptrdiff_t minRow, minCol;
    f_test.row(0).minCoeff(&minRow, &minCol);//返回函数值矩阵第一行中极小值对应的位置
    g_best = pos_mat(minRow, minCol);//初始化全局最优位置
}

void PSO()
{
    static std::mt19937 rng;
    static std::uniform_real_distribution<double> distribution(0, 1);
    for (int step = 1; step < iters; ++step)
    {
        for (int i = 0; i < p_num; ++i)
        {
            //更新速度向量和位置向量
            spd[i] = 0.5 * spd[i] + 2 * distribution(rng) * (p_best[i] - pos[i]) +
                2 * distribution(rng) * (g_best - pos[i]);
            pos[i] = pos[i] + spd[i];
            //如果越界则取边界值
            if (spd[i] < -2 || spd[i] > 4)
                spd[i] = 4;
            if (pos[i] < -1 || pos[i] > 2)
                pos[i] = -1;
            //更新位置矩阵
            pos_mat(step, i) = pos[i];
        }
        //更新函数值矩阵
        for (int i = 0; i < p_num; ++i)
        {
            auto temp = fun_test(pos[i]);
            f_test(step, i) = temp;
        }
        for (int i = 0; i < p_num; ++i)
        {
            MatrixXd temp_test;
            temp_test = f_test.col(i);//取函数值矩阵的每一列
            std::ptrdiff_t minRow, minCol;
            temp_test.minCoeff(&minRow, &minCol);//获取每一列的极小值对应的位置
            p_best[i] = pos_mat(minRow, i);//获取每一列的极小值，即每个粒子的历史最优位置
        }
        g_best = *min_element(p_best.begin(), p_best.end());//获取全局最优位置
    }
    cout << fun_test(g_best);
}

int main()
{
    init();
    PSO();
    system("pause");
    return 0;
}

PSO算法的优缺点：

优点：

1、不依赖于问题信息，采用实数求解，算法通用性强。
2、需要调整的参数少，原理简单，容易实现，这是PSO算法的最大优点。
3、协同搜索，同时利用个体局部信息和群体全局信息指导搜索。
4、收敛速度快， 算法对计算机内存和CPU要求不高。
5、更容易飞越局部最优信息。对于目标函数仅能提供极少搜索最优值的信息，在其他算法无法辨别搜索方向的情况下，PSO算法的粒子具有飞越性的特点使其能够跨过搜索平面上信息严重不足的障碍，飞抵全局最优目标值。比如Generalized Rosenbrock函数全局最小值在原占附近.但是此函数全局最优值与可到达的局部最优值之间右一条独长的山公，曲面山谷中点的最速下降方向几乎与到函数最小值的最佳方向垂直，找到全局最小值的可能性微乎其微， 但是PSO算法完全有可能找到全局最优值。

缺点：

1、算法局部搜索能力较差，搜索精度不够高。

2、算法不能绝对保证搜索到全局最优解，主要有两方面的原因：

①有时粒子群在俯冲过程中会错失全局最优解。粒子飞翔过程中的俯冲动作使搜索行为不够精细，不容易发现全局最优目标值，所以对粒子的最大飞翔速度进行限制既是为了使粒子不要冲出搜索区域的边界，同时也是为了使搜索行为不至于太粗糙。
②应用PSO算法处理高维复杂问题时，算法可能会早熟收敛，也就是粒子群在没有找到全局最优信息之前就陷入停顿状态，飞翔的动力不够，粒子群丧失了多样性，各粒子之间的抱合力增强，紧紧地聚集在一起，并且它们的飞翔速度几乎为零，虽然此时粒子距离全局最优解并不远，但是几乎为零的飞翔速度使其跳出停滞不前的状态，各个粒子力不从心。这些停滞不前的早熟点未必都是局部最优点，也可能是位于局部最优点邻域内的其他点，这一点与梯度搜索法不同，梯度搜索法如果出现早熟，通常只会陷人局部最优点，而不可能陷人局部最优点邻域内的其他点，这一点与梯度搜索算法不同，梯度搜索算法如果出现早熟，通常只会陷入局部最优点，而·不可能陷入局部最优点领域内的其它点。
3、算法搜索性能对参数具有一定的依赖性。对于特定的优化问题，如果用户经验不足，参数调整的确是个棘手的问题。参数值的大小直接影响到算法是否收敛以及求解结果的精度。
4、PSO算法是一种概率算法，算法理论不完善，缺乏独特性，理论成果偏少。从数学角度严格证明算法结果的正确性和可靠性还比较困难；缺少算法结构设计和参数选取的实用性指导原则，特别是全局收敛研究和大型多约束非线性规划的研究成果非常少。