1. 程式人生 > 程式設計 >python實現粒子群演算法

python實現粒子群演算法

粒子群演算法

粒子群演算法源於複雜適應系統(Complex Adaptive System,CAS)。CAS理論於1994年正式提出,CAS中的成員稱為主體。比如研究鳥群系統,每個鳥在這個系統中就稱為主體。主體有適應性,它能夠與環境及其他的主體進行交流,並且根據交流的過程“學習”或“積累經驗”改變自身結構與行為。整個系統的演變或進化包括:新層次的產生(小鳥的出生);分化和多樣性的出現(鳥群中的鳥分成許多小的群);新的主題的出現(鳥尋找食物過程中,不斷髮現新的食物)。

PSO初始化為一群隨機粒子(隨機解)。然後通過迭代找到最優解。在每一次的迭代中,粒子通過跟蹤兩個“極值”(pbest,gbest)來更新自己。

在找到這兩個最優值後,粒子通過下面的公式來更新自己的速度和位置。

python實現粒子群演算法

python實現粒子群演算法

i 表示第 i 個粒子, d 表示粒子的第 d 個維度。r1,r2 表示兩個位於 [0,1] 的隨機數(對於一個粒子的不同維度,r1,r2 的值不同)。pbest[i] 是指粒子取得最高(低)適應度時的位置,gbest[i] 指的是整個系統取得最高(低)適應度時的位置。

實踐

我們用 PSO 演算法求解如下函式的最小值

python實現粒子群演算法

可以在空間畫出影象

python實現粒子群演算法

下圖是使用 5 個粒子的收斂情況

python實現粒子群演算法

可以看到,fitness 在第 12 輪就幾乎收斂到 -10.0。

下面是完整程式碼

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


INF = 1e5

def plot_cost_func():
  """畫出適應度函式"""
  fig = plt.figure()
  ax = Axes3D(fig)
  X = np.arange(-4,4,0.25)
  Y = np.arange(-4,0.25)
  X,Y = np.meshgrid(X,Y)
  Z = (X**2 + Y**2) - 10
  ax.plot_surface(X,Y,Z,rstride=1,cstride=1,cmap='rainbow')
  plt.show()

def fitness(x):
  return x[0]**2 + x[1]**2 - 10

class PSOSolver(object):
  def __init__(self,n_iter,weight=0.5,c1=2,c2=2,n_particle=5):
    self.n_iter = n_iter
    self.weight = weight
    self.c1 = c1
    self.c2 = c2
    self.n_particle = n_particle
    self.gbest = np.random.rand(2)
    # gbest 對應的函式值
    self.gbest_fit = fitness(self.gbest)
    # 將位置初始化到 [-5,5]
    self.location = 10 * np.random.rand(n_particle,2) - 5
    # 將速度初始化到 [-1,1]
    self.velocity = 2 * np.random.rand(n_particle,2) - 1
    self.pbest_fit = np.tile(INF,n_particle)
    self.pbest = np.zeros((n_particle,2))
    # 記錄每一步的最優值
    self.best_fitness = []
  
  def new_velocity(self,i):
    r = np.random.rand(2,2)
    v = self.velocity[i]
    x = self.location[i]
    pbest = self.pbest[i]
    return self.weight * v + self.c1 * r[0] * (pbest - x) + \
        self.c2 * r[1] * (self.gbest - x)

  def solve(self):
    for it in range(self.n_iter):
      for i in range(self.n_particle):
        v = self.new_velocity(i)
        x = self.location[i] + v
        fit_i = fitness(x)
        if fit_i < self.pbest_fit[i]:
          self.pbest_fit[i] = fit_i
          self.pbest[i] = x
          if fit_i < self.gbest_fit:
            self.gbest_fit = fit_i
            self.gbest = x
        self.velocity[i] = v
        self.location[i] = x
      self.best_fitness.append(self.gbest_fit)

  
if __name__ == '__main__':
  plot_cost_func()
  n_iter = 20
  s = PSOSolver(n_iter)
  s.solve()
  print(s.gbest_fit)
  plt.title("Fitness Curve")
  plt.xlabel("iter")
  plt.ylabel("fitness")
  plt.plot(np.arange(n_iter),np.array(s.best_fitness))
  plt.show()

以上就是python實現粒子群演算法的詳細內容,更多關於python 粒子群演算法的資料請關注我們其它相關文章!