import matplotlib.pyplot as plt import numpy as np import random import pandas as pd

遺傳算法求函數最值

遺傳算法的特點有無須可導，可優化離散異構數據。

參考：

• 遺傳算法python實現
• 科學前沿：Genetic Algorithm - 遺傳算法究竟怎么回事

問題定義

如下圖所示，

在區間[-1, 2]上有很多極大值和極小值，如果要求其最大值或最小值，很多單點優化的方法（梯度下降等）就不適合，這種情況下就可以用遺傳算法（Genetic Algorithm）。def fun(x):return x * np.sin(10*np.pi*x) + 2Xs = np.linspace(-1, 2, 100) plt.plot(Xs, fun(Xs))

初始化原始種群

遺傳算法的一個特點是同時優化一批解（種群），例如下面在[-1, 2]范圍內隨機生成

個點：np.random.seed(0)# 初始化原始種群 def ori_popular(num, _min=-1, _max=2):return np.random.uniform(_min, _max, num) # 范圍[-1, 2)#__TEST__ population = ori_popular(10)for pop, fit in zip(population, fun(population)):print("x=%5.2f, fit=%.2f"%(pop, fit))plt.plot(Xs, fun(Xs)) plt.plot(population, fun(population), '*') plt.show()>>> x= 0.65, fit=2.64 x= 1.15, fit=0.87 x= 0.81, fit=2.21 x= 0.63, fit=2.56 x= 0.27, fit=2.21 x= 0.94, fit=1.13 x= 0.31, fit=1.88 x= 1.68, fit=3.17 x= 1.89, fit=2.53 x= 0.15, fit=1.85

上圖顯示我們隨機選的

個點離最大值（3.8左右）差距還挺遠，下面用GA算法看能否求到最優解。

編碼

編碼，也就是由表現型到基因型，性征到染色體。

二進制編碼的缺點：對于一些連續函數的優化問題，由于其隨機性使得其局部搜索能力較差，如對于一些高精度的問題，當解迫近于最優解后，由于其變異后表現型變化很大，不連續，所以會遠離最優解，達不到穩定。而格雷碼能有效地防止這類現象。

TODO:

• [ ] 用2**18擴展似乎會損失精度，準備用10000代替

下面分別將1, 10, 0轉成二進制編碼，注意浮點數0.1無法轉換：

print("1:", bin(1)) print("10:", bin(10)) print("0:", bin(0))try:print("0.1:", bin(0.1)) except Exception as E:print("Exception: {}".format(type(E).__name__), E)>>> 1: 0b1 10: 0b1010 0: 0b0 Exception: TypeError 'float' object cannot be interpreted as an integer

為了能編碼浮點數，需要擴大倍數轉成整數：

X = [1, 0.1, 0.002]print("2**18 =", 2**18, 'n') for x in X:tx = int(x * 2**18)print("%.4f => %6d => %s"%(x, tx, bin(tx)))2**18 = 262144 >>> 1.0000 => 262144 => 0b1000000000000000000 0.1000 => 26214 => 0b110011001100110 0.0020 => 524 => 0b1000001100for x in X:tx = int(x * 2**18)ex = bin(tx)dx = int(ex, 2) / 2**18print("%25s => %6d => %.4f"%(ex, tx, dx))>>> 0b1000000000000000000 => 262144 => 1.00000b110011001100110 => 26214 => 0.10000b1000001100 => 524 => 0.0020

需要注意的是，上面用%.4f進行了四舍五入，實際上用2**18=262144來放大編碼會有精度的損失，例如：

