快捷導(dǎo)航

Python代碼實(shí)現(xiàn)粒子群算法圖文詳解

更新時(shí)間：2021年08月26日 16:43:03 作者：大灰狼學(xué)編程

這篇文章主要介紹了Python編程實(shí)現(xiàn)粒子群算法詳解，涉及粒子群算法的原理，過(guò)程，以及實(shí)現(xiàn)代碼示例，具有一定參考價(jià)值，需要的朋友可以了解下

1.引言

粒子群優(yōu)化算法起源于對(duì)鳥群覓食活動(dòng)的分析。鳥群在覓食的時(shí)候通常會(huì)毫無(wú)征兆的聚攏，分散，以及改變飛行的軌跡，但是在不同個(gè)體之間會(huì)十分默契的保持距離。所以粒子群優(yōu)化算法模擬鳥類覓食的過(guò)程，將待求解問(wèn)題的搜索空間看作是鳥類飛行的空間，將每只鳥抽象成一個(gè)沒(méi)有質(zhì)量和大小的粒子，用這個(gè)粒子來(lái)表示待求解問(wèn)題的一個(gè)可行解。所以，尋找最優(yōu)解的過(guò)程就相當(dāng)于鳥類覓食的過(guò)程。

粒子群算法也是基于種群以及進(jìn)化的概念，通過(guò)個(gè)體間的競(jìng)爭(zhēng)與協(xié)作，實(shí)現(xiàn)復(fù)雜空間最優(yōu)解的求解。但是與遺傳算法不同的是，他不會(huì)對(duì)每個(gè)個(gè)體進(jìn)行“交叉”，“變異”等操作，而實(shí)以一定的規(guī)則，更新每個(gè)粒子的速度以及位置，使得每一個(gè)粒子向自身歷史最佳位置以及全局歷史最佳位置進(jìn)行移動(dòng)，從而實(shí)現(xiàn)整個(gè)種群向著最優(yōu)的方向進(jìn)化。

2.算法的具體描述：

2.1原理

在粒子群優(yōu)化算法中，粒子之間通過(guò)信息共享機(jī)制，獲得其它粒子的發(fā)現(xiàn)與飛行經(jīng)歷。粒子群算法中的信息共享機(jī)制實(shí)際上是一種合作共生的行為，在搜索最優(yōu)解的過(guò)程中，每個(gè)粒子能夠?qū)ψ约航?jīng)過(guò)的最佳的歷史位置進(jìn)行記憶，同時(shí)，每個(gè)粒子的行為有會(huì)受到群體中其他例子的影響，所以在搜索最優(yōu)解的過(guò)程中，粒子的行為既受其他粒子的影響，有受到自身經(jīng)驗(yàn)的指導(dǎo)。

粒子群優(yōu)化算法對(duì)于鳥群的模擬是按照如下的模式進(jìn)行的：假設(shè)一群鳥在空中搜索食物，所有鳥知道自己當(dāng)前距離食物有多遠(yuǎn)（這里的遠(yuǎn)近會(huì)用一個(gè)值來(lái)衡量，適應(yīng)度值），那么每只鳥最簡(jiǎn)單的搜索策略就是尋找距離目前距離食物最近的鳥的周圍空間。因此，在粒子群算法中，每個(gè)粒子都相當(dāng)于一只鳥，每個(gè)粒子有一個(gè)適應(yīng)度值，還有一個(gè)速度決定他們的飛行的距離與方向。所有的粒子追隨當(dāng)前最優(yōu)的粒子在解空間中搜索。每搜索一次，最優(yōu)的粒子會(huì)發(fā)生變化，其他的粒子又會(huì)追隨新的最優(yōu)粒子進(jìn)行搜索，如此反復(fù)迭代。

