1. 程式人生 > >Python遺傳和進化演算法框架(一)Geatpy快速入門

Python遺傳和進化演算法框架(一)Geatpy快速入門

Geatpy是一個高效能的Python遺傳演算法庫以及開放式進化演算法框架,由華南理工大學、華南農業大學、德州奧斯汀公立大學的學生聯合團隊開發。

它提供了許多已實現的遺傳和進化演算法相關運算元的庫函式,如初始化種群、選擇、交叉、變異、重插入、多種群遷移、多目標優化非支配排序等,並且提供開放式的進化演算法框架來實現多樣化的進化演算法。其執行效率高於Matlab遺傳演算法工具箱和Matlab第三方遺傳演算法工具箱Gaot、gatbx、GEATbx,學習成本低。

Geatpy支援二進位制/格雷碼編碼種群、實數值種群、整數值種群、排列編碼種群。支援輪盤賭選擇、隨機抽樣選擇、錦標賽選擇。提供單點交叉、兩點交叉、洗牌交叉、部分匹配交叉(PMX)、線性重組、離散重組、中間重組等重組運算元。提供簡單離散變異、實數值變異、整數值變異、互換變異等變異運算元。支援隨機重插入、精英重插入。支援awGA、rwGA、i-awGA、nsga2、快速非支配排序等多目標優化的庫函式、提供進化演算法框架下的常用進化演算法模板等。

關於遺傳演算法、進化演算法的學習資料,在官網中https://www.geatpy.com (repairing)有詳細講解以及相關的學術論文連結。同時網上也有很多資料。

閒話少說……下面分享一下近幾個星期以來的學習心得:

先說一下安裝方法

首先是要windows系統,Python要是3.5,3.6或3.7版本 ,並且安裝了pip。只需在控制檯執行

pip install geatpy

即可安裝成功。或者到github上下載:https://github.com/geatpy-dev/geatpy 推薦是直接用pip的方式安裝。因為這樣有利於後續的更新。我為了方便執行demo程式碼以及檢視原始碼,因此另外在github上也下載了(但仍用pip方式安裝)。

有些初學Python的讀者反映還是不知道怎麼安裝,或者安裝之後不知道怎麼寫程式碼。這裡推薦安裝Anaconda,它集成了Python的許多常用的執行庫,比如Numpy、Scipy等。其內建的Spyder風格跟Matlab類似,給人熟悉的感覺,更容易上手。

再說一下更新方法

(注意:專案改組後,Geatpy已於2018.09.20更新至1.0.6版本,新版本解決了一些版權問題。請及時更新。)

Geatpy在持續更新。可以通過以下命令使電腦上的版本與官方最新版保持一致:

pip install --upgrade geatpy

Geatpy提供2種方式實現遺傳演算法。先來講一下第一種最基本的實現方式:編寫程式設計指令碼。

1. 編寫指令碼實現遺傳演算法:

用過謝菲爾德大學的Matlab遺傳演算法庫Gatbx以及其升級版——GEATbx的朋友應該非常熟悉下面的Matlab指令碼:

%% matlab_gatbx_test.m

%遺傳演算法求f(x)=x*sin(10*pi*x)+2.0,在[-1,2]上的最大值
figure(1);
fplot(@(variable)variable.*sin(10*pi*variable)+2.0,[-1,2]);   %畫出函式曲線
tic %開始計時
%定義遺傳演算法引數
NIND=40;        %個體數目(Number of individuals)
MAXGEN=25;      %最大遺傳代數(Maximum number of generations)
PRECI=20;       %變數的二進位制位數(Precision of variables)
GGAP=0.9;       %代溝(Generation gap)說明子代與父代的重複率為0.1
trace=zeros(MAXGEN,2);                        %尋優結果的初始值
FieldD=[20;-1;2;1;0;1;1];                      %區域描述器(Build field descriptor),第2、3行為自變數的下界和上界
Chrom=crtbp(NIND, PRECI);                      %定義初始種群
gen=0;                                         %代計數器
variable=bs2rv(Chrom, FieldD);                 %計算初始種群的十進位制轉換
ObjV=shang(variable);        %計算目標函式值
while gen<MAXGEN %進化MAXGEN代
   FitnV=ranking(-ObjV);                                  %分配適應度值(Assign fitness values)ranking函式的功能就是目標值越小的分配值越大,
                                                          %本例求解最大值,應該要是他的適應度值更大,故必須使得ranking內的數越小,這樣分配的適應度值才能大
   SelCh=select('sus', Chrom, FitnV, GGAP);               %選擇,使用sus方式,也可以改用rws方式
   SelCh=recombin('xovsp', SelCh, 0.7);                   %重組,選用xovsp方式
   SelCh=mut(SelCh);                                      %變異
   variable=bs2rv(SelCh, FieldD);                         %子代個體的十進位制轉換,把染色體變為十進位制
   ObjVSel=shang(variable);             %計運算元代的目標函式值
   [Chrom ObjV]=reins(Chrom, SelCh, 1, 1, ObjV, ObjVSel); %重插入子代的新種群
   variable=bs2rv(Chrom, FieldD);                         %子代個體的十進位制轉換,轉為十進位制
   gen=gen+1;                                             %代計數器增加
   %輸出最優解及其序號,並在目標函式影象中標出,Y為最優解,I為種群的序號
   [Y, I]=max(ObjV);hold on;
   plot(variable(I), Y, 'ro');
   trace(gen,1)=max(ObjV);                                %遺傳演算法效能跟蹤,把當代的最優值放入trace矩陣的第一行第目前代數列
   trace(gen,2)=sum(ObjV)/length(ObjV);                   %把當代種群目標函式的均值,放入trace矩陣的第二行第目前代數列