x = 0.1 tx = int(x * 2**18) ex = bin(tx) int('0b110011001100110', 2) / (2**18)>>> 0.09999847412109375def encode(popular, _min=-1, _max=2, scale=2**18, bin_len=18): # popular應該是float類型的列表"""bin_len: 染色體序列長度"""norm_data = (popular-_min) / (_max-_min) * scalebin_data = np.array([np.binary_repr(x, width=bin_len) for x in norm_data.astype(int)]) # 轉成二進制編碼return bin_datachroms = encode(population) for pop, chrom, fit in zip(population, chroms, fun(population)):print("x=%.2f, chrom=%s, fit=%.2f"%(pop, chrom, fit))>>> x=0.65, chrom=100011000111111100, fit=2.64 x=1.15, chrom=101101110001011010, fit=0.87 x=0.81, chrom=100110100100111010, fit=2.21 x=0.63, chrom=100010110111110101, fit=2.56 x=0.27, chrom=011011000111010010, fit=2.21 x=0.94, chrom=101001010101100101, fit=1.13 x=0.31, chrom=011100000000010110, fit=1.88 x=1.68, chrom=111001000100101100, fit=3.17 x=1.89, chrom=111101101011001010, fit=2.53 x=0.15, chrom=011000100010100100, fit=1.85

解碼和適應度函數

通過基因（染色體）得到個體的適應度值，也就是評判當前種群每個個體（解）表現好壞，要對編碼解碼后才能代入函數。

注意，因為函數

的結果就是大小，所以直接用來當適應度函數。def decode(popular_gene, _min=-1, _max=2, scale=2**18):return np.array([(np.int(x, base=2)/scale*3)+_min for x in popular_gene])fitness = fun(decode(chroms))for pop, chrom, dechrom, fit in zip(population, chroms, decode(chroms), fitness):print("x=%5.2f, chrom=%s, dechrom=%.2f, fit=%.2f"%(pop, chrom, dechrom, fit))>>> x= 0.65, chrom=100011000111111100, dechrom=0.65, fit=2.64 x= 1.15, chrom=101101110001011010, dechrom=1.15, fit=0.87 x= 0.81, chrom=100110100100111010, dechrom=0.81, fit=2.21 x= 0.63, chrom=100010110111110101, dechrom=0.63, fit=2.56 x= 0.27, chrom=011011000111010010, dechrom=0.27, fit=2.21 x= 0.94, chrom=101001010101100101, dechrom=0.94, fit=1.13 x= 0.31, chrom=011100000000010110, dechrom=0.31, fit=1.88 x= 1.68, chrom=111001000100101100, dechrom=1.68, fit=3.17 x= 1.89, chrom=111101101011001010, dechrom=1.89, fit=2.53 x= 0.15, chrom=011000100010100100, dechrom=0.15, fit=1.85

這里將最小的適應度調整到0.000001，主要為了防止負數。

fitness = fitness - fitness.min() + 0.000001 # 保證所有的都為正

選擇與交叉

選擇就是淘汰不好的染色體（解），選擇方式有很多，這里用輪牌賭。

交叉就是讓染色體相互交換基因，方式有很多，這里在染色體中間進行交叉，交叉概率默認為0.66。

def SelectandCrossover(chroms, fitness, prob=0.6):probs = fitness/np.sum(fitness) # 各個個體被選擇的概率probs_cum = np.cumsum(probs) # 概率累加分布each_rand = np.random.uniform(size=len(fitness))#print(np.round(fitness, 2), np.round(probs, 2), np.round(probs_cum, 2), sep='n')#print([np.where(probs_cum > rand)[0][0] for rand in each_rand])# 根據隨機概率選擇出新的基因編碼newX = np.array([chroms[np.where(probs_cum > rand)[0][0]] for rand in each_rand])# 繁殖，隨機配對（概率為0.6)pairs = np.random.permutation(int(len(newX)*prob//2*2)).reshape(-1, 2)center = len(newX[0])//2for i, j in pairs:# 在中間位置交叉x, y = newX[i], newX[j]newX[i] = x[:center] + y[center:]newX[j] = y[:center] + x[center:]#print(x, y, newX[i], '', sep='n')return newX chroms = SelectandCrossover(chroms, fitness)dechroms = decode(chroms) fitness = fun(dechroms)for gene, dec, fit in zip(chroms, dechroms, fitness):print("chrom=%s, dec=%5.2f, fit=%.2f"%(gene, dec, fit))fig, (axs1, axs2) = plt.subplots(1, 2, figsize=(14, 5)) axs1.plot(Xs, fun(Xs)) axs1.plot(population, fun(population), 'o') axs2.plot(Xs, fun(Xs)) axs2.plot(dechroms, fitness, '*') plt.show()>>> chrom=111101101000010110, dec= 1.89, fit=2.64 chrom=011100000011001010, dec= 0.31, fit=1.86 chrom=011100000010100100, dec= 0.31, fit=1.86 chrom=011000100000010110, dec= 0.15, fit=1.85 chrom=100011000111111100, dec= 0.65, fit=2.64 chrom=100011000111111100, dec= 0.65, fit=2.64 chrom=100011000111111100, dec= 0.65, fit=2.64 chrom=111101101011001010, dec= 1.89, fit=2.53 chrom=111001000100101100, dec= 1.68, fit=3.17 chrom=111101101011001010, dec= 1.89, fit=2.53

