參考文檔連結：scikit-opt
本章繼續Python的優化算法系列。
優化算法，尤其是啟發式的仿生智能算法在最近很火，它适用于解決管理學，運籌學，統計學裡面的一些優化問題。比如線性規劃，整數規劃，動态規劃，非線性限制規劃，甚至是超參數搜尋等等方向的問題。
但是一般的優化算法還是matlab裡面用的多，Python相關代碼較少。部落客在參考了很多文章的代碼和子產品之後，決定學習 scikit-opt 這個子產品。這個優化算法子產品對新手很友好，代碼簡潔，上手簡單。而且代碼和官方文檔是中國人寫的，還有很多案例，學起來就沒什麼壓力...
缺點是包裝的算法種類目前還不算多，隻有七種：（差分進化算法、遺傳算法、粒子群算法、模拟退火算法、蟻群算法、魚群算法、免疫優化算法） /（其實已經夠用了）
本次帶來的是 數學模組化裡面經常使用的遺傳算法的使用示範。數學原理就不多說了
首先安裝子產品，在cmd裡面或者anaconda prompt裡面輸入：
pip install scikit-opt
這個包很小，很快就能裝好。
遺傳算法
首先定義我們要解決的優化的問題：
該函數有大量的局部最小值，具有很強的沖擊力，在(0,0) 處的全局最小值，值為 0
代碼
import numpy as np
def schaffer(p):
x1, x2 = p
x = np.square(x1) + np.square(x2)
return 0.5 + (np.square(np.sin(x)) - 0.5) / np.square(1 + 0.001 * x)
調用遺傳算法GA解決：
from sko.GA import GA
ga = GA(func=schaffer, n_dim=2, size_pop=50, max_iter=800, prob_mut=0.001, lb=[-1, -1], ub=[1, 1], precision=1e-7)
best_x, best_y = ga.run()
print('best_x:', best_x, '\n', 'best_y:', best_y)
可以看到基本找到了全局最小值和對應的x 。
畫出疊代次數的圖
import pandas as pd
import matplotlib.pyplot as plt
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()
參數詳解
上面用到了很多參數，不懂可以看看這裡
輸入參數
輸出參數
遺傳算法進行整數規劃
在多元優化時，想讓哪個變量限制為整數，就設定 precision 為 整數 即可。
例如，我想讓我的自定義函數 demo_func 的某些變量限制為整數+浮點數（分别是整數，整數，浮點數），那麼就設定 precision=[1, 1, 1e-7]
例子如下：
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,1,1e-7])
best_x, best_y = ga.run()
print('best_x:', best_x, '\n', 'best_y:', best_y)
可以看到第一個、第二個變量都是整數，第三個就是浮點數了 。
遺傳算法用于旅行商問題
商旅問題(TSP)就是路徑規劃的問題，比如有很多城市，你都要跑一遍，那麼先去哪個城市再去哪個城市可以讓你的總路程最小。
實際問題需要一個城市坐标資料，比如你的出發點位置為（0,0），第一個城市離位置為（x1,y1）,第二個為（x2,y2）.........這裡沒有實際資料，就直接随機生成了。
import numpy as np
from scipy import spatial
import matplotlib.pyplot as plt
num_points = 50
points_coordinate = np.random.rand(num_points, 2) # generate coordinate of points
points_coordinate
這裡定義的是50個城市，每個城市的坐标都在是上圖随機生成的矩陣。
然後我們把它變成類似相關系數裡面的矩陣
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric='euclidean')
distance_matrix.shape
這個矩陣就能得出每個城市之間的距離，算上自己和自己的距離（0），總共有2500個數。
定義問題
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)])
求解
from sko.GA import GA_TSP
ga_tsp = GA_TSP(func=cal_total_distance, n_dim=num_points, size_pop=50, max_iter=500, prob_mut=1)
best_points, best_distance = ga_tsp.run()
best_distance
得到的最小距離
畫圖檢視計算出來的路徑，還有疊代次數和y的關系
fig, ax = plt.subplots(1, 2,figsize=(10,4))
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')
ax[1].plot(ga_tsp.generation_best_Y)
plt.show()
使用遺傳算法進行曲線拟合
建構資料集
import numpy as np
import matplotlib.pyplot as plt
from sko.GA import GA
x_true = np.linspace(-1.2, 1.2, 30)
y_true = x_true ** 3 - x_true + 0.4 * np.random.rand(30)
plt.plot(x_true, y_true, 'o')
建構的資料是y=x^3-x+0.4,然後加上了随機擾動項。如圖
定義需要拟合的函數（三次函數），然後将殘差作為目标函數去求解
def f_fun(x, a, b, c, d):
return a * x ** 3 + b * x ** 2 + c * x + d #三次函數
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=0.919x^3-0.0244x^2-0.939x+0.2152.
拟合效果還行
畫出拟合線
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()