Simple Genetic Alogrithm是模擬生物進化過程而提出的一種優化算法。SGA采用隨機導向搜索全局最優解或者說近似全局最優解。傳統的爬山算法(例如梯度下降，牛頓法)一次只優化一個解，并且對于多峰的目標函數很容易陷入局部最優解，而SGA算法一次優化一個種群(即一次優化多個解)，SGA比傳統的爬山算法更容易收斂到全局最優解或者近似全局最優解。
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。接下來，要確定是那個染色體中哪個位置的基因發生了變異。將330按照10進制序號進行編碼即從0,1,2,.......229。隨機從330個數中選擇330*Pm個數，假設其中一個數時154，chromosomeIndex = 154/33 =4,
geneIndex = 154%33 = 22。由此確定了第154號位置的基因位于第4個染色體的第22個位置上，將此位置的基因值置反完成基本位變異操作。
9、以上步驟完成了一次迭代的所有操作。接下就是評估的過程。對變異后得到的最終的種群進行解碼，利用解碼值求得每個個體的適應度值，將最大的適應度值保存下來，對應的解碼后的決策變量的值也保存下來。
10、根據迭代次數，假設是500次，重復執行1-9的步驟，最終得到是一個500個數值的最優適應度取值的數組以及一個500*n的決策變量取值數組（假設有n個決策變量）。從500個值中找到最優的一個（最大或者最小，根據定義的適應度函數來選擇）以及對應的決策變量的取值。
對于以上流程不是很清楚的地方，在代碼中有詳細的注釋。也可以自行查找資料補充理論。本文重點是實現
本代碼實現的問題是: 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
初始種群的編碼結果如下圖所示:
初始種群的解碼結果如下圖所示:
適應度值如圖所示:
輪盤賭選擇后的種群如圖所示；
單點交叉后的種群如圖所示:
基本位變異后的種群如圖所示；
最終結果如下圖所示；
源代碼如下；
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
|
# !/usr/bin/env python # -*- coding:utf-8 -*- # Author: wsw # 簡單實現SGA算法 import numpy as np from scipy.optimize import fsolve, basinhopping import random import timeit # 根據解的精度確定染色體(chromosome)的長度 # 需要根據決策變量的上下邊界來確定 def getEncodedLength(delta=0.0001, boundarylist=[]): # 每個變量的編碼長度 lengths = [] for i in boundarylist: lower = i[0] upper = i[1] # lamnda 代表匿名函數f(x)=0,50代表搜索的初始解 res = fsolve(lambda x: ((upper - lower) * 1 / delta) - 2 ** x - 1, 50) length = int(np.floor(res[0])) lengths.append(length) return lengths pass # 隨機生成初始編碼種群 def getIntialPopulation(encodelength, populationSize): # 隨機化初始種群為0 chromosomes = np.zeros((populationSize, sum(encodelength)), dtype=np.uint8) for i in range(populationSize): chromosomes[i, :] = np.random.randint(0, 2, sum(encodelength)) # print('chromosomes shape:', chromosomes.shape) return chromosomes # 染色體解碼得到表現型的解 def decodedChromosome(encodelength, chromosomes, boundarylist, delta=0.0001): populations = chromosomes.shape[0] variables = len(encodelength) decodedvalues = np.zeros((populations, variables)) for k, chromosome in enumerate(chromosomes): chromosome = chromosome.tolist() start = 0 for index, length in enumerate(encodelength): # 將一個染色體進行拆分，得到染色體片段 power = length - 1 # 解碼得到的10進制數字 demical = 0 for i in range(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 return decodedvalues # 得到個體的適應度值及每個個體被選擇的累積概率 def getFitnessValue(func, chromosomesdecoded): # 得到種群規模和決策變量的個數 population, nums = chromosomesdecoded.shape # 初始化種群的適應度值為0 fitnessvalues = np.zeros((population, 1)) # 計算適應度值 for i in range(population): fitnessvalues[i, 0] = func(chromosomesdecoded[i, :]) # 計算每個染色體被選擇的概率 probability = fitnessvalues / np.sum(fitnessvalues) # 得到每個染色體被選中的累積概率 cum_probability = np.cumsum(probability) return fitnessvalues, cum_probability # 新種群選擇 def selectNewPopulation(chromosomes, cum_probability): m, n = chromosomes.shape newpopulation = np.zeros((m, n), dtype=np.uint8) # 隨機產生M個概率值 randoms = np.random.rand(m) for i, randoma in enumerate(randoms): logical = cum_probability >= randoma index = np.where(logical == 1) # index是tuple,tuple中元素是ndarray newpopulation[i, :] = chromosomes[index[0][0], :] return newpopulation pass # 新種群交叉 def crossover(population, Pc=0.8): """ :param population: 新種群 :param Pc: 交叉概率默認是0.8 :return: 交叉后得到的新種群 """ # 根據交叉概率計算需要進行交叉的個體個數 m, n = population.shape numbers = np.uint8(m * Pc) # 確保進行交叉的染色體個數是偶數個 if numbers % 2 != 0: numbers += 1 # 交叉后得到的新種群 updatepopulation = np.zeros((m, n), dtype=np.uint8) # 產生隨機索引 index = random.sample(range(m), numbers) # 不進行交叉的染色體進行復制 for i in range(m): if not index.__contains__(i): updatepopulation[i, :] = population[i, :] # crossover while len(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:] return updatepopulation pass # 染色體變異 def mutation(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) # 確定每個將要變異的基因在整個染色體中的基因座(即基因的具體位置) for gene in mutationGeneIndex: # 確定變異基因位于第幾個染色體 chromosomeIndex = gene // n # 確定變異基因位于當前染色體的第幾個基因位 geneIndex = gene % n # mutation if updatepopulation[chromosomeIndex, geneIndex] == 0: updatepopulation[chromosomeIndex, geneIndex] = 1 else: updatepopulation[chromosomeIndex, geneIndex] = 0 return updatepopulation pass # 定義適應度函數 def fitnessFunction(): return lambda x: 21.5 + x[0] * np.sin(4 * np.pi * x[0]) + x[1] * np.sin(20 * np.pi * x[1]) pass def main(max_iter=500): # 每次迭代得到的最優解 optimalSolutions = [] optimalValues = [] # 決策變量的取值范圍 decisionVariables = [[-3.0, 12.1], [4.1, 5.8]] # 得到染色體編碼長度 lengthEncode = getEncodedLength(boundarylist=decisionVariables) for iteration in range(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) # 搜索每次迭代的最優解，以及最優解對應的目標函數的取值 optimalValues.append(np.max(list(fitnessvalues))) index = np.where(fitnessvalues == max(list(fitnessvalues))) optimalSolutions.append(final_decoded[index[0][0], :]) # 搜索最優解 optimalValue = np.max(optimalValues) optimalIndex = np.where(optimalValues == optimalValue) optimalSolution = optimalSolutions[optimalIndex[0][0]] return optimalSolution, optimalValue solution, value = main() print('最優解: x1, x2') print(solution[0], solution[1]) print('最優目標函數值:', value) # 測量運行時間 elapsedtime = timeit.timeit(stmt=main, number=1) print('Searching Time Elapsed:(S)', elapsedtime) |
以上就是本文的全部內容，希望對大家的學習有所幫助，也希望大家多多支持服務器之家。
原文鏈接：http://blog.csdn.net/qq_30666517/article/details/78637255
