遺傳算法python實現

最近看了一下遺傳算法，使用輪盤賭選擇染色體，使用單點交叉，下面是代碼實現（python3）

 

  1 import numpy as np
  2 import random
  3 from scipy.optimize import fsolve
  4 import matplotlib.pyplot as plt
  5 import heapq
  6 
  7 # 求染色體長度
  8 def getEncodeLength(decisionvariables, delta):
  9     # 將每一個變量的編碼長度放入數組
 10     lengths = []
 11     for decisionvar in decisionvariables:
 12         uper = decisionvar[1]
 13         low = decisionvar[0]
 14         # res()返回一個數組
 15         res = fsolve(lambda x: ((uper - low) / delta - 2 ** x + 1), 30)
 16         # ceil()向上取整
 17         length = int(np.ceil(res[0]))
 18         lengths.append(length)
 19     # print("染色體長度:", lengths)
 20     return lengths
 21 
 22 
 23 # 隨機生成初始化種羣
 24 def getinitialPopulation(length, populationSize):
 25     chromsomes = np.zeros((populationSize, length), dtype=np.int)
 26     for popusize in range(populationSize):
 27         # np.random.randit()產生[0,2)之間的隨機整數，第三個參數表示隨機數的數量
 28         chromsomes[popusize, :] = np.random.randint(0, 2, length)
 29     return chromsomes
 30 
 31 
 32 # 染色體解碼獲得表現形的解
 33 def getDecode(population, encodelength, decisionvariables, delta):
 34     # 獲得population中有幾個元素
 35     populationsize = population.shape[0]
 36     length = len(encodelength)
 37     decodeVariables = np.zeros((populationsize, length), dtype=np.float)
 38     # 將染色體拆分添加到解碼數組decodeVariables中
 39     for i, populationchild in enumerate(population):
 40         # 設置起始點
 41         start = 0
 42         for j, lengthchild in enumerate(encodelength):
 43             power = lengthchild - 1
 44             decimal = 0
 45             for k in range(start, start + lengthchild):
 46                 # 二進制轉爲十進制
 47                 decimal += populationchild[k] * (2 ** power)
 48                 power = power - 1
 49             # 從下一個染色體開始
 50             start = lengthchild
 51             lower = decisionvariables[j][0]
 52             uper = decisionvariables[j][1]
 53             # 轉換爲表現形
 54             decodevalue = lower + decimal * (uper - lower) / (2 ** lengthchild - 1)
 55             # 將解添加到數組中
 56             decodeVariables[i][j] = decodevalue
 57     return decodeVariables
 58 
 59 
 60 # 獲得每一個個體的適應度值及累計機率
 61 def getFitnessValue(func, decode):
 62     # 獲得種羣的規模和決策變量的個數
 63     popusize, decisionvar = decode.shape
 64     # 初始化適應度值空間
 65     fitnessValue = np.zeros((popusize, 1))
 66     for popunum in range(popusize):
 67         fitnessValue[popunum][0] = func(decode[popunum][0], decode[popunum][1])
 68     # 獲得每一個個體被選擇的機率
 69     probability = fitnessValue / np.sum(fitnessValue)
 70     # 獲得每一個染色體被選中的累積機率，用於輪盤賭算子使用
 71     cum_probability = np.cumsum(probability)
 72     return fitnessValue, cum_probability
 73 
 74 
 75 # 選擇新的種羣
 76 def selectNewPopulation(decodepopu, cum_probability):
 77     # 獲取種羣的規模和
 78     m, n = decodepopu.shape
 79     # 初始化新種羣
 80     newPopulation = np.zeros((m, n))
 81     for i in range(m):
 82         # 產生一個0到1之間的隨機數
 83         randomnum = np.random.random()
 84         # 輪盤賭選擇
 85         for j in range(m):
 86             if (randomnum < cum_probability[j]):
 87                 newPopulation[i] = decodepopu[j]
 88                 break
 89     return newPopulation
 90 
 91 
 92 # 新種羣交叉
 93 def crossNewPopulation(newpopu, prob):
 94     m, n = newpopu.shape
 95     # uint8將數值轉換爲無符號整型
 96     numbers = np.uint8(m * prob)
 97     # 若是選擇的交叉數量爲奇數，則數量加1
 98     if numbers % 2 != 0:
 99         numbers = numbers + 1
