姓名：周小蓬 16019110037

轉載自：http://blog.csdn.net/qq_30666517/article/details/78637255

[嵌牛導讀]

本文用Python3完整實現了簡單遺傳算法(SGA)

Simple Genetic Alogrithm是模擬生物進化過程而提出的一種優(yōu)化算法。SGA采用隨機導向搜索全局最優(yōu)解或者說近似全局最優(yōu)解。傳統的爬山算法(例如梯度下降，牛頓法)一次只優(yōu)化一個解，并且對于多峰的目標函數很容易陷入局部最優(yōu)解，而SGA算法一次優(yōu)化一個種群(即一次優(yōu)化多個解)，SGA比傳統的爬山算法更容易收斂到全局最優(yōu)解或者近似全局最優(yōu)解。

[嵌牛鼻子]

遺傳算法

[嵌牛提問]

什么事遺傳算法，如何用python實現

[嵌牛正文]

SGA基本流程如下:

1、對問題的解進行二進制編碼。編碼涉及精度的問題，在本例中精度delta=0.0001,根據決策變量的上下界確定對應此決策變量的染色體基因的長度（m）。假設一個決策變量x0上界為upper,下界為lower,則精度delta = (upper-lower)/2^m-1。如果已知決策變量邊界和編碼精度，那么可以用下面的公式確定編碼決策變量x0所對應的染色體長度：

2^(length-1)<(upper-lower)/delta<=2^length-1

2、對染色體解碼得到表現形:

解碼后得到10進制的值；decoded = lower + binary2demical(chromosome)*delta。其中binary2demical為二進制轉10進制的函數，在代碼中有實現，chromosome是編碼后的染色體。

3、確定初始種群,初始種群隨機生成

4、根據解碼函數得到初始種群的10進制表現型的值

5、確定適應度函數，對于求最大值最小值問題，一般適應度函數就是目標函數。根據適應度函數確定每個個體的適應度值Fi=FitnessFunction(individual);然后確定每個個體被選擇的概率Pi=Fi/sum(Fi),sum(Fi)代表所有個體適應度之和。

6、根據輪盤賭選擇算子，選取適應度較大的個體。一次選取一個個體，選取n次,得到新的種群population

7、確定交叉概率Pc,對上一步得到的種群進行單點交叉。每次交叉點的位置隨機。

8、確定變異概率Pm,假設種群大小為10，每個個體染色體編碼長度為33，則一共有330個基因位，則變異的基因位數是330*Pm。接下來，要確定是那個染色體中哪個位置的基因發(fā)生了變異。將330按照10進制序號進行編碼即從0,1,2,.......229。隨機從330個數中選擇330*Pm個數，假設其中一個數時154，chromosomeIndex = 154/33 =4,

geneIndex = 154%33 = 22。由此確定了第154號位置的基因位于第4個染色體的第22個位置上，將此位置的基因值置反完成基本位變異操作。

9、以上步驟完成了一次迭代的所有操作。接下就是評估的過程。對變異后得到的最終的種群進行解碼，利用解碼值求得每個個體的適應度值，將最大的適應度值保存下來，對應的解碼后的決策變量的值也保存下來。

10、根據迭代次數，假設是500次，重復執(zhí)行1-9的步驟，最終得到是一個500個數值的最優(yōu)適應度取值的數組以及一個500*n的決策變量取值數組（假設有n個決策變量）。從500個值中找到最優(yōu)的一個（最大或者最小，根據定義的適應度函數來選擇）以及對應的決策變量的取值。

對于以上流程不是很清楚的地方，在代碼中有詳細的注釋。也可以自行查找資料補充理論。本文重點是實現

本代碼實現的問題是: maxf(x1,x2) = 21.5+x1*sin(4*pi*x1)+x2*sin(20*pi*x2)

s.t. -3.0<=x1<=12.1

4.1<=x2<=5.8

初始種群的編碼結果如下圖所示:

初始種群的解碼結果如下圖所示:

適應度值如圖所示:

輪盤賭選擇后的種群如圖所示；

單點交叉后的種群如圖所示:

基本位變異后的種群如圖所示；

最終結果如下圖所示；

源代碼如下；

[python]view plaincopy

#?!/usr/bin/env?python

#?-*-?coding:utf-8?-*-

#?Author:?wsw

