數學建模十大演算法之——遺傳演算法及Python實現

來自專欄 Eureka機器學習讀書筆記6 人贊了文章

專欄文章匯總

文章結構如下:

1: 遺傳演算法理論的由來

2: 遺傳演算法定義

3: 遺傳演算法具體步驟

3.1 初始化(編碼)

3.2 編碼補充

3.3 適應度函數

3.3 選擇

3.4 交叉

3.5 變異

4: 遺傳演算法應用

4.1 問題

4.2 創造染色體(編碼)

4.3 個體、種群與進化

4.4 python實現

5: 參考文獻


註:本文是根據參考文獻幾篇文章及個人理解整合而來


1: 遺傳演算法理論的由來

我們先從查爾斯·達爾文的一句名言開始:

能夠生存下來的往往不是最強大的物種,也不是最聰明的物種,而是最能適應環境的物種。

你也許在想:這句話和遺傳演算法有什麼關係?其實遺傳演算法的整個概念就基於這句話。

讓我們用一個基本例子來解釋 :

我們先假設一個情景,現在你是一國之王,為了讓你的國家免於災禍,你實施了一套法案:

  • 你選出所有的好人,要求其通過生育來擴大國民數量。
  • 這個過程持續進行了幾代。
  • 你將發現,你已經有了一整群的好人。

這個例子雖然不太可能,但是我用它是想幫助你理解概念。也就是說,我們改變了輸入值(比如:人口),就可以獲得更好的輸出值(比如:更好的國家)。現在,我假定你已經對這個概念有了大致理解,認為遺傳演算法的含義應該和生物學有關係。那麼我們就快速地看一些小概念,這樣便可以將其聯繫起來理解。


2: 遺傳演算法定義

首先我們回到前面討論的那個例子,並總結一下我們做過的事情。

  1. 首先,我們設定好了國民的初始人群大小。
  2. 然後,我們定義了一個函數,用它來區分好人和壞人。
  3. 再次,我們選擇出好人,並讓他們繁殖自己的後代。
  4. 最後,這些後代們從原來的國民中替代了部分壞人,並不斷重複這一過程。

遺傳演算法實際上就是這樣工作的,也就是說,它基本上儘力地在某種程度上模擬進化的過程。

因此,為了形式化定義一個遺傳演算法,我們可以將它看作一個優化方法,它可以嘗試找出某些輸入,憑藉這些輸入我們便可以得到最佳的輸出值或者是結果。遺傳演算法的工作方式也源自於生物學,具體流程見下圖:

遺傳演算法的基本特徵:

  1. 智能式搜索:依據適應度(目標函數)進行只能搜索
  2. 漸進式優化:利用複製、交換、突變等操作,使下一代結果優於上一代
  3. 全局最優解:採用交換和突變操作產生新個體,使得搜索得到的優化結果逼近全局最優解
  4. 黑箱式結構:根據問題的特性進行編碼(輸入)和確定適應度(輸出),具有隻考慮輸入與輸出關係的黑箱式就夠,並不深究輸入與輸出關係的原因
  5. 通用性強:不要求明確的數學表達式,只需要一些簡單的原則要求,可應用於解決離散問題、函數關係不明確的複雜問題
  6. 並行式運算:每次迭代計算都是對群體中所有個體同時進行運算,是並行式運算方式,搜索速度快

3: 遺傳演算法具體步驟

為了讓講解更為簡便,我們先來理解一下著名的組合優化問題「背包問題」。

比如,你準備要去野遊 1 個月,但是你只能背一個限重 30 公斤的背包。現在你有不同的必需物品,它們每一個都有自己的「生存點數」(具體在下表中已給出)。因此,你的目標是在有限的背包重量下,最大化你的「生存點數」。

3.1 初始化(編碼)

這裡我們用遺傳演算法來解決這個背包問題。第一步是定義我們的總體。總體中包含了個體,每個個體都有一套自己的染色體。

我們知道,染色體可表達為二進位數串,在這個問題中,1 代表接下來位置的基因存在,0 意味著丟失。(譯者註:作者這裡借用染色體、基因來解決前面的背包問題,所以特定位置上的基因代表了上方背包問題表格中的物品,比如第一個位置上是 Sleeping Bag,那麼此時反映在染色體的『基因』位置就是該染色體的第一個『基因』。)

