[python] 啟發式演算法庫scikit-opt使用指北

落痕的寒假發表於2024-07-30

scikit-opt是一個封裝了多種啟發式演算法的Python程式碼庫，可以用於解決最佳化問題。scikit-opt官方倉庫見：scikit-opt，scikit-opt官網文件見：scikit-opt-doc。

scikit-opt安裝程式碼如下：

pip install scikit-opt

# 呼叫scikit-opt並檢視版本
import sko
sko.__version__

'0.6.6'

0 背景介紹

啟發式演算法介紹

啟發式演算法，顧名思義，就是一種基於直覺或經驗來解決問題的演算法。它不像傳統演算法那樣一步一步地窮盡所有可能性，而是透過一些啟發式的規則或策略，快速找到一個可行的解。打個比方，若開車去一個陌生的地方，沒有導航儀。啟發式演算法就像問路一樣，可以向路人詢問，也可以根據路邊的標誌和指示牌來判斷方向。雖然這種方式不能保證找到最優路線，但通常能夠在較短時間內找到一個可行的路線。

啟發式演算法通常具有以下特點：

快速性: 啟發式演算法通常能夠在較短時間內找到一個可行的解，特別是在面對複雜問題時。
魯棒性: 啟發式演算法對問題的細節不敏感，即使問題輸入發生變化，也能找到一個合理的解。
有效性: 啟發式演算法通常能夠找到高質量的解，在許多實際應用中表現良好。

當然，啟發式演算法也存在一些缺點：

近似性: 啟發式演算法不能保證找到最優解，只能找到一個較好的近似解。
經驗性: 啟發式演算法的效能很大程度上依賴於經驗或規則的設計，不同的啟發式演算法可能具有不同的效果。
侷限性: 啟發式演算法通常針對特定的問題設計，難以應用於其他問題。

儘管存在一些缺點，啟發式演算法仍然是一種非常有效的解決問題的方法，特別是在面對複雜問題或需要實時決策的情況下。啟發式演算法在許多領域都有應用，例如：

人工智慧: 啟發式演算法被廣泛應用於人工智慧領域，例如機器學習、機器人、遊戲等。
運籌最佳化: 啟發式演算法被用於解決各種運籌最佳化問題，例如旅行商問題、排程問題、資源分配問題等。
計算機圖形學: 啟發式演算法被用於解決計算機圖形學中的各種問題，例如路徑規劃、影像分割、紋理生成等。

scikit-opt中的啟發式演算法

scikit-opt支援的啟發式演算法包括：

差分進化演算法 (Differential Evolution)：一種基於群體搜尋的最佳化演算法，透過模擬生物進化的過程來尋找最優解。
遺傳演算法 (Genetic Algorithm)：模擬自然選擇和遺傳機制，透過種群中個體的變異、交叉和選擇來最佳化問題。
粒子群演算法 (Particle Swarm Optimization)：模擬鳥群或魚群中個體的群體行為，透過個體間資訊共享來搜尋最優解。
模擬退火演算法 (Simulated Annealing)：受金屬退火過程啟發，透過接受狀態降低的解的機率來在搜尋空間中跳出區域性最優解。
蟻群演算法 (Ant Colony Optimization)：模擬螞蟻尋找食物的行為，透過資訊素的沉積和蒸發來尋找最佳化路徑。
免疫最佳化演算法(Immune Optimization Algorithm)：基於免疫系統的進化過程，透過模擬抗體的生成和免疫反應來最佳化問題。
人工魚群演算法(Artificial Fish-Swarm Algorithm)：模擬魚群覓食行為，透過個體之間的位置調整和資訊共享來搜尋最優解。

這些演算法都屬於啟發式演算法的範疇，適用於複雜的最佳化問題，如函式最佳化、引數調優等。scikit-opt提供了這些演算法的Python實現，並且通常還包括了對不同問題的適應性調整和最佳化引數的支援，使得使用者能夠更方便地應用這些演算法進行問題求解。

以下表格包含了這些演算法的優缺點和適用環境：

演算法名稱	優點	缺點	適用環境
差分進化演算法	魯棒性強，適用於多種最佳化問題，收斂速度快	引數設定對演算法效能影響較大，容易陷入區域性最優	連續最佳化問題，特別是高維空間的最佳化問題
遺傳演算法	具有全域性搜尋能力，適用於複雜的最佳化問題	引數設定複雜，收斂速度較慢	組合最佳化問題，如旅行商問題、排程問題等
粒子群演算法	實現簡單，收斂速度快，適用於多模態最佳化問題	易陷入區域性最優，引數設定對演算法效能影響較大	連續最佳化問題，特別是多模態最佳化問題
模擬退火演算法	具有跳出區域性最優的能力，適用於複雜最佳化問題	引數設定複雜，收斂速度較慢	組合最佳化問題，如電路設計、影像處理等
蟻群演算法	適用於離散最佳化問題，具有分散式計算的特點	引數設定複雜，收斂速度較慢	組合最佳化問題，如旅行商問題、車輛路徑問題等
免疫最佳化演算法	具有較強的魯棒性和自適應性，適用於動態最佳化問題	演算法複雜度較高，引數設定較困難	組合最佳化問題，如特徵選擇、資料分類等
魚群演算法	實現簡單，具有較好的並行性，適用於多目標最佳化問題	易陷入區域性最優，引數設定對演算法效能影響較大	連續最佳化問題，特別是多目標最佳化問題

1 演算法使用

1.1 差分進化演算法

差分進化演算法（Differential Evolution, DE）是一種基於種群的啟發式最佳化演算法，其核心思想是透過模擬生物進化過程中的變異和自然選擇機制，不斷探索和改進解，最終找到問題的最優或近似最優解。這個過程可以簡單比喻為一次探險旅行：

探險隊伍的組建：在探險開始前，我們會組建一支隊伍，這支隊伍由多個探險者組成，每個探險者都代表了一個可能的解決方案。這些探險者攜帶著地圖和裝備，準備出發探索未知的領域。
差分（Differential）：在探險過程中，每個探險者會觀察周圍的環境，並與其他探險者交流資訊。差分在這裡可以理解為探險者之間的資訊交流，他們透過比較彼此的地圖和裝備，發現哪些是更有效的策略。在DE演算法中，差分是指透過計算兩個隨機選擇的探險者（解）之間的差異來生成新的探索方向。
進化（Evolution）：隨著探險的進行，探險者會根據收集到的資訊不斷調整自己的路線和策略。在DE演算法中，進化是指透過交叉（crossover）操作將新的探索方向與當前的解結合，生成新的候選解。這個過程類似於自然選擇，優秀的策略會被保留下來，而不適應的策略則會被淘汰。
適應度評估：探險者在探索過程中會不斷評估自己的路線是否有效，比如是否更接近目標或者是否避開了危險區域。在DE演算法中，適應度評估是指根據最佳化問題的特定目標來評價每個候選解的效能。
迭代更新：探險者在每次評估後，會根據結果更新自己的策略。在DE演算法中，這意味著根據適應度評估的結果，選擇保留或替換當前的解，從而不斷進化和改進。
目標達成：最終，探險者會找到一條通往目標的最佳路徑。在DE演算法中，這意味著找到了一個接近最優解的候選解。

上面的比喻可能不完全準確，差分進化演算法詳細理解見：差分進化演算法的基本概念和原理。

scikit-opt庫中的DE類用於實現差分進化演算法演算法，DE類構建引數介紹如下：

引數	預設值	意義
func	-	需要最佳化的函式
n_dim	-	目標函式包含的自變數個數
size_pop	50	種群規模，種群規模越大，探索範圍越廣，但計算量也越大
max_iter	200	演算法執行的最大迭代次數
prob_mut	0.001	變異機率，變異機率越高，種群多樣性越大，但也可能導致陷入區域性解
F	0.5	變異係數，變異係數越大，變異幅度越大，種群的探索範圍也越大
lb	-1	每個自變數的最小值
ub	1	每個自變數的最大值
constraint_eq	空元組	等式約束
constraint_ueq	空元組	不等式約束，形式為小於等於0