end
variable=bs2rv(Chrom, FieldD);                            %最優個體的十進位制轉換,轉回十進位制,以便輸出
result = max(trace(:,1)) %輸出搜尋到的目標函式最大值
toc %結束計時
hold on;
grid on;
plot(variable,ObjV,'b*');
figure(2);
plot(trace(:,1)); %把trace矩陣的第一列畫出來(記錄的是每一代的最優值)
hold on;
plot(trace(:,2),'-.');grid %把trace矩陣的第2列畫出來(記錄的是每一代種群目標函式均值)
legend('解的變化','種群均值的變化')


function z=shang(x) % 目標函式的核心部分(即缺省了優化目標的純函式)
    z=x.*sin(10*pi*x)+2.0;
end

這是在Matlab的gatbx工具箱下實現簡單遺傳演算法搜尋f(x)=x*sin(10*pi*x)+2.0,在[-1,2]上的最大值的指令碼程式,執行結果如下:

result =

    3.8500

時間已過 0.081932 秒。

再看一下在Geatpy上如何編寫指令碼:

"""demo.py"""
import numpy as np
import geatpy as ga # 匯入geatpy庫
import matplotlib.pyplot as plt
import time

"""============================目標函式============================"""
def aim(x):                    # 傳入種群染色體矩陣解碼後的基因表現型矩陣
    return x * np.sin(10 * np.pi * x) + 2.0
x = np.linspace(-1, 2, 200)
plt.plot(x, aim(x)) # 繪製目標函式影象
start_time = time.time()       # 開始計時
"""============================變數設定============================"""
x1 = [-1, 2]                   # 自變數範圍
b1 = [1, 1]                    # 自變數邊界
codes = [1]                    # 變數的編碼方式,2個變數均使用格雷編碼
precisions =[5]                # 變數的精度
scales = [0]                   # 採用算術刻度
ranges=np.vstack([x1]).T       # 生成自變數的範圍矩陣
borders=np.vstack([b1]).T      # 生成自變數的邊界矩陣
"""========================遺傳演算法引數設定========================="""
NIND = 40;                     # 種群個體數目
MAXGEN = 25;                   # 最大遺傳代數
GGAP = 0.9;                    # 代溝:說明子代與父代的重複率為0.1                                      
"""=========================開始遺傳演算法進化========================"""
FieldD = ga.crtfld(ranges,borders,precisions,codes,scales) # 呼叫函式建立區域描述器
Lind = np.sum(FieldD[0, :])    # 計算編碼後的染色體長度
Chrom = ga.crtbp(NIND, Lind)   # 根據區域描述器生成二進位制種群
variable = ga.bs2rv(Chrom, FieldD) #對初始種群進行解碼
ObjV = aim(variable)           # 計算初始種群個體的目標函式值
pop_trace = (np.zeros((MAXGEN, 2)) * np.nan) # 定義進化記錄器,初始值為nan
ind_trace = (np.zeros((MAXGEN, Lind)) * np.nan) # 定義種群最優個體記錄器,記錄每一代最優個體的染色體,初始值為nan
# 開始進化!!
for gen in range(MAXGEN):
    FitnV = ga.ranking(-ObjV)  # 根據目標函式大小分配適應度值(由於遵循目標最小化約定,因此最大化問題要對目標函式值乘上-1)
    SelCh=ga.selecting('sus', Chrom, FitnV, GGAP) # 選擇
    SelCh=ga.recombin('xovsp', SelCh, 0.7) # 重組(採用單點交叉方式,交叉概率為0.7)
    SelCh=ga.mutbin(SelCh) # 二進位制種群變異
    variable = ga.bs2rv(SelCh, FieldD) # 對育種種群進行解碼(二進位制轉十進位制)
    ObjVSel = aim(variable) # 求育種個體的目標函式值
    [Chrom,ObjV] = ga.reins(Chrom,SelCh,SUBPOP,2,1,ObjV,ObjVSel) # 重插入得到新一代種群
    # 記錄
    best_ind = np.argmax(ObjV) # 計算當代最優個體的序號
    pop_trace[gen, 0] = ObjV[best_ind] # 記錄當代種群最優個體目標函式值
    pop_trace[gen, 1] = np.sum(ObjV) / ObjV.shape[0] # 記錄當代種群的目標函式均值
    ind_trace[gen, :] = Chrom[best_ind, :] # 記錄當代種群最優個體的變數值
