遺傳演算法（Genetic Algorithm, GA）是一種通用的優化演算法，屬於進化演算法簇中一個比較實用又有名的演算法。進化演算法融合了自然生物進化中共有的一些行為：繁殖、變異、競爭、選擇等。

基本流程

　　GA通過迭代來優化目標函式的引數，直到目標函式滿足一定條件時結束。迭代對目標函式的連續性並無要求，也就是說演算法的迭代並不基於目標函式在當前引數下的梯度等連續性質。演算法的轉化思想是將所有待優化引數看做生物染色體，單個引數看做染色體上對應位置的基因，而目標函式看做生物對生存環境的適應度。大致的迭代流程如下：

　　0、隨機生成一定數量的生物個體，它們的染色體各不相同，因此在生存環境中的適應度也不同。也就是說，這些隨機生成的引數對目標函式的滿足度各不相同。

　　1、根據個體對環境的適應度，適者生存。選擇適應度較高的一部分個體，生存下來，有資格繁殖下一代，剩下的都死在這一代。當然這裡也可以隨機少量選擇一些“幸運”的適應度較低的個體，讓它們也能參與繁殖，增加多樣性。

　　2、對選出的個體進行兩兩“交配”，交配就是讓它們DNA中對應位置的基因進行交叉和變異產生新的個體（不知道交叉和變異的回去讀高中生物書），得到下一代。具體的交配和交叉變異策略下面再詳細介紹。

　　3、在新的一代，對所有個體進行判斷，是否已經有個體滿足目標函式的終極條件。有就結束迭代，否則回到1。

演算法基本概念

　　串：表示生物個體的染色體，通常是二進位制串，或根據優化需要使用其它形式的串。

　　群體和群體大小：個體的集合稱為群體，串是群體的元素。群體中個體的數量稱為群體大小。

　　基因：基因是串中的元素，用來表示個體特徵。通常一個基因對應一個引數，或者在二進位制編碼時，一位就代表一個基因。

　　基因位置：基因在串中的位置，簡稱基因位。

　　串結構空間：基因任意組合所形成的串集合，也就是在當前優化背景下所有可能的串的集合。

　　適應度：表示某一個體對環境的適應程度，用一個函式表示，通常就是待優化的目標函式。

演算法組成部分

　　介紹完GA的大致流程和基本概念，你應該對其有了個大致瞭解。下面開始詳細講解演算法的各個組成部分，由四個部分組成：編碼機制、適應度函式、控制引數、遺傳運算元。

編碼機制

　　待優化的引數可以編碼為各種形式的串，以便進行上述的遺傳迭代。編碼要是可逆的，也就是說串與引數之間應該是一一對應的，從而每次迭代都能解碼並計算個體的適應度。並且GA的編碼可以有個十分廣泛的理解，不侷限於簡介中對待優化引數的編碼。比如，在分類問題，串可解釋為一個規則，即串的前半部為輸入或前件，後半部為輸出或後件、結論。因此，遺傳演算法可以對大量的問題進行建模，只要把問題的條件啊，約束啊，引數啊全都一股腦地編碼丟進去“進化”就好了，這也是為什麼它在數學建模比賽中如此熱門的原因。

　　編碼方式有二進位制編碼、格雷編碼、實數編碼等。二進位制編碼就是把引數都轉換為二進位制，然後排成一排，編碼的每一位都看做一個基因。實數編碼實際上就是不編碼，當引數都為實數而且很多時，就可以直接把這些實數引數作為染色體基因。下面介紹格雷碼。

格雷編碼

　　格雷編碼是對二進位制編碼進行變換後所得到的一種編碼方法，它要求兩個連續整數的編碼之間只能有一個碼位不同，其餘碼位都是完全相同的。因此，格雷碼解決了二進位制編碼的漢明懸崖問題。二進位制編碼與格雷碼之間的轉換如下：

　　設有二進位制串$b_1,b_2,...,b_n$，對應的格雷串為$a_1,a_2,…,a_n$，則從二進位制編碼到格雷編碼的變換為：