在迭代開始的時(shí)候，每個(gè)粒子通過(guò)隨機(jī)的方式初始化在空間中的速度和位置，然后在迭代過(guò)程中，粒子通過(guò)跟蹤兩個(gè)極值來(lái)自己在解空間中的位置和速度，一個(gè)極值是單個(gè)粒子自身在迭代的過(guò)程中的最優(yōu)位置（就是最優(yōu)適應(yīng)度值所對(duì)應(yīng)的空間解），這個(gè)稱之為粒子的個(gè)體極值。另一個(gè)極值是種群中所有的粒子在迭代過(guò)程中所找到的最優(yōu)位置，這個(gè)成為全局極值。如果粒子只是跟蹤一個(gè)極值的話，則算法稱為局部粒子群算法或者全局粒子群算法。

PSO從這種模型中得到啟示并用于解決優(yōu)化問(wèn)題。PSO 中，每個(gè)優(yōu)化問(wèn)題的潛在解都是搜索空間中的一只鳥，稱之為粒子。所有的粒子都有一個(gè)由被優(yōu)化的函數(shù)決定的適值( fitness value) ，每個(gè)粒子還有一個(gè)速度決定它們飛翔的方向和距離。然后粒子們就追隨當(dāng)前的最優(yōu)粒子在解空間中搜索。

PSO初始化為一群隨機(jī)粒子(隨機(jī)解)，然后通過(guò)迭代找到最優(yōu)解。在每一次迭代中，粒子通過(guò)跟蹤兩個(gè)極值來(lái)更新自己；第一個(gè)就是粒子本身所找到的最優(yōu)解，這個(gè)解稱為個(gè)體極值；另一個(gè)極值是整個(gè)種群目前找到的最優(yōu)解，這個(gè)極值是全局極值。另外也可以不用整個(gè)種群而只是用其中一部分作為粒子的鄰居，那么在所有鄰居中的極值就是局部極值。

圖解：

2.2標(biāo)準(zhǔn)粒子群算法流程

 算法的流程如下：

Step1:種群初始化:可以進(jìn)行隨機(jī)初始化或根據(jù)被優(yōu)化的問(wèn)題設(shè)計(jì)特定的初始化方法,包括群體規(guī)模，每個(gè)粒子的位置 X i X_{i} Xi 和速度 V i V_i Vi ，然后計(jì)算每個(gè)粒子的適應(yīng)度值，從而選擇出個(gè)體的局部最優(yōu)位置向量和種群的全局最優(yōu)位置向量。

Step2:迭代設(shè)置:設(shè)置迭代次數(shù) g m a x g_{max} gmax ，令當(dāng)前迭代次數(shù)g=1。

Step3:根據(jù)公式更新每個(gè)粒子的速度向量V。

Step4:根據(jù)公式更新每個(gè)粒子的位置向量X。

Step5:局部位置向量和全局位置向量更新:更新每個(gè)粒子的Pbest,和種群的Gbest。

Step6:終止條件判斷:判斷迭代次數(shù)時(shí)都達(dá)到 g m a x g_{max} gmax 或誤差是否足夠小，如果滿足則輸出Gbest.否則繼續(xù)進(jìn)行迭代，跳轉(zhuǎn)至步驟(3)。

對(duì)于粒子群優(yōu)化算法的實(shí)際應(yīng)用，因?yàn)橹饕菍?duì)速度和位置向量迭代算子的設(shè)計(jì)計(jì),選代算子是否合理將決定整個(gè)PSO算法性能的優(yōu)劣.，所以如何設(shè)計(jì) t pso的迭代算子是算法應(yīng)用的研究重點(diǎn)和難點(diǎn)。

3.代碼案例

3.1問(wèn)題

求解f(x,y)的最小值點(diǎn)

3.2繪圖

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D

# 生成X和Y的數(shù)據(jù)
X=np.arange(-5,5,0.1)
Y=np.arange(-5,5,0.1)
X,Y=np.meshgrid(X,Y)

# 目標(biāo)函數(shù)
Z=X**2+Y**2+X

# 繪圖
fig=plt.figure()
ax=Axes3D(fig)
surf=ax.plot_surface(X,Y,Z,cmap=cm.coolwarm)
plt.show()

3.3計(jì)算適應(yīng)度