差分進化演算法對引數設定比較敏感，需要根據具體問題進行調整。scikit-opt中DE類使用例項如下：

定義最佳化問題

最佳化問題如下：

min f(x1, x2, x3) = x1 + x2^2 + x3^3
s.t.
    x1*x2 >= 1
    x1*x2 <= 5
    x2 + x3 = 1
    0 <= x1, x2, x3 <= 5

這是一個帶有約束條件的非線性規劃問題。在這個問題中，目標是找到一組變數 \(x_1, x_2, x_3\) 的值，使得目標函式 \(f(x_1, x_2, x_3) = x_1 + x_2^2 + x_3^3\) 達到最小值，同時滿足一系列給定的約束條件。

目標函式（Objective Function）：

\(f(x_1, x_2, x_3) = x_1 + x_2^2 + x_3^3\)

約束條件（Constraints）：

\(x_1 \cdot x_2 \geq 1\) （這是一個非線性約束，確保 \(x_1\) 和 \(x_2\) 的乘積不小於1）
\(x_1 \cdot x_2 \leq 5\) （同樣是非線性約束，限制 \(x_1\) 和 \(x_2\) 的乘積不大於5）
\(x_2 + x_3 = 1\) （這是一個線性等式約束，表示 \(x_2\) 和 \(x_3\) 的和必須為1）
\(0 \leq x_1, x_2, x_3 \leq 5\) （這是變數邊界約束，限制 \(x_1, x_2, x_3\) 的取值範圍在0到5之間）

示例程式碼

# 定義約束最佳化問題
def obj_func(p):
    # 給出輸入
    x1, x2, x3 = p
    # 返回目標函式值，目的是最小化目標函式值
    return x1 + x2 ** 2 + x3 ** 3


# 定義第三個約束條件
constraint_eq = [
    lambda x: 1 - x[1] - x[2]
]

# 定義第一個和第二個約束條件
constraint_ueq = [
    lambda x: 1 - x[0] * x[1],
    lambda x: x[0] * x[1] - 5
]


# 呼叫差分進化演算法解決問題
from sko.DE import DE

# lb和ub定義第四個約束條件
de = DE(func=obj_func, n_dim=3, size_pop=50, max_iter=100, lb=[0, 0, 0], ub=[5, 5, 5],
        constraint_eq=constraint_eq, constraint_ueq=constraint_ueq)

# 將返回最優解的變數值和目標函式值
# 目的是最小化目標函式值
best_x, best_y = de.run()
print('best_x:', best_x, '\n', 'best_y:', best_y)

best_x: [1.60211276 0.73970843 0.26021027] 
 best_y: [2.16689999]

帶入最優解的變數值到目標函式中的計算結果和best_y相等：

best_x[0]+best_x[1]**2+best_x[2]**3

2.1668999924186294

功能函式

# 每次迭代的最優函式值
res = de.generation_best_Y 
# 長度為max_iter
len(res)

# 每次迭代的最優函式值對應的輸入值
res = de.generation_best_X
# 長度為max_iter
len(res)

# 每次迭代種群所有個體的函式值
res = de.all_history_Y
# 長度為max_iter
# len(res)
# 返回值為size_pop，也就是種群個數
len(res[0])

1.2 遺傳演算法

1.2.1 基礎遺傳演算法

遺傳演算法(Genetic Algorithm，GA)是一種模擬自然選擇和遺傳機制的搜尋和最佳化演算法。它透過模擬生物進化過程中的遺傳、突變、交叉和選擇等自然現象，來解決最佳化問題。透過挖寶這個例子來通俗地說明遺傳演算法的核心思想：

定義問題（寶藏在哪裡）：假設你在一個廣闊的區域內尋找寶藏。這個區域可以看作是一個巨大的地圖，寶藏可能藏在任何位置。
初始化種群（隨機選擇起點）：開始時，你不知道寶藏的具體位置，因此你隨機選擇一些起點，這些起點可以看作是“種群”中的“個體”。每個個體代表一個可能的寶藏位置。
適應度評估（評估寶藏位置）：每個個體（寶藏位置）需要被評估其“適應度”。在挖寶的例子中，適應度可以定義為距離寶藏的真實位置的遠近。越接近寶藏的位置，適應度越高。
選擇（淘汰不合適的個體）：根據適應度，選擇那些更接近寶藏的個體進入下一代。這就像是自然選擇中，更適應環境的生物更有可能生存並繁衍後代。
交叉（產生新的寶藏位置）：將選擇出的個體進行“交叉”（類似於生物的交配）。在遺傳演算法中，這通常透過交換個體之間的某些特徵來實現。例如，兩個寶藏位置的座標可以相互交換，產生新的寶藏位置。
變異（引入新的寶藏位置）：為了增加種群的多樣性，引入一些隨機變化（變異）。這就像是生物基因突變一樣，可能會產生新的寶藏位置。
重複迭代（不斷嘗試）：重複上述過程，直到找到寶藏或者達到一定的迭代次數。每一次迭代都會產生新的種群，逐漸接近寶藏的真實位置。
收斂（找到寶藏）：最終，透過不斷的迭代和最佳化，種群中的個體將越來越接近寶藏的真實位置，直到找到寶藏。

上面的比喻可能不完全準確，遺傳演算法詳細理解見：10分鐘搞懂遺傳演算法。

scikit-opt庫中的GA類用於實現遺傳演算法，GA類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式
n_dim	-	目標函式的維度
size_pop	50	種群規模
max_iter	200	最大迭代次數
prob_mut	0.001	變異機率
lb	-1	每個自變數的最小值
ub	1	每個自變數的最大值
constraint_eq	空元組	等式約束
constraint_ueq	空元組	不等式約束
precision	1e-7	精準度，int/float或者它們組成的列表

案例 1

如下程式碼展示了利用遺傳演算法求解schaffer函式：

import numpy as np

# schaffer函式是一種用於最佳化演算法測試的標準基準函式
def schaffer(p):
    '''
    該函式有許多區域性最小值，並且存在強烈的震盪。
    全域性最小值出現在 (0,0) 處，其函式值為 0。
    '''
    # 解包傳入的引數 p，其中 p 是一個包含兩個元素的元組或陣列
    x1, x2 = p
    
    # 計算 x1 和 x2 的平方和
    x = np.square(x1) + np.square(x2)
    
    # 計算schaffer函式的值
    # 這裡使用了正弦函式和平方函式，並透過一個分母項來增加函式的複雜性
    return 0.5 + (np.square(np.sin(x)) - 0.5) / np.square(1 + 0.001 * x)

# 匯入遺傳演算法模組
from sko.GA import GA

# 初始化遺傳演算法
ga = GA(
    func=schaffer,  # 需要最佳化的目標函式，返回
    n_dim=2,        # 決策變數的維度
    size_pop=30,    # 種群大小
    max_iter=100,   # 最大迭代次數
    prob_mut=0.001, # 突變機率
    lb=[-1, -1],    # 決策變數的下界
    ub=[1, 1],      # 決策變數的上界
    precision=1e-7  # 演算法精度
)

# 將返回最優解的變數值和目標函式值
# 目的是最小化目標函式值
best_x, best_y = ga.run()

# 輸出最優解的決策變數和目標函式值
# 可以看到best_x中兩個變數值接近0
print('best_x:', best_x, '\n', 'best_y:', best_y)

best_x: [-0.0001221  0.0074937] 
 best_y: [5.9325642e-08]

# 匯入pandas和matplotlib.pyplot庫，用於資料處理和繪圖
import pandas as pd
import matplotlib.pyplot as plt

# 將遺傳演算法過程中的所有個體的函式值歷史記錄轉換為DataFrame
# ga.all_history_Y為每次迭代種群所有個體的函式值
# Y_history形狀為(max_iter,size_pop)
Y_history = pd.DataFrame(ga.all_history_Y)