$ a_i = \left\{ \begin{aligned} &b_i,\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;i=1\\ &b_{i-1}\oplus b_i,\;\;\;\;i>1 \end{aligned} \right. $

　　其中$\oplus$為異或運算。格雷碼到二進位制編碼的變換為：

 $b_i = (\sum\limits_{j=1}^ia_i)mod\; 2$

　　儘管轉換式子有了，但是這樣看不出格雷碼的特點，下面通過列舉格雷碼來總結它的特點。

　　可以發現，格雷碼每增加一位，低位（僅除了最高位外）編碼即為所有位數少於它的格雷碼的反序排列。比如上面的四位格雷碼，低三位分別是100、101、111、110、010、011、001、000，觀察小於四位的格雷碼，從7到0，格雷碼同樣也是100、101、111、110、010、011、001、000。

適應度函式

　　適應度函式用於判斷個體對環境的適應度、評估個體，定出優劣程度。

　　原始適應度函式直接將待求解問題的目標函式$f(x)$定義為演算法的適應度函式。比如求解極值問題

$\max\limits_{x\in [a,b]}f(x)$

　　直接將$f(x)$作為原始適應度函式。

　　還有其他的定義方式，但是因為對於不同問題，定義方式各不相同，所以不再舉例，只要適應度函式是向著優化目標定義的即可（物種適應度越大，越接近你的優化目標）。

控制引數

　　在GA的實際操作時，需適當確定某些引數的值以提高選優的效果。這些引數包含：

　　字串所含字元的個數，即串長。這一長度為常數，記為$L$。

　　每一代群體的大小，也稱群體的容量，記為$n$。

　　交換率，即施行交換運算元的概率，記為$P_c$。

　　突變率，即施行突變運算元的概率，記為$P_m$。

遺傳運算元

　　遺傳運算元就是執行選擇、交叉、變異等操作的方法。GA中最常用的運算元有三種：選擇運算元、交叉運算元、變異運算元。

選擇運算元

　　選擇運算元從群體中按某一概率成對選擇個體，某個體$x_i$被選擇的概率$P_i$與其適應度值成正比。最通常的實現方法是輪盤賭型。

　　所謂輪盤賭，就是根據適應度函式來確定個體被選擇的概率。如果適應度函式$f(x)>0$，則可以這樣定義：

$\displaystyle P(x_i) = \frac{f(x_i)}{\sum\limits_{j=1}^Nf(x_j)}$

　　否則定義成softmax也行（我覺得）：

$\displaystyle P(x_i) = \frac{\exp f(x_i)}{\sum\limits_{j=1}^N\exp f(x_j)}$

　　然後根據概率選擇個體。這就是最簡單的概率選擇，不知道為什麼要扯出一個所謂輪盤賭來，是怕別人理解不了嗎。

交叉運算元

　　交叉運算元將被選中的兩個個體的基因鏈按概率$P_c$進行交叉，生成兩個新的個體。其中$P_c$是一個系統引數，通常看群體的大小，群體越小，$P_c$可以取得更大一些，但一般小於0.5。這是因為大於0.5等於把大部分基因都交換了，而這和交換剩餘的小部分基因等價。具體來說，交叉有所謂的單點交叉、兩點交叉、多點交叉、均勻交叉等。這些分類感覺沒有任何意義，直接就按$P_c$進行交叉就行了。過程如下：

　　首先生成長度為$L$的二進位制模板串，其中每位為1的概率為$P_c$。然後按照這個串對選擇出來的兩個父串進行交叉，為1的位置交叉，為0不交叉。比如對於兩個父串

$X = 1011010$

$Y = 1101001$

　　模板串為

$T = 0101010$

　　則交叉後變為

$X'=1111000$

$Y'=1001011$

變異運算元

　　變異運算元將新個體的基因鏈的各位按概率$P_m$進行變異，$P_m$通常取得很小比如0.0001。對二值基因鏈(0,1編碼)來說即是取反，對於實數基因來說，則可以在這個數的範圍內取隨機值，或是在某個固定範圍內取隨機值。另外，還可以進行基因之間的數字的調換、迴圈移動等。

總結

　　以上這些操作對於具體的優化任務而言，並沒有什麼原理可以講解，也沒有理論依據表明交叉和變異的確可以讓個體更優。但是從巨集觀上來看，群體的適應度的確會在選擇操作下不斷向好發展，最終達到一個良好的點。歸根結底，我認為這個演算法就是在保持當前較好結果的前提下的隨機迭代。它是十分隨機的，有時可能迭代幾步就搞定了，而有時可能幾千步都出不來，所以會浪費大量的計算資源。而且對於很複雜的非凸問題，它同樣也可能會陷入區域性最優。所以能用傳統方法求解的問題，儘量不用這種啟發式演算法。

實驗

　　實驗使用遺傳演算法求解以下函式的最大值：

$f(x,y) = 2\exp\left[\frac{-(x+3)^2-(y-3)^2}{10}\right] + 1.2\exp\left[\frac{-(x-3)^2-(y+3)^2}{10}\right] + \frac{1}{5}\exp\left[-\cos(3x)-\sin(3y)\right]$

　　$x,y$的範圍都是$(-100,100)$，函式三維影象如下：

　　看錶達式和影象可以判斷最大值在$(-3,3)$附近。下面先以均勻分佈產生1000個點，然後進行40次迭代。1000個點的迭代過程如下：

　　每次迭代的最大值和座標：

　　這裡交叉率是0.2，變異率是0.3，變異使用的正態分佈標準差是10，因此比較容易跳出這個函式的區域性最優點。比如函式在$(3,-3)$附近的區域性最優點，這點和$(-3,3)$之間差了$6\sqrt{2}$左右，剛好是10左右。另外，在$100\times 100$的範圍內，把種群的大小設定為1000實際上過於飽和了，隨隨便便就找到了最優解，這和直接網格暴力搜尋都有的一拼了。現在把種群大小調為200，變異標準差調為20，交叉率調為0.4，變異率調為0.1，看看迭代幾次可以達到最優。以下是迭代過程圖：

　　迭代了150次，也到了$(-3,3)$附近，時間上消耗少了一些。如果要再精確可能還要更多次的迭代，但因為它的變異方差較大，所以比較難。

　　以下是實驗程式碼：

#%%定義函式
import time
import numpy as np 
 
def func(X,Y):
  return np.exp((-(X+3)**2-(Y-3)**2)/10)*2 + np.exp((-(X-3)**2-(Y+3)**2)/10)*1.2 + np.exp(-np.cos(X*3)-np.sin(Y*3))/5
 
#%%遺傳演算法
time_start=time.time()
dim = 2
iterate_times = 150
individual_n = 200
p_c = 0.4
p_m = 0.15
std_var = 20
num_range = [-100,100]
def initialize(individual_n,num_range):
  return np.random.uniform(num_range[0],num_range[1],size = [individual_n,dim]) 
def overlapping(a,b,pc,dim):
  T = np.random.choice(a=2, size=dim,p=[1-pc,pc])
  for i in range(dim):
    if T[i] == 1:
      a[i],b[i] = b[i],a[i]
def variation(a,pm,dim,num_range):
  T = np.random.choice(a=2, size=dim,p=[1-pm,pm]) 
  for i in range(dim):
    if T[i] == 1:
      a[i] = np.clip(np.random.normal(loc = a[i], scale=std_var),num_range[0],num_range[1])
def get_select_p(v):
  exp_v = np.exp(v)
  sum_exp_v = np.sum(exp_v)
  return exp_v/sum_exp_v
def select_and_reproduction(population,num_range): 
  values = func(population[:,0],population[:,1])
  select_p = get_select_p(values)

  choices = np.random.choice(a=len(population), size=len(population), replace=True, p=select_p)
  new_population = population[choices]
  for i in range(int(len(new_population)/2)):
    overlapping(new_population[2*i],new_population[2*i+1],p_c,dim) 
    variation(new_population[2*i],p_m,dim,num_range)
    variation(new_population[2*i+1],p_m,dim,num_range)
  argmax_value = np.argmax(values)
  return values[argmax_value],population[argmax_value],new_population

population = initialize(individual_n,num_range)
times_population = np.empty([iterate_times,individual_n,dim])
max_value = 0
max_value_coordinate = []
for i in range(iterate_times):
  times_population[i] = population
  last_max_value,last_max_value_pos,population = select_and_reproduction(population,num_range)
  if last_max_value>max_value:
    max_value = last_max_value
    max_value_coordinate = last_max_value_pos
  print(i,":  ",last_max_value) 
  
time_end=time.time()
print('總共耗時：',time_end-time_start)
print('迭代最大值：',max_value)
print('對應座標：',max_value_coordinate)
#%%畫出種群的迭代過程 
import matplotlib.pyplot as plt
import matplotlib.animation as animation

fig, ax = plt.subplots()
xdata, ydata = [], []
ln, = plt.plot([], [],'ro',animated=True,color='black',markersize = 2) 
def init():
  ax.set_xlim(num_range[0],num_range[1])
  ax.set_ylim(num_range[0],num_range[1])
  ax.set_xlabel('X')
  ax.set_ylabel('Y')
  return ln, 
def update(frame):  
  frame = int(frame) 
  xdata = times_population[frame,:,0]
  ydata = times_population[frame,:,1]
  ln.set_data(xdata, ydata) 
  return ln, 

anim = animation.FuncAnimation(fig, update, frames=np.linspace(0,iterate_times-1,iterate_times),interval=100,
                    init_func=init,blit=True)
plt.show()