#?簡單實現SGA算法

importnumpy?as?np

fromscipy.optimizeimportfsolve,?basinhopping

importrandom

importtimeit

#?根據解的精度確定染色體(chromosome)的長度

#?需要根據決策變量的上下邊界來確定

defgetEncodedLength(delta=0.0001,?boundarylist=[]):

#?每個變量的編碼長度

lengths?=?[]

foriinboundarylist:

lower?=?i[0]

upper?=?i[1]

#?lamnda?代表匿名函數f(x)=0,50代表搜索的初始解

res?=?fsolve(lambdax:?((upper?-?lower)?*1/?delta)?-2**?x?-1,50)

length?=?int(np.floor(res[0]))

lengths.append(length)

returnlengths

pass

#?隨機生成初始編碼種群

defgetIntialPopulation(encodelength,?populationSize):

#?隨機化初始種群為0

chromosomes?=?np.zeros((populationSize,?sum(encodelength)),?dtype=np.uint8)

foriinrange(populationSize):

chromosomes[i,?:]?=?np.random.randint(0,2,?sum(encodelength))

#?print('chromosomes?shape:',?chromosomes.shape)

returnchromosomes

#?染色體解碼得到表現型的解

defdecodedChromosome(encodelength,?chromosomes,?boundarylist,?delta=0.0001):

populations?=?chromosomes.shape[0]

variables?=?len(encodelength)

decodedvalues?=?np.zeros((populations,?variables))

fork,?chromosomeinenumerate(chromosomes):

chromosome?=?chromosome.tolist()

start?=0

forindex,?lengthinenumerate(encodelength):

#?將一個染色體進行拆分，得到染色體片段

power?=?length?-1

#?解碼得到的10進制數字

demical?=0

foriinrange(start,?length?+?start):

demical?+=?chromosome[i]?*?(2**?power)

power?-=1

lower?=?boundarylist[index][0]

upper?=?boundarylist[index][1]

decodedvalue?=?lower?+?demical?*?(upper?-?lower)?/?(2**?length?-1)

decodedvalues[k,?index]?=?decodedvalue

#?開始去下一段染色體的編碼

start?=?length

returndecodedvalues

#?得到個體的適應度值及每個個體被選擇的累積概率

defgetFitnessValue(func,?chromosomesdecoded):

#?得到種群規(guī)模和決策變量的個數

population,?nums?=?chromosomesdecoded.shape

#?初始化種群的適應度值為0

fitnessvalues?=?np.zeros((population,1))

#?計算適應度值

foriinrange(population):

fitnessvalues[i,0]?=?func(chromosomesdecoded[i,?:])

#?計算每個染色體被選擇的概率

probability?=?fitnessvalues?/?np.sum(fitnessvalues)

#?得到每個染色體被選中的累積概率

cum_probability?=?np.cumsum(probability)

returnfitnessvalues,?cum_probability

#?新種群選擇

defselectNewPopulation(chromosomes,?cum_probability):

m,?n?=?chromosomes.shape

newpopulation?=?np.zeros((m,?n),?dtype=np.uint8)

#?隨機產生M個概率值

randoms?=?np.random.rand(m)

fori,?randomainenumerate(randoms):

logical?=?cum_probability?>=?randoma

index?=?np.where(logical?==1)

#?index是tuple,tuple中元素是ndarray

newpopulation[i,?:]?=?chromosomes[index[0][0],?:]

returnnewpopulation

pass

#?新種群交叉

defcrossover(population,?Pc=0.8):

"""

:param?population:?新種群

:param?Pc:?交叉概率默認是0.8

:return:?交叉后得到的新種群

"""

#?根據交叉概率計算需要進行交叉的個體個數

m,?n?=?population.shape

numbers?=?np.uint8(m?*?Pc)

#?確保進行交叉的染色體個數是偶數個

ifnumbers?%2!=0:

numbers?+=1

#?交叉后得到的新種群

updatepopulation?=?np.zeros((m,?n),?dtype=np.uint8)

#?產生隨機索引

index?=?random.sample(range(m),?numbers)

#?不進行交叉的染色體進行復制

foriinrange(m):

ifnotindex.__contains__(i):

updatepopulation[i,?:]?=?population[i,?:]

#?crossover

whilelen(index)?>0:

a?=?index.pop()