# 建立一個包含兩個子圖的圖形
fig, ax = plt.subplots(2, 1)

# 在第一個子圖中繪製函式值的歷史記錄
ax[0].plot(Y_history.index, Y_history.values, '.', color='red')

# 在第二個子圖中繪製每次迭代的各種群的函式最小值
Y_history.min(axis=1).cummin().plot(kind='line')

# 顯示圖形
plt.show()

png

GA類還提供了其他的功能函式：

ga.generation_best_Y # 每一代的最優函式值
ga.generation_best_X # 每一代的最優函式值對應的輸入值
ga.all_history_FitV # 每一代的每個個體的適應度
ga.all_history_Y # 每一代每個個體的函式值

案例 2

如下程式碼展示了利用precision引數控制變數的精度。當precision引數被設定為整數時，系統自動啟用整數規劃模式，使得變數值嚴格遵循整數約束。在整數規劃模式下，為了達到最佳的收斂速度與效果，推薦變數的可能取值數量儘量為\(2^n\)的形式：

from sko.GA import GA

demo_func = lambda x: (x[0] - 1) ** 2 + (x[1] - 0.05) ** 2 + x[2] ** 2
ga = GA(func=demo_func, n_dim=3, max_iter=500, lb=[-1, -1, -1], ub=[5, 1, 1], precision=[1, 2, 1e-7])
best_x, best_y = ga.run()
print('best_x:', best_x, '\n', 'best_y:', best_y)

best_x: [1.00000000e+00 1.00000000e+00 2.98023233e-08] 
 best_y: [0.9025]

案例 3

以下程式碼展示了使用遺傳演算法進行曲線擬合：

import numpy as np
# 匯入matplotlib.pyplot用於繪圖
import matplotlib.pyplot as plt
# 匯入sko庫中的遺傳演算法模組
from sko.GA import GA

# 建立資料點x_true，範圍從-1.2到1.2，共20個點
x_true = np.linspace(-1.2, 1.2, 20)
# 計算對應的y_true值，為三次多項式的值加上隨機噪聲
y_true = x_true ** 3 - x_true + 0.4 * np.random.rand(20)

# 定義三次多項式函式f_fun，引數包括x值和係數a, b, c, d
def f_fun(x, a, b, c, d):
    return a * x ** 3 + b * x ** 2 + c * x + d

# 定義遺傳演算法的適應度函式obj_fun，用於計算多項式擬合的殘差平方和
def obj_fun(p):
    a, b, c, d = p  # 解包引數
    # 計算擬合多項式與原始資料點的殘差平方和
    residuals = np.square(f_fun(x_true, a, b, c, d) - y_true).sum()
    return residuals  # 返回殘差平方和

# 初始化遺傳演算法，設定適應度函式、維度、種群大小、迭代次數、引數界限
ga = GA(func=obj_fun, n_dim=4, size_pop=100, max_iter=500,
        lb=[-2] * 4, ub=[2] * 4)

# 執行遺傳演算法，獲取最優引數和對應的殘差
best_params, residuals = ga.run()
# 列印最優引數和殘差
print('best_x:', best_params, '\n', 'best_y:', residuals)

# 使用最優引數計算預測的y值
y_predict = f_fun(x_true, *best_params)

# 建立繪圖視窗和座標軸
fig, ax = plt.subplots()

# 在同一個座標軸上繪製原始資料點和預測的曲線
ax.plot(x_true, y_true, 'o')
ax.plot(x_true, y_predict, '-')

# 顯示繪圖視窗
plt.show()

best_x: [ 0.82131277 -0.01955447 -0.86175541  0.16721198] 
 best_y: [0.21969302]

png

1.2.2 用於解決旅行商問題的遺傳演算法

scikit-opt提供了GA_TSP模組以專門為解決旅行商問題（Traveling Salesman Problem, TSP）而設計的遺傳演算法。它透過對問題的特定部分進行過載，比如交叉（crossover）和變異（mutation）操作，來適應TSP問題的特點。使用GA_TSP時，首先需要定義問題的座標和距離矩陣，然後建立一個GA_TSP物件並呼叫run方法來執行演算法，最終得到最短路徑和對應的總距離。

所謂是TSP問題一個經典的組合最佳化問題，通常描述為：給定一組城市和各城市之間的距離，求解一條訪問每個城市一次並且回到起始城市的最短路徑。這個問題在電腦科學和運籌學中非常著名，也是NP-hard問題的代表之一，意味著隨著城市數量的增加，找到最優解的計算複雜度呈指數增長。TSP的解決方案不僅限於商業旅行推銷員的路線最佳化，也可以應用於諸如電路板生產、DNA測序等各種現實生活中的路徑規劃問題。

GA_TSP的輸入引數如下，可控制引數比GA演算法少：

引數	預設值	意義
func	-	目標函式
n_dim	-	城市個數
size_pop	50	種群規模
max_iter	200	最大迭代次數
prob_mut	0.001	變異機率

案例 1

以下是使用scikit-opt中的GA_TSP類來求解TSP問題的簡單示例：

import numpy as np 
from scipy import spatial  # 匯入scipy庫中的spatial模組，用於計算空間距離  
import matplotlib.pyplot as plt  
  
# 設定點的數量  
num_points = 20  
  
# 生成隨機點座標  
points_coordinate = np.random.rand(num_points, 2)  # 生成一個num_points x 2的陣列，包含隨機座標  
  
# 計算點之間的距離矩陣，使用歐幾里得距離  
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric='euclidean')  
  
# 定義目標函式，計算給定路線（routine）的總距離  
def cal_total_distance(routine):  
    '''  
    目標函式。輸入一個路線（routine），返回一個總距離。  
    route是一個一維陣列，表示訪問點的順序
    '''  
    num_points, = routine.shape  # 獲取路線陣列的長度  
    # 使用列表推導式和距離矩陣計算總距離  
    # 遍歷路線中的每個點對，計算相鄰點之間的距離並求和  
    return sum([distance_matrix[routine[i % num_points], routine[(i + 1) % num_points]] for i in range(num_points)])  
  
# 從sko庫匯入遺傳演算法（GA）解決旅行商問題（TSP）的類  
from sko.GA import GA_TSP  
  
# 初始化GA_TSP物件  
ga_tsp = GA_TSP(func=cal_total_distance, n_dim=num_points, size_pop=50, max_iter=100, prob_mut=1)  

# 執行遺傳演算法，得到最優解  
best_points, best_distance = ga_tsp.run()  
print(best_points)
print(best_distance) # 最小距離

[13 15  1 11 16  7  3  9  5 12 14 19  6 17  2 10  8 18  4  0]
[3.91036766]

視覺化程式碼如下：

# 繪圖展示結果  
fig, ax = plt.subplots(1, 2)  # 建立一個包含兩個子圖的圖形  
# 繪製最優路線  
best_points_ = np.concatenate([best_points, [best_points[0]]])  # 將起點新增到終點，形成閉環  
best_points_coordinate = points_coordinate[best_points_, :]  # 獲取最優路線上的點座標  
ax[0].plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1], 'o-r')  # 繪製最優路線  
# 新增箭頭表示路線方向
for i in range(num_points):
    ax[0].annotate('', xy=best_points_coordinate[i + 1], xytext=best_points_coordinate[i],
                   arrowprops=dict(arrowstyle='->', color='green'))
ax[0].plot(best_points_coordinate[0, 0], best_points_coordinate[0, 1], 'gs')  # 使用綠色方塊標記起點    

# 繪製演算法迭代過程中的最優解變化  
ax[1].plot(ga_tsp.generation_best_Y)  # generation_best_Y儲存了每代的最優解  
plt.show()  # 顯示圖形

png

案例 2

下面的示例展示了遺傳TSP問題中固定起點和終點的方法。假設起點和終點的座標分別指定為(0, 0)和(1, 1)，考慮共有n+2個點，最佳化目標是中間的n個點，而起點和終點則固定不參與最佳化。目標函式即為實際路徑的總距離：