變異

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

def Mutate(chroms:np.array):prob = 0.1clen = len(chroms[0])m = {'0':'1', '1':'0'}newchroms = []each_prob = np.random.uniform(size=len(chroms))for i, chrom in enumerate(chroms):if each_prob[i] < prob:pos = np.random.randint(clen) chrom = chrom[:pos] + m[chrom[pos]] + chrom[pos+1:]newchroms.append(chrom)return np.array(newchroms)newchroms = Mutate(chroms)def PltTwoChroms(chroms1, chroms2, fitfun):Xs = np.linspace(-1, 2, 100)fig, (axs1, axs2) = plt.subplots(1, 2, figsize=(14, 5))dechroms = decode(chroms1)fitness = fitfun(dechroms)axs1.plot(Xs, fitfun(Xs))axs1.plot(dechroms, fitness, 'o')dechroms = decode(chroms2)fitness = fitfun(dechroms)axs2.plot(Xs, fitfun(Xs))axs2.plot(dechroms, fitness, '*')plt.show()PltTwoChroms(chroms, newchroms, fun)

完整運行

上面我們分解了每一步的運算，下面我們用GA算法完整地迭代求出最優解，這里種群個數為100，迭代次數為1000：

np.random.seed(0) population = ori_popular(100) chroms = encode(population)for i in range(1000):fitness = fun(decode(chroms))fitness = fitness - fitness.min() + 0.000001 # 保證所有的都為正newchroms = Mutate(SelectandCrossover(chroms, fitness))if i % 300 == 1:PltTwoChroms(chroms, newchroms, fun)chroms = newchromsPltTwoChroms(chroms, newchroms, fun)

深入思考

其實GA每一步都有各種各樣的算法，包括超參數（種群個數、交叉概率、變異概率），沒有一套組合適合所有的情況，需要根據實際情況來調整。例如我們在進行“選擇與交叉”時，可以考慮保留5%最優秀的（因為有一定概率被淘汰），甚至這5%都不去交叉。

背包問題

背包問題是一個典型的優化組合問題：有一個容量固定的背包，已知

個物品的體積及其價值，我們該選擇哪些物品裝入背包使得在其容量約束限制之內所裝物品的價值最大？比如：背包容量限制為，有個物品的體積為，其價值則分別是。

背包問題在計算理論中屬于NP-完全問題，傳統的算法是動態規劃法，可求出最優總價值為

，其時間復雜度對于輸入是指數級的。而如果用貪心算法，計算物品性價比依次為，依次選擇物品和，總價值只能達到，這并非最優解。

現在我們試用遺傳算法求解這個問題。

編碼

將待求解的各量表示成長度為n的二進制串。例如：

代表一個解，它表示將第、號物體放入背包中，其它的物體則不放入。

生成初始種群

本例中，群體規模的大小取為