現在,我們將圖中的 4 條染色體看作我們的總體初始值。

3.2 編碼補充

二進位編碼

二進位編碼由二進位符號0和1所組成的二值符號集。它有以下一些優點: 1. 編碼、解碼操作簡單易行 2. 交叉、變異等遺傳操作便於實現 3. 合最小字符集編碼原則 4. 利用模式定理對演算法進行理論分析。

二進位編碼的缺點是:對於一些連續函數的優化問題,由於其隨機性使得其局部搜索能力較差,如對於一些高精度的問題,當解迫近於最優解後,由於其變異後表現型變化很大,不連續,所以會遠離最優解,達不到穩定。

格雷碼

格雷碼編碼是其連續的兩個整數所對應的編碼之間只有一個碼位是不同的,其餘碼位完全相同。

二進位碼轉為格雷碼:異或運算:同則為0,異則為1。公式如下:

g_m = b_m, g_i = b_{i+1}igoplus b_i ,i=m-1,m-2,cdots 2,1

由格雷碼轉二進位的轉換公式為:

b_m = g_m, b_i = b_{i+1}igoplus g_i, i=m-1,m-2,cdots 2,1

優點:增強遺傳演算法的局部搜索能力,便於對連續函數進行局部空間搜索。使用非常廣泛。

浮點編碼法

二進位編碼雖然簡單直觀,但明顯地。但是存在著連續函數離散化時的映射誤差。個體長度較短時,可能達不到精度要求,而個體編碼長度較長時,雖然能提高精度,但增加了解碼的難度,使遺傳演算法的搜索空間急劇擴大。

所謂浮點法,是指個體的每個基因值用某一範圍內的一個浮點數來表示。編碼長度等於決策變數的個數。 在浮點數編碼方法中,必須保證基因值在給定的區間限制範圍內,遺傳演算法中所使用的交叉、變異等遺傳運算元也必須保證其運算結果所產生的新個體的基因值也在這個區間限制範圍內。

優點:

  1. 適用於在遺傳演算法中表示範圍較大的數
  2. 適用於精度要求較高的遺傳演算法
  3. 便於較大空間的遺傳搜索
  4. 改善了遺傳演算法的計算複雜性,提高了運算交率
  5. 便於遺傳演算法與經典優化方法的混合使用
  6. 便於設計針對問題的專門知識的知識型遺傳運算元
  7. 便於處理複雜的決策變數約束條件

符號編碼法

符號編碼法是指個體染色體編碼串中的基因值取自一個無數值含義、而只有代碼含義的符號集如{A,B,C…}。

優點:

  1. 符合有意義積術塊編碼原則
  2. 便於在遺傳演算法中利用所求解問題的專門知識
  3. 便於遺傳演算法與相關近似演算法之間的混合使用。

3.3 適應度函數

接下來,讓我們來計算一下前兩條染色體的適應度分數。對於 A1 染色體 [100110] 而言,有:

類似地,對於 A2 染色體 [001110] 來說,有:

對於這個問題,我們認為,當染色體包含更多生存分數時,也就意味著它的適應性更強。因此,由圖可知,染色體 1 適應性強於染色體 2。

3.4 選擇

現在,我們可以開始從總體中選擇適合的染色體,來讓它們互相『交配』,產生自己的下一代了。這個是進行選擇操作的大致想法,但是這樣將會導致染色體在幾代之後相互差異減小,失去了多樣性。因此,我們一般會進行「輪盤賭選擇法」(Roulette Wheel Selection method)。

想像有一個輪盤,現在我們將它分割成 m 個部分,這裡的 m 代表我們總體中染色體的個數。每條染色體在輪盤上佔有的區域面積將根據適應度分數成比例表達出來。

基於上圖中的值,我們建立如下「輪盤」。

現在,這個輪盤開始旋轉,我們將被圖中固定的指針(fixed point)指到的那片區域選為第一個親本。然後,對於第二個親本,我們進行同樣的操作。有時候我們也會在途中標註兩個固定指針,如下圖:

通過這種方法,我們可以在一輪中就獲得兩個親本。我們將這種方法成為「隨機普遍選擇法」(Stochastic Universal Selection method)。

3.5 交叉

在上一個步驟中,我們已經選擇出了可以產生後代的親本染色體。那麼用生物學的話說,所謂「交叉」,其實就是指的繁殖。現在我們來對染色體 1 和 4(在上一個步驟中選出來的)進行「交叉」,見下圖:

這是交叉最基本的形式,我們稱其為「單點交叉」。這裡我們隨機選擇一個交叉點,然後,將交叉點前後的染色體部分進行染色體間的交叉對調,於是就產生了新的後代。

如果你設置兩個交叉點,那麼這種方法被成為「多點交叉」,見下圖:

3.5 變異

如果現在我們從生物學的角度來看這個問題,那麼請問:由上述過程產生的後代是否有和其父母一樣的性狀呢?答案是否。在後代的生長過程中,它們體內的基因會發生一些變化,使得它們與父母不同。這個過程我們稱為「變異」,它可以被定義為染色體上發生的隨機變化,正是因為變異,種群中才會存在多樣性。

下圖為變異的一個簡單示例:

變異完成之後,我們就得到了新為個體,進化也就完成了,整個過程如下圖:

在進行完一輪「遺傳變異」之後,我們用適應度函數對這些新的後代進行驗證,如果函數判定它們適應度足夠,那麼就會用它們從總體中替代掉那些適應度不夠的染色體。這裡有個問題,我們最終應該以什麼標準來判斷後代達到了最佳適應度水平呢?

一般來說,有如下幾個終止條件: 1. 在進行 X 次迭代之後,總體沒有什麼太大改變。 2. 我們事先為演算法定義好了進化的次數。 3. 當我們的適應度函數已經達到了預先定義的值。


4: 遺傳演算法應用

這節將講述如何利用遺傳演算法解決一個二元函數的最大值求解問題。

4.1 問題

二元函數如下: f(x,y)=0.5-frac{sin^2sqrt{x^2+y^2}-0.5}{1+0.001(x^2+y^2)^2}

# 畫出圖像如下from mpl_toolkits.mplot3d import Axes3Dimport numpy as npfrom matplotlib import pyplot as pltfig = plt.figure(figsize=(10,6))ax = Axes3D(fig)x = np.arange(-10, 10, 0.1)y = np.arange(-10, 10, 0.1)X, Y = np.meshgrid(x, y) Z = 0.5 - (np.sin(np.sqrt(X**2+Y**2))**2 - 0.5)/(1 + 0.001*(x**2 + y**2)**2)plt.xlabel(x)plt.ylabel(y)ax.set_zlim([-1,5])ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=rainbow)plt.show()

我們任務是找到 xin [-10,10],y in [-10, 10] 範圍之內的最大值。

4.2 創造染色體(編碼)

我們嘗試為上文所述的函數 f(x, y) 的最大值所對應的 xy 的值構造染色體。也就是說,要在一組二進位位中存儲函數 f(x, y) 的定義域中的數值信息。

假如設定求解的精度為小數點後6位,可以將區間[-10,10]劃分為 20	imes 10^6 個子區間。若用一組二進位位形式的染色體來表示這個數值集合, 2^{24}<20	imes 10^6<2^{25} ,需要25位二進位數來表示這些解。換句話說,一個解的編碼就是一個25位的二進位串。

現在,我們已經創建了一種 25 位長度的二進位位類型的染色體,那麼對於任意一個這樣的染色體,我們如何將其復原為[-10, 10]這個區間中的數值呢?這個很簡單,只需要使用下面的公式即可:

f(c)=-10 +ccdot frac{10-(-10)}{2^{25}-1}

例如 0000 0000 0000 0000 0000 0000 0000 0 和 1111 1111 1111 1111 1111 1111 1 這兩個二進位數,將其化為 10 進位數,代入上式,可得 -10.0 和 10.0。這意味著長度為 25 位的二進位數總是可以通過上式轉化為[-10, 10]區間中的數。