import numpy as np 
from scipy import spatial  # 匯入scipy庫中的spatial模組，用於計算空間距離  
import matplotlib.pyplot as plt  
  
# 設定點的數量  
num_points = 15  
  
# 生成隨機點座標  
points_coordinate = np.random.rand(num_points, 2)  # 生成一個num_points x 2的陣列，包含隨機座標  
  
start_point=[[0,0]] # 起始點
end_point=[[1,1]] # 結束點
points_coordinate = np.concatenate([points_coordinate,start_point,end_point])
# 計算點之間的距離矩陣，使用歐幾里得距離  
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric='euclidean') 
  
# 定義目標函式，計算給定路線（routine）的總距離  
def cal_total_distance(routine):  
    '''  
    目標函式。輸入一個路線（routine），返回一個總距離。  
    route是一個一維陣列，表示訪問點的順序
    '''  
    num_points, = routine.shape
    # start_point,end_point 本身不參與最佳化。給一個固定的值，參與計算總路徑
    # num_points，num_points+1為start_point,end_point的標號
    routine = np.concatenate([[num_points], routine, [num_points+1]])
    # 遍歷路線中的每個點對，計算相鄰點之間的距離並求和  
    return sum([distance_matrix[routine[i], routine[i + 1]] for i in range(num_points+2-1)])

# 從sko庫匯入遺傳演算法（GA）解決旅行商問題（TSP）的類  
from sko.GA import GA_TSP  
  
# 初始化GA_TSP物件  
ga_tsp = GA_TSP(func=cal_total_distance, n_dim=num_points, size_pop=50, max_iter=100, prob_mut=1)  

# 執行遺傳演算法，得到最優解  
best_points, best_distance = ga_tsp.run()  
  
# 繪圖展示結果  
fig, ax = plt.subplots(1, 2)  # 建立一個包含兩個子圖的圖形  
# 繪製最優路線  
best_points_ = np.concatenate([best_points, [best_points[0]]])  # 將起點新增到終點，形成閉環  
best_points_coordinate = points_coordinate[best_points_, :]  # 獲取最優路線上的點座標  
ax[0].plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1], 'o-r')  # 繪製最優路線  
# 新增箭頭表示路線方向
for i in range(num_points):
    ax[0].annotate('', xy=best_points_coordinate[i + 1], xytext=best_points_coordinate[i],
                   arrowprops=dict(arrowstyle='->', color='green'))
ax[0].plot(best_points_coordinate[0, 0], best_points_coordinate[0, 1], 'gs')  # 使用綠色方塊標記起點    

# 繪製演算法迭代過程中的最優解變化  
ax[1].plot(ga_tsp.generation_best_Y)  # generation_best_Y儲存了每代的最優解  
plt.show()  # 顯示圖形

png

1.3 粒子群演算法

粒子群最佳化演算法（Particle Swarm Optimization, PSO）是一種模擬鳥群覓食行為的最佳化演算法。想象一群鳥在廣闊的天空中尋找食物，它們透過觀察彼此的位置和移動方向來尋找食物最多的地方。粒子群演算法也是基於類似的原理，下面透過探險的例子來通俗解釋其核心理念：

探險隊員：在粒子群演算法中，每個“粒子”可以想象成一名探險隊員。每個隊員都有一定的位置，代表了一個可能的解決方案。
位置和速度：每個探險隊員都有自己當前的位置（解決方案）和速度（向目標移動的方向和速度）。
個人經驗：每個探險隊員都會記得自己曾經到達過的最好位置，也就是他們個人找到的最豐富的食物源。
群體經驗：除了個人經驗外，探險隊員們還會共享整個隊伍找到的最好位置，這是整個群體的共識，代表了一個更優的解決方案。
探索與利用：探險隊員們在移動時會考慮個人經驗（自己找到的最好位置）和群體經驗（整個隊伍找到的最好位置）。他們會根據這兩個因素調整自己的速度和方向，既探索新的地方，也利用已知的資訊。
更新位置：根據速度和方向，探險隊員們更新自己的位置，繼續尋找更好的解決方案。
迭代過程：這個過程會不斷重複，每次迭代都是一次新的探索，探險隊員們會根據新的位置更新個人經驗和群體經驗。
收斂：隨著時間的推移，探險隊員們會逐漸聚集到寶藏最多的地方，這代表演算法找到了最優解或接近最優解的位置。

上面的比喻可能不完全準確，粒子群演算法詳細理解見：粒子群最佳化演算法的詳細解讀。

scikit-opt庫中的PSO類用於實現粒子群最佳化演算法，PSO類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式
n_dim	-	目標函式的維度
size_pop	50	種群規模
max_iter	200	最大迭代次數
lb	None	每個引數的最小值
ub	None	每個引數的最大值
w	0.8	慣性權重，控制粒子運動速度對上次速度的影響程度。w值越大，粒子對歷史速度的記憶力越強，全域性搜尋能力越強
c1	0.5	個體記憶，控制粒子向自身歷史最優位置移動的程度。c1值越大，粒子越傾向於向自身歷史最優位置移動
c2	0.5	集體記憶，控制粒子向全域性最優位置移動的程度。c2值越大，粒子越傾向於向全域性最優位置移動
constraint_ueq	空元組	不等式約束，形式為小於等於0

案例 1

以下程式碼展示了PSO類簡單使用示例：

# 定義問題
def demo_func(x):
    """
    目標函式。
    輸入：
        x: 一個包含三個元素的列表，代表問題的三個變數。
    輸出：
        目標函式值：最佳化目標尉最小化該值
    """
    x1, x2, x3 = x  # 將列表x解包為三個變數
    return x1 ** 2 + (x2 - 0.01) ** 2 + x3 ** 2  # 計算目標函式值

# 呼叫粒子群演算法
from sko.PSO import PSO

pso = PSO(
    func=demo_func,  # 目標函式
    n_dim=3,  # 問題的維度（變數個數）
    pop=40,  # 種群大小（粒子個數）
    max_iter=150,  # 最大迭代次數
    lb=[0, -1, 0.2],  # 變數的下界
    ub=[1, 1, 1],  # 變數的上界
    w=0.8,  # 慣性權重
    c1=0.5,  # 個體學習因子
    c2=0.5  # 集體學習因子
)
pso.run()  # 執行粒子群演算法

print('best_x is ', pso.gbest_x, 'best_y is', pso.gbest_y)  # 列印最優解

import matplotlib.pyplot as plt

plt.plot(pso.gbest_y_hist)  # 繪製目標函式值隨迭代次數的變化曲線
plt.show()

best_x is  [0.   0.01 0.2 ] best_y is [0.04]

png

案例 2

以下程式碼展示了帶非線性約束的PSO類使用示例及視覺化程式碼（注意該程式碼執行時間較長）：

import numpy as np  
from sko.PSO import PSO  
  
def demo_func(x):
    x1, x2 = x
    return x1 ** 2 + (x2 - 0.05) ** 2 
  
# 定義非線性約束
constraint_ueq = (  
  lambda x: (x[0] - 1) ** 2 + (x[1] - 0) ** 2 - 0.5 ** 2  
  ,  
)  
  
# 設定粒子群最佳化演算法的最大迭代次數  
max_iter = 50  
# 建立PSO例項，設定目標函式、維度、種群大小、迭代次數、變數範圍以及非線性約束
pso = PSO(func=demo_func, n_dim=2, pop=40, max_iter=max_iter, lb=[-2, -2], ub=[2, 2],  
          constraint_ueq=constraint_ueq)  
# 開啟記錄模式，以便後續繪製動態圖  
pso.record_mode = True  
# 執行PSO演算法  
pso.run()  
# 列印最優解的位置和值  
print('best_x is ', pso.gbest_x, 'best_y is', pso.gbest_y)  
  
