python實現GA求二元函數最大值(來自知乎)

原文鏈接：https://zhuanlan.zhihu.com/p/43546261

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

問題

二元函數如下： 

# 畫出圖像如下
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
from matplotlib import pyplot as plt

fig = 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()

 

 

我們任務是找到  范圍之內的最大值。

 

創造染色體（編碼）

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

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

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

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

 

個體、種群與進化

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

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

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

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

 

python實現

 

import math, random

class 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()