# 進化完成
end_time = time.time() # 結束計時
"""============================輸出結果及繪圖================================"""
print('目標函式最大值:',np.max(pop_trace[:, 0])) # 輸出目標函式最大值
variable = ga.bs2rv(ind_trace, FieldD) # 解碼得到表現型
print('用時:', end_time - start_time)
plt.plot(variable, aim(variable),'bo')

執行結果如下:

目標函式最大值: 3.850273756279405 用時: 0.0509946346282959

對比發現,Geatpy的執行效率要高於Matlab,而且結果較好。

對比Matlab程式碼和Python程式碼,我們會發現Geatpy提供風格極為相似的庫函式,有Matlab相關程式設計經驗的基本上可以無縫轉移到Python上利用Geatpy進行遺傳演算法程式開發。

Geatpy提供了詳盡的API文件,比如要檢視上面程式碼中的"ranking"函式是幹什麼的,可以在python中執行

import geatpy as ga
help(ga.ranking)

即可看到"ranking"函式的相關使用方法。另外也可以參見github上面的文件:

另外官網上也有更多詳盡的Geatpy教程,Geatpy官網http:www.geatpy.com將於九月初重新上線。

2. 利用框架實現遺傳演算法。

Geatpy提供開放的進化演算法框架。即“函式介面”+“進化演算法模板”。對於一些複雜的進化演算法,如多目標進化優化、改進的遺傳演算法等,利用上面所說的編寫指令碼是非常麻煩的,改用框架的方法可以極大提高程式設計效率

這裡給出一個利用框架實現NSGA-II演算法求多目標優化函式ZDT-1的帕累託前沿面的例子:

首先編寫函式介面檔案:

""" aimfuc.py """

# ZDT1
def aimfuc(Chrom):
    ObjV1 = Chrom[:, 0]
    gx = 1 + (9 / 29) * np.sum(Chrom[:, 1:30], 1)
    hx = 1 - np.sqrt(ObjV1 / gx)
    ObjV2 = gx * hx
    return np.array([ObjV1, ObjV2]).T

然後編寫指令碼,使用Geatpy提供的nsga2演算法的進化演算法模板(nsga2_templet):

"""main.py"""
import numpy as np
import geatpy as ga # 匯入geatpy庫

# 獲取函式介面地址
AIM_M = __import__('aimfuc')
AIM_F = 'aimfuc'

"""============================變數設定============================"""
ranges = np.vstack([np.zeros((1,30)), np.ones((1,30))])   # 生成自變數的範圍矩陣
borders = np.vstack([np.ones((1,30)), np.ones((1,30))])   # 生成自變數的邊界矩陣
precisions = [4] * 30              # 自變數的編碼精度
"""========================遺傳演算法引數設定========================="""
NIND = 25                 # 種群規模
MAXGEN = 1000             # 最大遺傳代數
GGAP = 1;                 # 代溝:子代與父代的重複率為(1-GGAP),由於後面使用NSGA2演算法,因此該引數無用
selectStyle = 'tour'      # 遺傳演算法的選擇方式
recombinStyle = 'xovdprs' # 遺傳演算法的重組方式,設為兩點交叉
recopt = 0.9              # 交叉概率
pm = 0.1                  # 變異概率
SUBPOP = 1                # 設定種群數為1f
maxormin = 1              # 設定標記表明這是最小化目標
MAXSIZE = 1000            # 帕累托最優集最大個數
FieldDR = ga.crtfld(ranges, borders, precisions) # 生成區域描述器
"""=======================呼叫進化演算法模板進行種群進化==================="""
# 得到帕累托最優解集NDSet以及解集對應的目標函式值NDSetObjV
[ObjV, NDSet, NDSetObjV, times] = ga.nsga2_templet(AIM_M, AIM_F, None, None, FieldDR, 'R', maxormin, MAXGEN, MAXSIZE, NIND, SUBPOP, GGAP, selectStyle, recombinStyle, recopt, pm, distribute = True, drawing = 1)

執行結果如下:

用時: 7.359716176986694 秒 帕累託前沿點個數: 479 個 單位時間找到帕累託前沿點個數: 65 個

可以改用q_sorted_templet快速非支配排序的多目標優化進化演算法模板,可以得到更好的效率和更好的結果:

進化演算法的核心演算法邏輯是寫在進化演算法模板內部的,程式碼是開源的,我們可以參考Geatpy進化演算法模板的原始碼來自定義演算法模板,以實現豐富多樣的進化演算法,如差分進化演算法、改進的遺傳演算法等:

後面的部落格將深入理解Geatpy的庫函式用法,以及探討框架的核心——進化演算法模板的實現。還會講一些使用Geatpy解決問題的案例。歡迎繼續跟進~感謝!

下一篇:Python遺傳和進化演算法框架(二)Geatpy庫函式和資料結構