# 匯入matplotlib庫進行繪圖  
import matplotlib.pyplot as plt  
from matplotlib.animation import FuncAnimation  
  
# 從PSO記錄中獲取資料  
record_value = pso.record_value  
X_list, V_list = record_value['X'], record_value['V']  
  
# 建立圖形和座標軸  
fig, ax = plt.subplots(1, 1)  
ax.set_title('title', loc='center')  # 暫時設定標題，後續更新  
line = ax.plot([], [], 'b.')  # 繪製初始空點  
  
# 建立網格並計算目標函式值  
X_grid, Y_grid = np.meshgrid(np.linspace(-2.0, 2.0, 40), np.linspace(-2.0, 2.0, 40))  
Z_grid = demo_func((X_grid, Y_grid))  
ax.contour(X_grid, Y_grid, Z_grid, 30)  # 繪製等高線圖  
  
# 設定座標軸範圍  
ax.set_xlim(-2, 2)  
ax.set_ylim(-2, 2)  
  
# 繪製圓形約束的邊界  
t = np.linspace(0, 2 * np.pi, 40)  
ax.plot(0.5 * np.cos(t) + 1, 0.5 * np.sin(t), color='r')  
  
# 定義更新散點圖的函式  
def update_scatter(frame):  
  # 計算當前幀對應的迭代次數和子幀數  
  i, j = frame // 10, frame % 10  
  # 更新標題  
  ax.set_title('iter = ' + str(i))  
  # 計算當前粒子位置（考慮到速度和迭代）  
  X_tmp = X_list[i] + V_list[i] * j / 10.0  
  # 更新散點圖的資料  
  plt.setp(line, 'xdata', X_tmp[:, 0], 'ydata', X_tmp[:, 1])  
  return line  
  
# 建立動畫  
ani = FuncAnimation(fig, update_scatter, blit=True, interval=25, frames=max_iter * 10)  
# 儲存動畫為gif檔案  
ani.save('pso.gif', writer='pillow')

best_x is  [0.50070802 0.02494922] best_y is [0.25133606]

png

1.4 模擬退火演算法

1.4.1 基礎模擬退火演算法

模擬退火演算法 (Simulated Annealing，SA)來源於冶金學中的退火過程，這是一種透過加熱和緩慢冷卻金屬來減少其內部缺陷的方法。在演算法中，這個過程被用來尋找問題的最優解。核心理念可以用探險來通俗說明：

隨機出發：想象你是一個探險家，開始時你站在一個未知的森林中，不知道寶藏在哪裡。
探索周邊：你開始隨機走動，探索周圍的區域，這就像是演算法在解空間中隨機選擇候選解。
接受好壞：在探險中，你可能會遇到好地方（離寶藏更近）或壞地方（離寶藏更遠）。模擬退火演算法允許你在早期階段接受好壞解，這就像是即使知道某個方向可能不是最佳路徑，但為了探索未知，你還是決定去嘗試。
溫度降低：隨著探險的進行，你會逐漸減少探索的範圍，變得更加謹慎。在演算法中，這相當於降低“溫度引數”，意味著演算法在選擇新解時變得更加挑剔。
最佳解：最終你可能會找到一個寶藏，或者至少是一個相對較好的地方。演算法透過這種方式逐步收斂到問題的最優解或近似最優解。
結束探險：當溫度降到足夠低，你幾乎不再移動，演算法也就停止搜尋，此時你所在的地點就是演算法找到的解。

模擬退火演算法的關鍵特點在於它能夠在搜尋過程中接受次優解，這有助於演算法跳出區域性最優解，從而有更大機會找到全域性最優解。隨著“溫度”的降低，演算法逐漸變得更加專注於尋找更好的解，直到最終收斂。

上面的比喻可能不完全準確，模擬退火演算法演算法詳細理解見：模擬退火演算法詳解。scikit-opt庫中的SA類用於實現模擬退火演算法，SA類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式
x0	-	迭代初始點，演算法開始搜尋時的初始自變數值
T_max	100	最大溫度，模擬退火演算法中的初始溫度，通常設定為一個較高的值以允許較大的解空間探索
T_min	1e-7	最小溫度度，演算法停止時的溫度閾值，當溫度降至此值以下時，演算法停止
L	300	鏈長，模擬退火演算法中鏈的步數或長度
max_stay_counter	150	冷卻耗時，當達到此次數時，無論是否找到更好的解，都將降低溫度並繼續搜尋
lb		每個自變數的最小值
ub		每個自變數的最大值

案例 1

以下程式碼展示瞭如何將模擬退火演算法用於多元函式最佳化：

demo_func = lambda x: x[0] ** 2 + (x[1] - 0.02) ** 2 + x[2] ** 2

from sko.SA import SA

sa = SA(func=demo_func, x0=[1, 1, 1], T_max=1, T_min=1e-9, L=300, max_stay_counter=10)
best_x, best_y = sa.run()
print('best_x:', best_x, 'best_y', best_y)

import matplotlib.pyplot as plt
import pandas as pd

# de.generation_best_Y 每一代的最優函式值
# de.generation_best_X 每一代的最優函式值對應的輸入值
plt.plot(pd.DataFrame(sa.best_y_history).cummin(axis=0))
plt.show()

best_x: [7.21355253e-06 2.00035659e-02 1.69597500e-05] best_y 3.5238375871831515e-10

png

案例 2

scikit-opt庫的SA模組還包括了"Fast"、"Cauchy"和"Boltzmann"三種不同的溫度下降策略或機率分佈。關於這些SA演算法的詳細介紹，請參考：從模擬退火到退火進化演算法。其中SA模組的SA類預設使用"Fast"版本演算法。

Fast、Cauchy和Boltzmann各自有不同的優缺點：

Fast（快速冷卻）：

優點：冷卻速度快，可以在較短的時間內達到較低的溫度，從而快速收斂到一個解。
缺點：由於冷卻速度過快，可能導致演算法過早地陷入區域性最優解，而無法探索到全域性最優解。

Cauchy（柯西分佈）：

優點：使用柯西分佈作為機率分佈來選擇新解，可以更好地避免區域性最優，因為柯西分佈的尾部較重，有助於探索更遠的解空間。
缺點：由於柯西分佈的特性，演算法可能會在某些情況下過於集中在某些區域，導致搜尋效率降低。

Boltzmann（玻爾茲曼分佈）：

優點：使用玻爾茲曼分佈作為機率分佈，可以平衡探索和利用之間的關係，隨著溫度的降低，演算法逐漸偏向於選擇更優的解。
缺點：溫度下降策略需要精心設計，以避免過快或過慢的冷卻速度，這可能會影響演算法的效能。

以上三種策略都有其適用的場景和問題型別。例如，對於需要快速得到一個解的問題，Fast策略可能更合適；而對於需要避免陷入區域性最優的問題，Cauchy或Boltzmann策略可能更有優勢。在實際應用中，選擇哪種策略往往需要根據具體問題的特性和要求來決定。示例程式碼如下：

# 定義一個匿名函式，用於計算目標函式的值
# 該函式接受一個列表 x 作為引數，返回 x 中元素的平方和
demo_func = lambda x: x[0] ** 2 + (x[1] - 0.02) ** 2 + x[2] ** 2

# 從 sko.SA 模組匯入 SAFast 類，用於快速模擬退火演算法
from sko.SA import SAFast
# 建立 SAFast 類的例項 sa_fast，用於執行模擬退火演算法
# 引數說明：
# func: 目標函式
# x0: 初始解
# T_max: 最大溫度
# T_min: 最小溫度
# q: 溫度衰減因子
# L: 鄰域搜尋次數
# lb: 下界列表
# ub: 上界列表
# max_stay_counter: 最大停滯次數
# SA等同於呼叫SAFast函式
sa_fast = SAFast(func=demo_func, x0=[1, 1, 1], T_max=1, T_min=1e-9, q=0.99, L=300, max_stay_counter=150,
                 lb=[-1, 1, -1], ub=[2, 3, 4])