4.3 個體、種群與進化

染色體表達了某種特徵,這種特徵的載體,可以稱為「個體」。例如,我本人就是一個「個體」,我身上載有 23 對染色體,也許我的相貌、性別、性格等因素主要取決於它們。眾多個體便構成「種群」。

對於所要解決的二元函數最大值求解問題,個體可採用上一節所構造的染色體表示,並且數量為2個,其含義可理解為函數f(x, y)定義域內的一個點的坐標。許多這樣的個體便構成了一個種群,其含義為一個二維點集,包含於對角定點為(-10.0, -10.0)和(10.0, 10.0)的正方形區域。

也許有這樣一個種群,它所包含的個體對應的函數值會比其他個體更接近於函數f(x, y)的理論最大值,但是它一開始的時候可能並不比其他個體優秀,它之所以優秀是因為它選擇了不斷的進化,每一次的進化都要盡量保留種群中的優秀個體,淘汰掉不理想的個體,並且在優秀個體之間進行染色體交叉,有些個體還可能出現變異。種群的每一次進化後,必定會產生一個最優秀的個體。種群所有世代中的那個最優個體也許就是函數f(x, y)的最大值對應的定義域中的點。如果種群不休止的進化,它總是能夠找到最好的解。但是,由於我們的時間是有限的,有可能等不及種群的最優進化結果,通常是在得到了一個看上去還不錯的解時,便終止了種群的進化。

那麼,對於一個給定的種群,通常上述講過的選擇、交叉、變異來進行進化。

4.4 python實現