# 速度
# Vi+1 = w*Vi + c1 * r1 * (pbest_i - Xi) + c2 * r2 * (gbest_i - Xi)
# 位置
# Xi+1 = Xi + Vi+1
# vi, xi 分別表示粒子第i維的速度和位置
# pbest_i, gbest_i 分別表示某個(gè)粒子最好位置第i維的值、整個(gè)種群最好位置第i維的值

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl

mpl.rcParams['font.sans-serif'] = ['SimHei']  # 指定默認(rèn)字體
mpl.rcParams['axes.unicode_minus'] = False  # 解決保存圖像是負(fù)號(hào)'-'顯示為方塊的問(wèn)題

def fitness_func(X):
    """計(jì)算粒子的的適應(yīng)度值，也就是目標(biāo)函數(shù)值，X 的維度是 size * 2 """
    A = 10
    pi = np.pi
    x = X[:, 0]
    y = X[:, 1]
    return x**2+y**2+x

3.4更新速度

def velocity_update(V, X, pbest, gbest, c1, c2, w, max_val):
    """
    根據(jù)速度更新公式更新每個(gè)粒子的速度
    :param V: 粒子當(dāng)前的速度矩陣，20*2 的矩陣
    :param X: 粒子當(dāng)前的位置矩陣，20*2 的矩陣
    :param pbest: 每個(gè)粒子歷史最優(yōu)位置，20*2 的矩陣
    :param gbest: 種群歷史最優(yōu)位置，1*2 的矩陣
    """
    size = X.shape[0]
    r1 = np.random.random((size, 1))
    r2 = np.random.random((size, 1))
    V = w*V+c1*r1*(pbest-X)+c2*r2*(gbest-X)
    # 防止越界處理
    V[V < -max_val] = -max_val
    V[V > -max_val] = max_val
    return V

3.5更新粒子位置

def position_update(X, V):
    """
    根據(jù)公式更新粒子的位置
    :param X: 粒子當(dāng)前的位置矩陣，維度是 20*2
    :param V: 粒子當(dāng)前的速度舉著，維度是 20*2
    """
    return X+V

3.6主要算法過(guò)程

def pos():
    w = 1
    c1 = 2
    c2 = 2
    r1 = None
    r2 = None
    dim = 2
    size = 20
    iter_num = 1000
    max_val = 0.5
    best_fitness = float(9e10)
    fitness_val_list = []
    # 初始化種群各個(gè)粒子的位置
    X = np.random.uniform(-5, 5, size=(size, dim))
    # 初始化各個(gè)粒子的速度
    V = np.random.uniform(-0.5, 0.5, size=(size, dim))
    # print(X)
    p_fitness = fitness_func(X)
    g_fitness = p_fitness.min()
    fitness_val_list.append(g_fitness)

    # 初始化的個(gè)體最優(yōu)位置和種群最優(yōu)位置
    pbest = X
    gbest = X[p_fitness.argmin()]
    # 迭代計(jì)算
    for i in range(1, iter_num):
        V = velocity_update(V, X, pbest, gbest, c1, c2, w, max_val)
        X = position_update(X, V)
        p_fitness2 = fitness_func(X)
        g_fitness2 = p_fitness2.min()

        # 更新每個(gè)粒子的歷史最優(yōu)位置
        for j in range(size):
            if p_fitness[j] > p_fitness2[j]:
                pbest[j] = X[j]
                p_fitness[j] = p_fitness2[j]
            # 更新群體的最優(yōu)位置
            if g_fitness > g_fitness2:
                gbest = X[p_fitness2.argmin()]
                g_fitness = g_fitness2
            # 記錄最優(yōu)迭代記錄
            fitness_val_list.append(g_fitness)
            i += 1

    # 輸出迭代結(jié)果
    print("最優(yōu)值是：%.5f" % fitness_val_list[-1])
    print("最優(yōu)解是：x=%.5f,y=%.5f" % (gbest[0], gbest[1]))

    # 繪圖
    plt.plot(fitness_val_list, color='r')
    plt.title('迭代過(guò)程')
    plt.show()

pos()