# 執行模擬退火演算法
sa_fast.run()
# 列印演算法找到的最佳解和對應的目標函式值
print('Fast Simulated Annealing: best_x is ', sa_fast.best_x, 'best_y is ', sa_fast.best_y)

# 從 sko.SA 模組匯入 SABoltzmann 類，用於 Boltzmann 模擬退火演算法
from sko.SA import SABoltzmann
# 建立 SABoltzmann 類的例項 sa_boltzmann，用於執行 Boltzmann 模擬退火演算法
# 其餘引數與 SAFast 類似，但使用了不同的溫度更新策略
sa_boltzmann = SABoltzmann(func=demo_func, x0=[1, 1, 1], T_max=1, T_min=1e-9, q=0.99, L=300, max_stay_counter=150,
                           lb=-1, ub=[2, 3, 4])
# 列印演算法找到的最佳解和對應的目標函式值，注意這裡列印的 best_y 應該是 sa_boltzmann.best_y 而不是 sa_fast.best_y
print('Boltzmann Simulated Annealing: best_x is ', sa_boltzmann.best_x, 'best_y is ', sa_boltzmann.best_y)  

# 從 sko.SA 模組匯入 SACauchy 類，用於 Cauchy 模擬退火演算法
from sko.SA import SACauchy
# 建立 SACauchy 類的例項 sa_cauchy，用於執行 Cauchy 模擬退火演算法
# 其餘引數與 SABoltzmann 類似，但使用了不同的接受準則
sa_cauchy = SACauchy(func=demo_func, x0=[1, 1, 1], T_max=1, T_min=1e-9, q=0.99, L=300, max_stay_counter=150,
                     lb=[-1, 1, -1], ub=[2, 3, 4])
sa_cauchy.run()
# 列印演算法找到的最佳解和對應的目標函式值
print('Cauchy Simulated Annealing: best_x is ', sa_cauchy.best_x, 'best_y is ', sa_cauchy.best_y)

Fast Simulated Annealing: best_x is  [2.58236441e-06 1.00000000e+00 2.91945906e-06] best_y is  0.9604000000151918
Boltzmann Simulated Annealing: best_x is  [1 1 1] best_y is  2.9604
Cauchy Simulated Annealing: best_x is  [-7.75482132e-04  1.00000000e+00  4.41811491e-05] best_y is  0.9604006033245103

1.4.2 用於解決旅行商問題的模擬退火演算法

scikit-opt提供了SA_TSP模組以專門為解決旅行商問題（Traveling Salesman Problem, TSP）而設計的遺傳演算法。SA_TSP的輸入引數如下，可控制引數比SA類少：

引數	預設值	意義
func	-	目標函式
x0	-	迭代初始點
T_max	100	最大溫度
T_min	1e-7	最小溫度
L	300	鏈長
max_stay_counter	150	冷卻耗時

示例程式碼如下：

import numpy as np 
# 匯入scipy庫中的spatial模組，用於計算空間距離
from scipy import spatial  
import matplotlib.pyplot as plt  

# 設定點的數量
num_points = 15  

# 生成隨機點座標，生成一個num_points x 2的陣列，包含隨機的x和y座標
points_coordinate = np.random.rand(num_points, 2)  

# 定義起始點和結束點的座標
start_point = [[0,0]]  # 起始點座標
end_point = [[1,1]]    # 結束點座標

# 將起始點、隨機點座標和結束點合併為一個陣列
points_coordinate = np.concatenate([points_coordinate, start_point, end_point])

# 計算點之間的距離矩陣，使用歐幾里得距離
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric='euclidean') 

# 定義目標函式，計算給定路線（routine）的總距離
def cal_total_distance(routine):  
    '''  
    目標函式。輸入一個路線（routine），返回一個總距離。
    route是一個一維陣列，表示訪問點的順序。
    '''
    num_points, = routine.shape  # 獲取路線中點的數量
    # 將起始點和結束點的索引新增到路線陣列中
    routine = np.concatenate([[num_points], routine, [num_points+1]])
    # 計算路線中相鄰點之間的距離，並累加得到總距離
    return sum([distance_matrix[routine[i], routine[i + 1]] for i in range(num_points + 1)])

# 匯入模擬退火演算法的TSP求解器
from sko.SA import SA_TSP

# 初始化模擬退火TSP求解器
sa_tsp = SA_TSP(func=cal_total_distance, x0=range(num_points), T_max=100, T_min=1, L=10 * num_points)

# 執行求解器，找到最佳路線和對應的總距離
best_points, best_distance = sa_tsp.run()

# 列印最佳路線和對應的總距離
print(best_points, best_distance, cal_total_distance(best_points))

# 準備繪圖
from matplotlib.ticker import FormatStrFormatter

# 建立一個1x2的子圖
fig, ax = plt.subplots(1, 2, figsize=(12,6))

# 將最佳路線的點索引擴充套件，包括起始點和結束點
best_points_ = np.concatenate([best_points, [best_points[0]]])
# 獲取最佳路線對應的座標點
best_points_coordinate = points_coordinate[best_points_, :]

# 在第一個子圖中繪製模擬退火演算法的最優解歷史
ax[0].plot(sa_tsp.best_y_history)
ax[0].set_xlabel("Iteration")
ax[0].set_ylabel("Distance")

# 在第二個子圖中繪製最佳路線
ax[1].plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1],
           marker='o', markerfacecolor='white', color='red', linestyle='-')

# 新增箭頭表示路線方向
for i in range(num_points):
    ax[1].annotate('', xy=best_points_coordinate[i + 1], xytext=best_points_coordinate[i],
                   arrowprops=dict(arrowstyle='->', color='green'))

# 使用綠色方塊標記起點
ax[1].plot(best_points_coordinate[0, 0], best_points_coordinate[0, 1], 'gs')

# 設定座標軸的格式
ax[1].xaxis.set_major_formatter(FormatStrFormatter('%.3f'))
ax[1].yaxis.set_major_formatter(FormatStrFormatter('%.3f'))

# 設定座標軸的標籤
ax[1].set_xlabel("Longitude")
ax[1].set_ylabel("Latitude")

# 顯示圖形
plt.show()

[ 8 12  4  2  3  6  9 13  1  7 11  5  0 14 10] 5.6269966435102825 5.6269966435102825

png

1.5 蟻群演算法

蟻群演算法(Ant Colony Algorithm，ACA)是一種模擬螞蟻尋找食物的行為的演算法，用於解決路徑最佳化問題。靈感來源於螞蟻尋找最短路徑到達食物源的現象。下面用探險來解釋蟻群演算法的核心理念：

探險準備：設想一群探險者要穿越一片未知的森林尋找寶藏。他們不知道寶藏的確切位置，只能透過試錯來探索路徑。
資訊素的發現：探險者們在走動時，會留下一種特殊的“資訊素”。這種資訊素會隨著時間的推移逐漸消失，但它會吸引其他探險者朝這個方向前進。
選擇路徑：當探險者站在一個交叉路口時，他們傾向於選擇資訊素濃度較高的路徑，因為這些路徑可能已經被其他探險者探索並確認為較優的路徑。
區域性探索與全域性最佳化：探險者在探索過程中會進行區域性探索，即跟隨資訊素濃度高的路徑前進，但同時演算法也允許他們進行隨機探索，以避免陷入區域性最優解。
更新資訊素：每當探險者找到一條新的路徑到達寶藏，他們就會在這條路徑上留下更多的資訊素。這樣，其他探險者在遇到這個交叉路口時，更可能選擇這條新的、被更多資訊素標記的路徑。
資訊素的更新規則：資訊素的更新不僅取決於找到新路徑的探險者，還取決於他們找到寶藏的速度。如果探險者找到了更快的路徑，那麼他們留下的資訊素會更多，從而更快地影響其他探險者的選擇。
收斂與穩定：隨著時間的推移，探險者們會逐漸收斂到最優或近似最優的路徑上，因為資訊素的不斷更新使得所有探險者都傾向於選擇越來越短的路徑。
演算法的終止：探險結束的條件可以是找到寶藏，或者所有探險者都集中在一條路徑上，或者達到一定的探索次數。