import math, randomclass Population: # 種群的設計 def __init__(self, size, chrom_size, cp, mp, gen_max): # 種群信息合 self.individuals = [] # 個體集合 self.fitness = [] # 個體適應度集 self.selector_probability = [] # 個體選擇概率集合 self.new_individuals = [] # 新一代個體集合 self.elitist = {chromosome:[0, 0], fitness:0, age:0} # 最佳個體的信息 self.size = size # 種群所包含的個體數 self.chromosome_size = chrom_size # 個體的染色體長度 self.crossover_probability = cp # 個體之間的交叉概率 self.mutation_probability = mp # 個體之間的變異概率 self.generation_max = gen_max # 種群進化的最大世代數 self.age = 0 # 種群當前所處世代 # 隨機產生初始個體集,並將新一代個體、適應度、選擇概率等集合以 0 值進行初始化 v = 2 ** self.chromosome_size - 1 for i in range(self.size): self.individuals.append([random.randint(0, v), random.randint(0, v)]) self.new_individuals.append([0, 0]) self.fitness.append(0) self.selector_probability.append(0) # 基於輪盤賭博機的選擇 def decode(self, interval, chromosome): 將一個染色體 chromosome 映射為區間 interval 之內的數值 d = interval[1] - interval[0] n = float (2 ** self.chromosome_size -1) return (interval[0] + chromosome * d / n) def fitness_func(self, chrom1, chrom2): 適應度函數,可以根據個體的兩個染色體計算出該個體的適應度 interval = [-10.0, 10.0] (x, y) = (self.decode(interval, chrom1), self.decode(interval, chrom2)) n = lambda x, y: math.sin(math.sqrt(x*x + y*y)) ** 2 - 0.5 d = lambda x, y: (1 + 0.001 * (x*x + y*y)) ** 2 func = lambda x, y: 0.5 - n(x, y)/d(x, y) return func(x, y) def evaluate(self): 用於評估種群中的個體集合 self.individuals 中各個個體的適應度 sp = self.selector_probability for i in range (self.size): self.fitness[i] = self.fitness_func (self.individuals[i][0], # 將計算結果保存在 self.fitness 列表中 self.individuals[i][1]) ft_sum = sum (self.fitness) for i in range (self.size): sp[i] = self.fitness[i] / float (ft_sum) # 得到各個個體的生存概率 for i in range (1, self.size): sp[i] = sp[i] + sp[i-1] # 需要將個體的生存概率進行疊加,從而計算出各個個體的選擇概率 # 輪盤賭博機(選擇) def select(self): (t, i) = (random.random(), 0) for p in self.selector_probability: if p > t: break i = i + 1 return i # 交叉 def cross(self, chrom1, chrom2): p = random.random() # 隨機概率 n = 2 ** self.chromosome_size -1 if chrom1 != chrom2 and p < self.crossover_probability: t = random.randint(1, self.chromosome_size - 1) # 隨機選擇一點(單點交叉) mask = n << t # << 左移運算符 (r1, r2) = (chrom1 & mask, chrom2 & mask) # & 按位與運算符:參與運算的兩個值,如果兩個相應位都為1,則該位的結果為1,否則為0 mask = n >> (self.chromosome_size - t) (l1, l2) = (chrom1 & mask, chrom2 & mask) (chrom1, chrom2) = (r1 + l2, r2 + l1) return (chrom1, chrom2) # 變異 def mutate(self, chrom): p = random.random () if p < self.mutation_probability: t = random.randint (1, self.chromosome_size) mask1 = 1 << (t - 1) mask2 = chrom & mask1 if mask2 > 0: chrom = chrom & (~mask2) # ~ 按位取反運算符:對數據的每個二進位位取反,即把1變為0,把0變為1 else: chrom = chrom ^ mask1 # ^ 按位異或運算符:當兩對應的二進位相異時,結果為1 return chrom # 保留最佳個體 def reproduct_elitist (self): # 與當前種群進行適應度比較,更新最佳個體 j = -1 for i in range (self.size): if self.elitist[fitness] < self.fitness[i]: j = i self.elitist[fitness] = self.fitness[i] if (j >= 0): self.elitist[chromosome][0] = self.individuals[j][0] self.elitist[chromosome][1] = self.individuals[j][1] self.elitist[age] = self.age # 進化過程 def evolve(self): indvs = self.individuals new_indvs = self.new_individuals # 計算適應度及選擇概率 self.evaluate() # 進化操作 i = 0 while True: # 選擇兩個個體,進行交叉與變異,產生新的種群 idv1 = self.select() idv2 = self.select() # 交叉 (idv1_x, idv1_y) = (indvs[idv1][0], indvs[idv1][1]) (idv2_x, idv2_y) = (indvs[idv2][0], indvs[idv2][1]) (idv1_x, idv2_x) = self.cross(idv1_x, idv2_x) (idv1_y, idv2_y) = self.cross(idv1_y, idv2_y) # 變異 (idv1_x, idv1_y) = (self.mutate(idv1_x), self.mutate(idv1_y)) (idv2_x, idv2_y) = (self.mutate(idv2_x), self.mutate(idv2_y)) (new_indvs[i][0], new_indvs[i][1]) = (idv1_x, idv1_y) # 將計算結果保存於新的個體集合self.new_individuals中 (new_indvs[i+1][0], new_indvs[i+1][1]) = (idv2_x, idv2_y) # 判斷進化過程是否結束 i = i + 2 # 循環self.size/2次,每次從self.individuals 中選出2個 if i >= self.size: break # 最佳個體保留 # 如果在選擇之前保留當前最佳個體,最終能收斂到全局最優解。 self.reproduct_elitist() # 更新換代:用種群進化生成的新個體集合 self.new_individuals 替換當前個體集合 for i in range (self.size): self.individuals[i][0] = self.new_individuals[i][0] self.individuals[i][1] = self.new_individuals[i][1] def run(self): 根據種群最大進化世代數設定了一個循環。 在循環過程中,調用 evolve 函數進行種群進化計算,並輸出種群的每一代的個體適應度最大值、平均值和最小值。 for i in range (self.generation_max): self.evolve () print (i, max (self.fitness), sum (self.fitness)/self.size, min (self.fitness))if __name__ == __main__: # 種群的個體數量為 50,染色體長度為 25,交叉概率為 0.8,變異概率為 0.1,進化最大世代數為 150 pop = Population (50, 24, 0.8, 0.1, 150) pop.run()


5: 參考文獻

一文讀懂遺傳演算法工作原理(附Python實現)

# emerge -e world

【演算法】超詳細的遺傳演算法(Genetic Algorithm)解析


推薦閱讀:
相关文章