100     # 初始化新的交叉種羣
101     updatepopulation = np.zeros((m, n), dtype=np.uint8)
102     # 隨機生成須要交叉的染色體的索引號
103     index = random.sample(range(m), numbers)
104     # 不須要交叉的染色體直接複製到新的種羣中
105     for i in range(m):
106         if not index.__contains__(i):
107             updatepopulation[i] = newpopu[i]
108     # 交叉操做
109     j = 0
110     while j < numbers:
111         # 隨機生成一個交叉點，np.random.randint()返回的是一個列表
112         crosspoint = np.random.randint(0, n, 1)
113         crossPoint = crosspoint[0]
114         # a = index[j]
115         # b = index[j+1]
116         updatepopulation[index[j]][0:crossPoint] = newpopu[index[j]][0:crossPoint]
117         updatepopulation[index[j]][crossPoint:] = newpopu[index[j + 1]][crossPoint:]
118         updatepopulation[index[j + 1]][0:crossPoint] = newpopu[j + 1][0:crossPoint]
119         updatepopulation[index[j + 1]][crossPoint:] = newpopu[index[j]][crossPoint:]
120         j = j + 2
121     return updatepopulation
122 
123 
124 # 變異操做
125 def mutation(crosspopulation, mutaprob):
126     # 初始化變異種羣
127     mutationpopu = np.copy(crosspopulation)
128     m, n = crosspopulation.shape
129     # 計算須要變異的基因數量
130     mutationnums = np.uint8(m * n * mutaprob)
131     # 隨機生成變異基因的位置
132     mutationindex = random.sample(range(m * n), mutationnums)
133     # 變異操做
134     for geneindex in mutationindex:
135         # np.floor()向下取整返回的是float型
136         row = np.uint8(np.floor(geneindex / n))
137         colume = geneindex % n
138         if mutationpopu[row][colume] == 0:
139             mutationpopu[row][colume] = 1
140         else:
141             mutationpopu[row][colume] = 0
142     return mutationpopu
143 
144 # 找到從新生成的種羣中適應度值最大的染色體生成新種羣
145 def findMaxPopulation(population, maxevaluation, maxSize):
146     #將數組轉換爲列表
147     maxevalue = maxevaluation.flatten()
148     maxevaluelist = maxevalue.tolist()
149     # 找到前100個適應度最大的染色體的索引
150     maxIndex = map(maxevaluelist.index, heapq.nlargest(100, maxevaluelist))
151     index = list(maxIndex)
152     colume = population.shape[1]
153     # 根據索引生成新的種羣
154     maxPopulation = np.zeros((maxSize, colume))
155     i = 0
156     for ind in index:
157         maxPopulation[i] = population[ind]
158         i = i + 1
159     return maxPopulation
160 
161 # 適應度函數，使用lambda能夠不用在函數總傳遞參數
162 def fitnessFunction():
163     return lambda a, b: 21.5 + a * np.sin(4 * np.pi * a) + b * np.sin(20 * np.pi * b)
164 
165 
166 def main():
167     optimalvalue = []
168     optimalvariables = []
169 
170     # 兩個決策變量的上下界，多維數組之間必須加逗號
171     decisionVariables = [[-3.0, 12.1], [4.1, 5.8]]
172     # 精度
173     delta = 0.0001
174     # 獲取染色體長度
175     EncodeLength = getEncodeLength(decisionVariables, delta)
176     # 種羣數量
177     initialPopuSize = 100
178     # 初始生成100個種羣
179     population = getinitialPopulation(sum(EncodeLength), initialPopuSize)
180     # 最大進化代數
181     maxgeneration = 100
182     # 交叉機率
183     prob = 0.8
184     # 變異機率
185     mutationprob = 0.01
186     # 新生成的種羣數量
187     maxPopuSize = 100
188 
189     for generation in range(maxgeneration):
190         # 對種羣解碼獲得表現形
191         decode = getDecode(population, EncodeLength, decisionVariables, delta)
192         # 獲得適應度值和累計機率值
193         evaluation, cum_proba = getFitnessValue(fitnessFunction(), decode)
194         # 選擇新的種羣
195         newpopulations = selectNewPopulation(population, cum_proba)
196         # 新種羣交叉
197         crossPopulations = crossNewPopulation(newpopulations, prob)
198         # 變異操做
199         mutationpopulation = mutation(crossPopulations, mutationprob)
200         # 將父母和子女合併爲新的種羣
201         totalpopulation = np.vstack((population, mutationpopulation))
202         # 最終解碼
203         final_decode = getDecode(totalpopulation, EncodeLength, decisionVariables, delta)
204         # 適應度評估
205         final_evaluation, final_cumprob = getFitnessValue(fitnessFunction(), final_decode)
206         #選出適應度最大的100個從新生成種羣
207         population = findMaxPopulation(totalpopulation, final_evaluation, maxPopuSize)
208         # 找到本輪中適應度最大的值
209         optimalvalue.append(np.max(final_evaluation))
210         index = np.where(final_evaluation == max(final_evaluation))
211         optimalvariables.append(list(final_decode[index[0][0]]))
212 
213     x = [i for i in range(maxgeneration)]
214     y = [optimalvalue[i] for i in range(maxgeneration)]
215     plt.plot(x, y)
216     plt.show()
217 
218     optimalval = np.max(optimalvalue)
219     index = np.where(optimalvalue == max(optimalvalue))
220     optimalvar = optimalvariables[index[0][0]]
221     return optimalval, optimalvar
222 
223 
224 if __name__ == "__main__":
225     optval, optvar = main()
226 
227     print("f(x1,x2) = 21.5+x1*sin(4*pi*x1)+x2*sin(20*pi*x2)")
228     print("x1:", optvar[0])
229     print("X2:", optvar[1])
230     print("maxValue:", optval)