上面的比喻可能不完全準確，蟻群演算法詳細理解見：一文搞懂什麼是蟻群最佳化演算法。
透過這種模擬螞蟻尋找食物的策略，蟻群演算法能夠有效地解決旅行商問題（TSP）和其他最佳化問題，尤其是在問題規模較大且求解空間複雜時。因此在scikit-opt庫中的蟻群演算法主要用於求解旅行商問題。ACA_TSP類用於實現蟻群演算法演算法，ACA_TSP類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式
n_dim	-	城市個數
size_pop	10	螞蟻數量，更多的螞蟻可能會提高找到全域性最優解的機率，但計算成本也會增加
max_iter	20	最大迭代次數
distance_matrix	-	城市之間的距離矩陣，其中每個元素代表兩個城市之間的距離，用於計算資訊素的揮發
alpha	1	資訊素重要程度，決定了資訊素在螞蟻選擇路徑時的影響。較高的alpha值意味著資訊素的影響更大
beta	2	適應度的重要程度，決定了路徑適應度（如距離或成本）在選擇路徑時的影響。較高的beta值意味著適應度的影響更大
rho	0.1	資訊素揮發速度，控制資訊素隨時間的揮發速度。較高的rho值意味著資訊素揮發得更快，從而減少資訊素對路徑選擇的影響

以下案例展示瞭如何利用ACA_TSP類求解旅行商問題：

# 匯入所需的庫
import numpy as np
from scipy import spatial
import pandas as pd
import matplotlib.pyplot as plt

# 設定點的數量
num_points = 20

# 生成點的座標
points_coordinate = np.random.rand(num_points, 2)  

# 計算點與點之間的歐幾里得距離矩陣
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric='euclidean')

# 定義計算旅行商問題總距離的函式
def cal_total_distance(routine):
    # 獲取路徑長度
    num_points, = routine.shape
    # 計算並返回路徑的總距離
    return sum([distance_matrix[routine[i % num_points], routine[(i + 1) % num_points]] for i in range(num_points)])

# 從 sko 庫匯入 ACA_TSP 演算法用於解決 TSP 問題
from sko.ACA import ACA_TSP

# 初始化 ACA_TSP 演算法，設定問題引數
# 注意直接執行以下程式碼會報錯：module 'numpy' has no attribute 'int'
# 需要將ACA_TSP函式報錯位置改為self.Table = np.zeros((size_pop, n_dim)).astype(np.int32)
aca = ACA_TSP(func=cal_total_distance, n_dim=num_points,
              size_pop=50, max_iter=200,
              distance_matrix=distance_matrix)

# 執行演算法並獲取最優解
best_x, best_y = aca.run()

# 建立繪圖視窗和軸
fig, ax = plt.subplots(1, 2)

# 將最優路徑的點連線起來，包括起點和終點
best_points_ = np.concatenate([best_x, [best_x[0]]])
best_points_coordinate = points_coordinate[best_points_, :]

# 在第一個子圖中繪製最優路徑
ax[0].plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1], 'o-r')

# 在最優路徑上新增箭頭，表示行進方向
for i in range(num_points):
    ax[0].annotate('', xy=best_points_coordinate[i + 1], xytext=best_points_coordinate[i],
                   arrowprops=dict(arrowstyle='->', color='green'))

# 用綠色方塊標記最優路徑的起點
ax[0].plot(best_points_coordinate[0, 0], best_points_coordinate[0, 1], 'gs')

# 在第二個子圖中繪製演算法效能隨迭代次數的變化
pd.DataFrame(aca.y_best_history).cummin().plot(ax=ax[1])

# 顯示圖表
plt.show()

png

1.6 免疫最佳化演算法

免疫最佳化演算法 (Immune Optimization Algorithm)是一種受人體免疫系統啟發的計算方法，它模仿了生物免疫系統中識別和消滅外來病原體的過程。這種演算法通常用於解決最佳化問題，尤其是在那些傳統演算法難以解決的複雜和高維度問題上。以下是免疫最佳化演算法核心理念的探險式解釋：

免疫系統的啟發：人體的免疫系統能夠識別並消滅入侵的病原體，如細菌和病毒。這種能力是透過識別特定的抗原（病原體的特定部分）來實現的。
抗體的多樣性：免疫系統產生大量不同的抗體，以覆蓋廣泛的抗原型別。在最佳化演算法中，這相當於生成多種可能的解決方案。
克隆選擇：當抗體成功識別並結合到抗原上時，免疫系統會克隆這些抗體，以增加其數量，從而更有效地對抗病原體。在演算法中，這相當於對好的解決方案進行復制和增強。
記憶細胞：免疫系統會保留一些成功的抗體作為記憶細胞，以便在將來遇到相同或相似的病原體時快速響應。在演算法中，這意味著保留並更新最佳解決方案。
負選擇：免疫系統透過負選擇過程去除那些可能攻擊自身組織的抗體。在演算法中，這相當於避免找到的解決方案與問題約束條件相沖突。
親和力成熟：抗體透過突變和選擇過程提高其與抗原的結合能力。在演算法中，這相當於透過區域性搜尋來改進解決方案的質量。
多樣性維持：免疫系統透過引入新的抗體來保持其多樣性，以應對不斷變化的病原體。在演算法中，這涉及到引入新的解決方案，以避免陷入區域性最優。
並行處理：免疫系統同時處理多個抗原，這在演算法中相當於並行評估多個解決方案。

上面的比喻可能不完全準確，免疫最佳化演算法演算法詳細理解見：萬字長文了解免疫演算法原理。在scikit-opt庫中的免疫最佳化演算法主要用於求解旅行商問題。提供IA_TSP類用於實現免疫最佳化演算法，IA_TSP類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式，用於評估解的適應度或質量
n_dim	-	城市個數，表示問題規模，即最佳化問題中的變數數量
size_pop	50	種群規模，表示每次迭代中個體的數量
max_iter	200	最大迭代次數，演算法執行的迭代次數上限
prob_mut	0.001	變異機率，表示個體在迭代過程中發生變異的可能性
T	0.7	抗體與抗體之間的親和度閾值，用於確定個體間的相似性，大於這個閾值認為個體間親和，否則認為不親和
alpha	0.95	多樣性評價指數，用於平衡抗體的適應度和多樣性，即抗體和抗原的重要性與抗體濃度的重要性

以下案例展示瞭如何利用IA_TSP求解旅行商問題：

import numpy as np

# 從 sko 庫的 demo_func 模組中匯入用於生成TSP問題的函式
from sko.demo_func import function_for_TSP

# num_points 是點的數量，points_coordinate 是點的座標，distance_matrix 是點之間的距離矩陣
# cal_total_distance 是一個函式，用於計算給定路徑的總距離
num_points, points_coordinate, distance_matrix, cal_total_distance = function_for_TSP(num_points=8)

from sko.IA import IA_TSP
# 定義最佳化函式
ia_tsp = IA_TSP(func=cal_total_distance, n_dim=num_points, size_pop=500, max_iter=800, prob_mut=0.2,
                T=0.7, alpha=0.95)

# best_points 是最優路徑的點的索引序列
# best_distance 是最優路徑的總距離
# 請注意，該段程式碼執行速度可能較慢，因為它涉及到大量的迭代計算
best_points, best_distance = ia_tsp.run()
# ia.generation_best_Y 每一代的最優函式值
# ia.generation_best_X 每一代的最優函式值對應的輸入值
# ia.all_history_FitV 每一代的每個個體的適應度
# ia.all_history_Y 每一代每個個體的函式值
# ia.best_y 最優函式值
# ia.best_x 最優函式值對應的輸入值
print('best routine:', best_points, 'best_distance:', best_distance)