b?=?index.pop()

#?隨機產生一個交叉點

crossoverPoint?=?random.sample(range(1,?n),1)

crossoverPoint?=?crossoverPoint[0]

#?one-single-point?crossover

updatepopulation[a,0:crossoverPoint]?=?population[a,0:crossoverPoint]

updatepopulation[a,?crossoverPoint:]?=?population[b,?crossoverPoint:]

updatepopulation[b,0:crossoverPoint]?=?population[b,0:crossoverPoint]

updatepopulation[b,?crossoverPoint:]?=?population[a,?crossoverPoint:]

returnupdatepopulation

pass

#?染色體變異

defmutation(population,?Pm=0.01):

"""

:param?population:?經交叉后得到的種群

:param?Pm:?變異概率默認是0.01

:return:?經變異操作后的新種群

"""

updatepopulation?=?np.copy(population)

m,?n?=?population.shape

#?計算需要變異的基因個數

gene_num?=?np.uint8(m?*?n?*?Pm)

#?將所有的基因按照序號進行10進制編碼，則共有m*n個基因

#?隨機抽取gene_num個基因進行基本位變異

mutationGeneIndex?=?random.sample(range(0,?m?*?n),?gene_num)

#?確定每個將要變異的基因在整個染色體中的基因座(即基因的具體位置)

forgeneinmutationGeneIndex:

#?確定變異基因位于第幾個染色體

chromosomeIndex?=?gene?//?n

#?確定變異基因位于當前染色體的第幾個基因位

geneIndex?=?gene?%?n

#?mutation

ifupdatepopulation[chromosomeIndex,?geneIndex]?==0:

updatepopulation[chromosomeIndex,?geneIndex]?=1

else:

updatepopulation[chromosomeIndex,?geneIndex]?=0

returnupdatepopulation

pass

#?定義適應度函數

deffitnessFunction():

returnlambdax:21.5+?x[0]?*?np.sin(4*?np.pi?*?x[0])?+?x[1]?*?np.sin(20*?np.pi?*?x[1])

pass

defmain(max_iter=500):

#?每次迭代得到的最優(yōu)解

optimalSolutions?=?[]

optimalValues?=?[]

#?決策變量的取值范圍

decisionVariables?=?[[-3.0,12.1],?[4.1,5.8]]

#?得到染色體編碼長度

lengthEncode?=?getEncodedLength(boundarylist=decisionVariables)

foriterationinrange(max_iter):

#?得到初始種群編碼

chromosomesEncoded?=?getIntialPopulation(lengthEncode,10)

#?種群解碼

decoded?=?decodedChromosome(lengthEncode,?chromosomesEncoded,?decisionVariables)

#?得到個體適應度值和個體的累積概率

evalvalues,?cum_proba?=?getFitnessValue(fitnessFunction(),?decoded)

#?選擇新的種群

newpopulations?=?selectNewPopulation(chromosomesEncoded,?cum_proba)

#?進行交叉操作

crossoverpopulation?=?crossover(newpopulations)

#?mutation

mutationpopulation?=?mutation(crossoverpopulation)

#?將變異后的種群解碼，得到每輪迭代最終的種群

final_decoded?=?decodedChromosome(lengthEncode,?mutationpopulation,?decisionVariables)

#?適應度評價

fitnessvalues,?cum_individual_proba?=?getFitnessValue(fitnessFunction(),?final_decoded)

#?搜索每次迭代的最優(yōu)解，以及最優(yōu)解對應的目標函數的取值

optimalValues.append(np.max(list(fitnessvalues)))

index?=?np.where(fitnessvalues?==?max(list(fitnessvalues)))

optimalSolutions.append(final_decoded[index[0][0],?:])

#?搜索最優(yōu)解

optimalValue?=?np.max(optimalValues)

optimalIndex?=?np.where(optimalValues?==?optimalValue)

optimalSolution?=?optimalSolutions[optimalIndex[0][0]]

returnoptimalSolution,?optimalValue

solution,?value?=?main()

print('最優(yōu)解:?x1,?x2')

print(solution[0],?solution[1])

print('最優(yōu)目標函數值:',?value)

#?測量運行時間

elapsedtime?=?timeit.timeit(stmt=main,?number=1)

print('Searching?Time?Elapsed:(S)',?elapsedtime)