，即群體由個個體組成，每個個體可通過隨機的方法產生。例如：np.random.seed(0) X = np.random.randint(0, 2, size=(4,4)) X>>> array([[0, 1, 1, 0],[1, 1, 1, 1],[1, 1, 1, 0],[0, 1, 0, 0]])

種群適應度計算

對每個個體，計算所有物體的體積totalSize，所有物體的價值totalValue，以及個體適應度fit。若所有物體的總體積小于背包容量C，fit就是totalValue；否則應進行懲罰，fit應減去alpha?(totalSize?C)，這里我們取alpha=2，得到如下結果：

each_size = np.array([2, 3, 4, 5]) each_value = np.array([3, 4, 5, 7]) C = 9 alpha = 2each_fit = np.zeros_like(each_size) for i, x in enumerate(X):total_size = (x*each_size).sum()total_value = (x*each_value).sum()fit = total_value if (total_size <= C) else total_value - alpha*(total_size - C)each_fit[i] = fitprint("x=%s, total_size=%d, total_value=%d, fit=%d"%(x, total_size, total_value, fit))>>> x=[0 1 1 0], total_size=7, total_value=9, fit=9 x=[1 1 1 1], total_size=14, total_value=19, fit=9 x=[1 1 1 0], total_size=9, total_value=12, fit=12 x=[0 1 0 0], total_size=3, total_value=4, fit=4

選擇

根據染色體適應度，我們可以采用輪盤賭進行染色體的選擇，即與適應度成正比的概率來確定每個個體復制到下一代群體中的數量。具體操作過程如下：

先計算出群體中所有個體的適應度的總和。

然后用每個個體適應度除以總和，結果就是個體被遺傳到下一代群體中的概率。

每個概率值組成一個區間，全部概率值之和為。

最后再產生一個[0, 1]之間的隨機數，依據該隨機數出現在上述哪一個區間內來確定各個個體被選中的次數。

fit_sum = each_fit.sum() probs = each_fit / fit_sum probs_cum = probs.cumsum()for i, x in enumerate(X):print('probility of %s is %.2f%%'%(x, probs[i]*100)) print('cumulative sum of the probilities:', np.round(probs_cum, 4)*100)>>> probility of [0 1 1 0] is 26.47% probility of [1 1 1 1] is 26.47% probility of [1 1 1 0] is 35.29% probility of [0 1 0 0] is 11.76% cumulative sum of the probilities: [ 26.47 52.94 88.24 100. ]new_num = 4 each_rand = np.random.uniform(size=4)np.round(each_rand, 2)>>> array([0.96, 0.38, 0.79, 0.53])

我們隨機產生了

個[0-1)之間的數each_rand，用于在probs_cum中符合的區間上選中對應的基因編碼，分別選中了第88.24~100、26.47~52.94、52.94~88.24和26.47~52.94上的區間，所以新的基因編碼為：newX = np.array([ X[np.where(probs_cum > rand)[0][0]] for rand in each_rand])newX>>> array([[0, 1, 0, 0],[1, 1, 1, 1],[1, 1, 1, 0],[1, 1, 1, 1]])

沒想到概率最大的[1 1 1 0]反而被淘汰了。

交叉

采用單點交叉的方法，先對群體進行隨機配對，其次隨機設置交叉點位置，最后再相互交換配對染色體之間的部分基因。

突變

我們采用基本位變異的方法來進行變異運算，首先確定出各個個體的基因變異位置，然后依照某一概率將變異點的原有基因值取反。若突變之后最大的適應值個體已經符合期望的價值則停止，否則返回第3步繼續迭代。這里我們設突變的概率為

：mutation_prob = 0.1 for i, x in enumerate(newX):if np.random.uniform() < mutation_prob:pos = np.random.randint(0, 4)newX[i][pos] = abs(x[pos] - 1) newX>>> array([[0, 1, 0, 0],[1, 1, 1, 1],[1, 0, 1, 0],[1, 1, 1, 1]])

總結

以上是生活随笔為你收集整理的python 求函数最大值_遗传算法与Python图解的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。