# 匯入 matplotlib 庫的 pyplot 模組，用於繪圖
import matplotlib.pyplot as plt

# 建立一個圖形和座標軸物件
fig, ax = plt.subplots(1, 1)

# 將最優路徑的點連同起點一起構成閉環路徑
best_points_ = np.concatenate([best_points, [best_points[0]]])

# 根據最優路徑的點索引獲取對應的座標
best_points_coordinate = points_coordinate[best_points_, :]

# 繪製最優路徑的點，並用紅色線條連線
ax.plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1], 'o-r')

# 在最優路徑上新增箭頭，表示行進方向
for i in range(num_points):
    # 使用 annotate 函式在兩點之間新增箭頭
    ax.annotate('', xy=best_points_coordinate[i + 1], xytext=best_points_coordinate[i],
                   arrowprops=dict(arrowstyle='->', color='green'))

# 用綠色方塊標記最優路徑的起點
ax.plot(best_points_coordinate[0, 0], best_points_coordinate[0, 1], 'gs')

# 顯示圖形
plt.show()

best routine: [0 3 4 6 1 5 7 2] best_distance: [2.49060449]

png

1.7 人工魚群演算法

人工魚群演算法(Artificial Fish-Swarm Algorithm, AFSA)是一種模擬自然界中魚群行為的最佳化演算法。它的核心理念可以從以下幾個方面來解釋，就像探險一樣：

探索與開發：探險家需要在未知的領域進行探索，同時也需要開發已知的資源。類似地，AFSA中的個體（模擬魚）在搜尋空間中進行探索，尋找最優解，同時也對當前找到的較好解進行開發，以期進一步最佳化。
群體行為：探險隊中的成員通常會相互協作，共享資訊。魚群演算法中，魚群透過群體行為，如追隨、避障等，來共同尋找食物資源。在演算法中，個體會根據群體中其他個體的位置和資訊來調整自己的搜尋策略。
隨機性與適應性：探險過程中充滿了不確定性，探險家需要根據環境的變化靈活調整策略。人工魚群演算法中的個體也會根據環境的反饋進行自適應調整，比如改變搜尋方向或速度。
生存競爭：在探險中，資源有限，探險家之間可能存在競爭。魚群演算法中，個體之間也存在競爭，透過競爭可以促進演算法的多樣性，避免早熟收斂。
學習和記憶：探險家會根據以往的經驗來指導未來的行動。在AFSA中，個體會學習並記憶先前的成功經驗，利用這些資訊來指導當前的搜尋過程。
環境互動：探險家需要與環境互動，瞭解環境特徵。魚群演算法中的個體也會與搜尋空間的環境進行互動，根據環境反饋調整搜尋策略。
目標導向：探險的最終目的是達到某個特定的目標。AFSA的目標是找到問題的最優解，所有的行為和策略都是為了這一目標服務。

上面的比喻可能不完全準確，人工魚群演算法演算法詳細理解見：人工魚群演算法超詳細解析。scikit-opt庫中的AFSA類用於實現人工魚群演算法，AFSA類的構造引數介紹如下：

引數	預設值	意義
func	-	目標函式
n_dim	-	目標函式的維度
size_pop	50	種群規模
max_iter	300	最大迭代次數
max_try_num	100	最大嘗試捕食次數
step	0.5	每一步的最大位移比例
visual	0.3	定義了個體在搜尋過程中能夠感知到其他個體的範圍
q	0.98	魚的感知範圍衰減係數
delta	0.5	擁擠度閾值，較大值可能導致過度聚集和區域性搜尋

以下程式碼展示瞭如何將人工魚群演算法用於多元函式最佳化：

def func(x):
    x1, x2, x3 = x
    return (x1 - x2) ** 2 + (x2 - 0.01) ** 2 + x3 ** 2

from sko.AFSA import AFSA

afsa = AFSA(func, n_dim=3, size_pop=50, max_iter=100,
            max_try_num=100, step=0.5, visual=0.3,
            q=0.98, delta=0.5)
best_x, best_y = afsa.run()
print(best_x, best_y)

[ 0.00593808  0.00693321 -0.00039819] 1.0554063699211183e-05

2 參考

scikit-opt
scikit-opt-doc
差分進化演算法的基本概念和原理
10分鐘搞懂遺傳演算法
粒子群最佳化演算法的詳細解讀
模擬退火演算法詳解
從模擬退火到退火進化演算法
一文搞懂什麼是蟻群最佳化演算法
萬字長文了解免疫演算法原理
人工魚群演算法超詳細解析

【scikit-opt】七大啟發式演算法的使用
2024-08-25
演算法
[python] Python平行計算庫Joblib使用指北
2024-08-10
Python
[python] Python日誌記錄庫loguru使用指北
2024-06-30
Python
元啟發式演算法庫 MEALPY 初體驗-遺傳演算法為例
2024-04-11
演算法
[機器學習] 低程式碼機器學習工具PyCaret庫使用指北
2024-06-01
機器學習
Protocol Buffer 使用指北
2019-02-24
Protocol
UIAppearcance 使用指北
2018-10-27
UIAPP
proxychains 使用指北
2018-04-17
AI
前端開發 ngrok 指北
2019-04-21
前端
curl 使用指北&實戰
2019-09-21
[深度學習] 時間序列分析工具TSLiB庫使用指北
2024-08-28
深度學習
laravel-echo前端使用指北
2021-02-14
Laravel前端
tar命令基本、進階使用指北
2020-05-03
Android-ViewModel 使用指北
2018-12-03
AndroidView
Performance API不完全使用指北
2023-02-06
ORMAPI
啟發式合併
2024-10-08
Flutter 外掛 webview_flutter 使用指北
2019-09-12
FlutterWebView
Python 簡單入門指北（試讀版）
2017-11-03
Python
運籌優化（十四）--離散優化的啟發式演算法
2019-01-19
優化演算法
js版九宮格拼圖與啟發式搜尋（A*演算法）
2018-03-30
JS演算法
RUI：受 SwiftUI啟發的實驗性宣告式Rust UI庫
2022-03-02
UISwiftRust
[深度學習] 計算機視覺低程式碼工具Supervision庫使用指北
2024-03-18
深度學習計算機視覺
msf stagers開發不完全指北(二)
2020-06-28
Python演算法庫的安裝
2018-11-08
Python演算法
Android啟發式尋路
2021-09-09
Android
樹上啟發式合併
2024-05-02
python-函式使用
2020-12-25
Python函式
『併發包入坑指北』之阻塞佇列
2019-04-10
佇列
【機器學習演算法實現】logistic迴歸__基於Python和Numpy函式庫
2016-05-10
機器學習演算法Python函式
React Hooks 指北
2022-02-21
ReactHook
Python curses庫如何使用
2021-09-11
Python
21.8 Python 使用BeautifulSoup庫
2023-10-27
Python
【Python3網路爬蟲開發實戰】3-基本庫的使用-3正規表示式
2018-03-19
Python爬蟲
TSP問題（換位表達，啟發式交叉，啟發式變異）C++實現
2020-11-20
C++
輕鬆開發React元件和造輪子指北
2019-03-04
React元件
python函式的基本使用
2021-03-04
Python函式
寫給自己看的在 Vue 下使用 Typescript 指北
2019-02-26
VueTypeScript
python演算法常用技巧與內建庫
2020-10-17
Python演算法

[python] 啟發式演算法庫scikit-opt使用指北

0 背景介紹

1 演算法使用

1.1 差分進化演算法

1.2 遺傳演算法

1.2.1 基礎遺傳演算法

1.2.2 用於解決旅行商問題的遺傳演算法

1.3 粒子群演算法

1.4 模擬退火演算法

1.4.1 基礎模擬退火演算法

1.4.2 用於解決旅行商問題的模擬退火演算法

1.5 蟻群演算法

1.6 免疫最佳化演算法

1.7 人工魚群演算法

2 參考

相關文章