python據說功能強大，觸角伸到各個領域，網上搜了一下其科學計算和工程計算能力也相當強，具備各種第三方包，除了性能軟肋外，其他無可指摘，甚至可以同matlab等專業工具一較高下。
從網上找了一個使用遺傳算法實現數據擬合的例子學習了一下，確實Python相當貼合自然語言，終於編程語言也能說人話了，代碼整體簡潔、優雅。。
代碼功能：給出一個隱藏函數 例如 z=x^2+y^2,生成200個數據，利用這200個數據，使用遺傳算法猜測這些數據是什么公式生成的。 （說的太直白，一點都不高大上）
代碼如下：
1
#
coding=utf-8
2 from random import random, randint, choice,uniform
3 from copy import deepcopy
4 import numpy as np
5 import matplotlib.pyplot as plt
6
7 from random import random, randint, choice
8 from copy import deepcopy
9 import numpy as np
10
11 # 運算類
12 class fwrapper:
13 def __init__(self, function, childcount, name):
14 self.function = function
15 self.childcount = childcount
16 self.name = name
17
18 # 節點類
19 class node:
20 def __init__(self, fw, children):
21 self.function = fw.function
22 self.name = fw.name
23 self.children = children
24 # 將inp指定的運算符作用到子節點上
25 def evaluate(self, inp):
26 # 循環調用子節點的子節點的子節點....的evaluate方法
27 results = [n.evaluate(inp) for n in self.children]
28 # 返回運算結果
29 return self.function(results)
30 # 打印本節點及所屬節點的操作運算符
31 def display(self, indent=0):
32 print( ' ' * indent) + self.name
33 for c in self.children:
34 c.display(indent + 1)
35
36 # 參數節點類，x+y 其中x，y都是參數節點
37 class paramnode:
38 def __init__(self, idx):
39 self.idx = idx
40 # evaluate方法返回paramnode節點值本身
41 def evaluate(self, inp):
42 return inp[self.idx]
43
44 def display(self, indent=0):
45 print ' %sp%d ' % ( ' ' * indent, self.idx)
46
47 # 常數節點
48 class constnode:
49 def __init__(self, v):
50 self.v = v
51
52 def evaluate(self, inp):
53 return self.v
54
55 def display(self, indent=0):
56 print ' %s%d ' % ( ' ' * indent, self.v)
57
58 # 操作運算符類
59 class opera:
60 # 使用前面定義的fwrapper類生產常用的加減乘除運算,第一個參數是本運算執行方式，第二個參數是本運算接受的參數個數，第三個參數是本運算名稱
61 addw = fwrapper( lambda l: l[0] + l[1], 2, ' add ')
62 subw = fwrapper( lambda l: l[0] - l[1], 2, ' subtract ')
63 mulw = fwrapper( lambda l: l[0] * l[1], 2, ' multiply ')
64
65 def iffunc(l):
66 if l[0] > 0:
67 return l[1]
68 else:
69 return l[2]
70 # 定義if運算
71 ifw = fwrapper(iffunc, 3, ' if ')
72
73 def isgreater(l):
74 if l[0] > l[1]:
75 return 1
76 else:
77 return 0
78 # 定義greater運算
79 gtw = fwrapper(isgreater, 2, ' isgreater ')
80 # 構建運算符集合
81 flist = [addw, mulw, ifw, gtw, subw]
82
83 # 使用node，paramnode，fwrapper構建一個example
84 def exampletree(self):
85 return node(self.ifw, [node(self.gtw, [paramnode(0), constnode(3)]), node(self.addw, [paramnode(1), constnode(5)]),
86 node(self.subw, [paramnode(1), constnode(2)]), ])
87
88
89 # 構建一顆隨機運算數，pc為參數（分叉）個數，maxdepth為樹的深度，fpr為運算符個數在運算符加節點總數中所占比例，ppr為參數個數在參數加常數個數總數中所占的比例
90 def makerandomtree(self,pc, maxdepth=4, fpr=0.5, ppr=0.6):
91 if random() < fpr and maxdepth > 0:
92 f = choice(self.flist)
93 # 遞歸調用makerandomtree實現子節點的創建
94 children = [self.makerandomtree(pc, maxdepth - 1, fpr, ppr) for i in range(f.childcount)]
95 return node(f, children)
96 elif random() < ppr:
97 return paramnode(randint(0, pc - 1))
98 else:
99 return constnode(randint(0, 10))
100
101
102 # 變異，變異概率probchange=0.1
103 def mutate(self,t, pc, probchange=0.1):
104 # 變異后返回一顆隨機樹
105 if random() < probchange:
106 return self.makerandomtree(pc)
107 else:
108 result = deepcopy(t)
109 # 遞歸調用，給其子節點變異的機會
110 if isinstance(t, node):
111 result.children = [self.mutate(c, pc, probchange) for c in t.children]
112 return result
113
114 # 交叉
115 def crossover(self,t1, t2, probswap=0.7, top=1):
116 # 如果符合交叉概率，就將t2的值返回，實現交叉；
117 if random() < probswap and not top:
118 return deepcopy(t2)
119 else:
120 # 如果本層節點未實現交配，遞歸詢問子節點是否符合交配條件
121 # 首先使用deepcopy保存本節點
122 result = deepcopy(t1)
123 if isinstance(t1, node) and isinstance(t2, node):
124 # 依次遞歸詢問t1下的各子孫節點交配情況，交配對象為t2的各子孫；t1,t2家族同輩交配
125 result.children = [self.crossover(c, choice(t2.children), probswap, 0) for c in t1.children]
126 return result
127
128 # random2.display()
129 # muttree=mutate(random2,2)
130 # muttree.display()
131 # cross=crossover(random1,random2)
132 # cross.display()
133
134 # 設置一個隱藏函數，也就是原始函數
135 def hiddenfunction(self,x, y):
136 return x ** 2+ y**2
137
138 # 依照隱藏函數，生成坐標數據
139 def buildhiddenset(self):
140 rows = []
141 for i in range(200):
142 x = randint(0, 10)
143 x=uniform(-1,1)
144 y = randint(0, 10)
145 y=uniform(-1,1)
146 rows.append([x, y, self.hiddenfunction(x, y)])
147 print( " rows: ",rows)
148 return rows
149
150 # 擬合成績函數，判定擬合函數（實際是一顆圖靈樹）貼近原始函數的程度
151 def scorefunction(self,tree, s):
152 dif = 0
153 # print("tree:",tree)
154 # print("s:",s)
155 for data in s:
156 # print("data[0]:",data[0])
157 # print("data[1]:",data[1])
158 v = tree.evaluate([data[0],data[1]])
159 # 累加每個數據的絕對值偏差
160 dif += abs(v - data[2])
161 return dif
162
163 # 返回一個成績判定函數rankfunction的句柄
164 def getrankfunction(self,dataset):
165 # 此函數調用擬合成績函數，並對成績排序，返回各個種群的成績
166 def rankfunction(population):
167 scores = [(self.scorefunction(t, dataset), t) for t in population]
168 scores.sort()
169 return scores
170 return rankfunction
171
172 # hiddenset=buildhiddenset()
173 # scorefunction(random2,hiddenset)
174 # scorefunction(random1,hiddenset)
175
176 def evolve(self,pc, popsize, rankfunction, maxgen=500, mutationrate=0.1, breedingrate=0.4, pexp=0.7, pnew=0.05):
177 # 輪盤算法
178 def selectindex():
179 return int(np.log(random()) / np.log(pexp))
180 # 使用隨機樹生成第一代各種群
181 population = [self.makerandomtree(pc) for i in range(popsize)]
182 # 計算每一代各種群的成績，
183 for i in range(maxgen):
184 scores = rankfunction(population)
185 # 打印歷代最好成績
186 print( ' the best score in genneration ',i, ' : ',scores[0][0])
187 # 成績完全吻合原函數的話，退出函數
188 if scores[0][0] == 0:
189 break
190 # 創建新一代各種群，成績前兩名的直接進入下一代
191 newpop = [scores[0][1], scores[1][1]]
192 while len(newpop) < popsize:
193 # pnew為引進隨機種群概率，未達此概率的，使用原種群的交配、變異生成新種群
194 if random() > pnew:
195 newpop.append(
196 self.mutate(self.crossover(scores[selectindex()][1], scores[selectindex()][1], probswap=breedingrate), pc,
197 probchange=mutationrate))
198 # 引入隨機種群
199 else:
200 newpop.append(self.makerandomtree(pc))
201 population = newpop
202 # 打印歷代最好種群
203 # scores[0][1].display()
204 return scores[0][1]
205
206
207
208 def main(argv):
209 e=opera()
210 def exampletree():
211 return node(e.ifw,[node(e.gtw,[paramnode(0),constnode(3)]),node(e.addw,[paramnode(1),constnode(5)]),node(e.subw,[paramnode(1),constnode(2)])])
212
213
214 # random1=e.makerandomtree(2)
215 # random1.evaluate([7,1])
216 # random1.evaluate([2,4])
217 # random2=e.makerandomtree(2)
218 # random2.evaluate([5,3])
219 # random2.evaluate([5,20])
220 # random1.display()
221 # random2.display()
222
223 # exampletree = e.exampletree()
224 # exampletree.display()
225 # print(exampletree.evaluate([6, 1]))
226 # print('after evaluate:')
227 # exampletree.display()
228 # exampletree.evaluate([2, 3])
229 # exampletree.evaluate([5, 3])
230 # exampletree.display()
231
232 a=opera()
233 row2=a.buildhiddenset()
234 # fig=plt.figure()
235 # ax=fig.add_subplot(1,1,1)
236 # plt.plot(np.random.randn(1000).cumsum())
237 # plt.show()
238
239
240
241 from mpl_toolkits.mplot3d import Axes3D
242 fig = plt.figure()
243 ax = fig.add_subplot(111, projection= ' 3d ')
244 X = [1, 1, 2, 2]
245 Y = [3, 4, 4, 3]
246 Z = [1, 2, 1, 1]
247 rx=[]
248 ry=[]
249 rz=[]
250 for i in row2:
251 rx.append(i[0])
252 ry.append(i[1])
253 rz.append(i[2])
254
255 ax.plot_trisurf(rx, ry, rz)
256 rz2=[]
257 rf = a.getrankfunction(row2)
258 final = a.evolve(2, 100, rf, mutationrate=0.2, breedingrate=0.1, pexp=0.7, pnew=0.1,maxgen=500)
259 print( ' __________________is it?_________________________ ')
260 final.display()
261 for j in range(0,len(rx)):
262 rz2.append(final.evaluate([rx[j],ry[j]]))
263 fig2 = plt.figure()
264 ax2 = fig2.add_subplot(111, projection= ' 3d ')
265 ax2.plot_trisurf(rx, ry, rz2)
266
267 plt.show()
268
269
270 # print(rf)
271 # final = a.evolve(2, 500, rf, mutationrate=0.2, breedingrate=0.1, pexp=0.7, pnew=0.1)
272 # print("final:",final)
273 # print(final.evaluate([1,8]))
274 # print(final.evaluate([2,9]))
275
276
277
278
279
280 if __name__== " __main__ ":
2 from random import random, randint, choice,uniform
3 from copy import deepcopy
4 import numpy as np
5 import matplotlib.pyplot as plt
6
7 from random import random, randint, choice
8 from copy import deepcopy
9 import numpy as np
10
11 # 運算類
12 class fwrapper:
13 def __init__(self, function, childcount, name):
14 self.function = function
15 self.childcount = childcount
16 self.name = name
17
18 # 節點類
19 class node:
20 def __init__(self, fw, children):
21 self.function = fw.function
22 self.name = fw.name
23 self.children = children
24 # 將inp指定的運算符作用到子節點上
25 def evaluate(self, inp):
26 # 循環調用子節點的子節點的子節點....的evaluate方法
27 results = [n.evaluate(inp) for n in self.children]
28 # 返回運算結果
29 return self.function(results)
30 # 打印本節點及所屬節點的操作運算符
31 def display(self, indent=0):
32 print( ' ' * indent) + self.name
33 for c in self.children:
34 c.display(indent + 1)
35
36 # 參數節點類，x+y 其中x，y都是參數節點
37 class paramnode:
38 def __init__(self, idx):
39 self.idx = idx
40 # evaluate方法返回paramnode節點值本身
41 def evaluate(self, inp):
42 return inp[self.idx]
43
44 def display(self, indent=0):
45 print ' %sp%d ' % ( ' ' * indent, self.idx)
46
47 # 常數節點
48 class constnode:
49 def __init__(self, v):
50 self.v = v
51
52 def evaluate(self, inp):
53 return self.v
54
55 def display(self, indent=0):
56 print ' %s%d ' % ( ' ' * indent, self.v)
57
58 # 操作運算符類
59 class opera:
60 # 使用前面定義的fwrapper類生產常用的加減乘除運算,第一個參數是本運算執行方式，第二個參數是本運算接受的參數個數，第三個參數是本運算名稱
61 addw = fwrapper( lambda l: l[0] + l[1], 2, ' add ')
62 subw = fwrapper( lambda l: l[0] - l[1], 2, ' subtract ')
63 mulw = fwrapper( lambda l: l[0] * l[1], 2, ' multiply ')
64
65 def iffunc(l):
66 if l[0] > 0:
67 return l[1]
68 else:
69 return l[2]
70 # 定義if運算
71 ifw = fwrapper(iffunc, 3, ' if ')
72
73 def isgreater(l):
74 if l[0] > l[1]:
75 return 1
76 else:
77 return 0
78 # 定義greater運算
79 gtw = fwrapper(isgreater, 2, ' isgreater ')
80 # 構建運算符集合
81 flist = [addw, mulw, ifw, gtw, subw]
82
83 # 使用node，paramnode，fwrapper構建一個example
84 def exampletree(self):
85 return node(self.ifw, [node(self.gtw, [paramnode(0), constnode(3)]), node(self.addw, [paramnode(1), constnode(5)]),
86 node(self.subw, [paramnode(1), constnode(2)]), ])
87
88
89 # 構建一顆隨機運算數，pc為參數（分叉）個數，maxdepth為樹的深度，fpr為運算符個數在運算符加節點總數中所占比例，ppr為參數個數在參數加常數個數總數中所占的比例
90 def makerandomtree(self,pc, maxdepth=4, fpr=0.5, ppr=0.6):
91 if random() < fpr and maxdepth > 0:
92 f = choice(self.flist)
93 # 遞歸調用makerandomtree實現子節點的創建
94 children = [self.makerandomtree(pc, maxdepth - 1, fpr, ppr) for i in range(f.childcount)]
95 return node(f, children)
96 elif random() < ppr:
97 return paramnode(randint(0, pc - 1))
98 else:
99 return constnode(randint(0, 10))
100
101
102 # 變異，變異概率probchange=0.1
103 def mutate(self,t, pc, probchange=0.1):
104 # 變異后返回一顆隨機樹
105 if random() < probchange:
106 return self.makerandomtree(pc)
107 else:
108 result = deepcopy(t)
109 # 遞歸調用，給其子節點變異的機會
110 if isinstance(t, node):
111 result.children = [self.mutate(c, pc, probchange) for c in t.children]
112 return result
113
114 # 交叉
115 def crossover(self,t1, t2, probswap=0.7, top=1):
116 # 如果符合交叉概率，就將t2的值返回，實現交叉；
117 if random() < probswap and not top:
118 return deepcopy(t2)
119 else:
120 # 如果本層節點未實現交配，遞歸詢問子節點是否符合交配條件
121 # 首先使用deepcopy保存本節點
122 result = deepcopy(t1)
123 if isinstance(t1, node) and isinstance(t2, node):
124 # 依次遞歸詢問t1下的各子孫節點交配情況，交配對象為t2的各子孫；t1,t2家族同輩交配
125 result.children = [self.crossover(c, choice(t2.children), probswap, 0) for c in t1.children]
126 return result
127
128 # random2.display()
129 # muttree=mutate(random2,2)
130 # muttree.display()
131 # cross=crossover(random1,random2)
132 # cross.display()
133
134 # 設置一個隱藏函數，也就是原始函數
135 def hiddenfunction(self,x, y):
136 return x ** 2+ y**2
137
138 # 依照隱藏函數，生成坐標數據
139 def buildhiddenset(self):
140 rows = []
141 for i in range(200):
142 x = randint(0, 10)
143 x=uniform(-1,1)
144 y = randint(0, 10)
145 y=uniform(-1,1)
146 rows.append([x, y, self.hiddenfunction(x, y)])
147 print( " rows: ",rows)
148 return rows
149
150 # 擬合成績函數，判定擬合函數（實際是一顆圖靈樹）貼近原始函數的程度
151 def scorefunction(self,tree, s):
152 dif = 0
153 # print("tree:",tree)
154 # print("s:",s)
155 for data in s:
156 # print("data[0]:",data[0])
157 # print("data[1]:",data[1])
158 v = tree.evaluate([data[0],data[1]])
159 # 累加每個數據的絕對值偏差
160 dif += abs(v - data[2])
161 return dif
162
163 # 返回一個成績判定函數rankfunction的句柄
164 def getrankfunction(self,dataset):
165 # 此函數調用擬合成績函數，並對成績排序，返回各個種群的成績
166 def rankfunction(population):
167 scores = [(self.scorefunction(t, dataset), t) for t in population]
168 scores.sort()
169 return scores
170 return rankfunction
171
172 # hiddenset=buildhiddenset()
173 # scorefunction(random2,hiddenset)
174 # scorefunction(random1,hiddenset)
175
176 def evolve(self,pc, popsize, rankfunction, maxgen=500, mutationrate=0.1, breedingrate=0.4, pexp=0.7, pnew=0.05):
177 # 輪盤算法
178 def selectindex():
179 return int(np.log(random()) / np.log(pexp))
180 # 使用隨機樹生成第一代各種群
181 population = [self.makerandomtree(pc) for i in range(popsize)]
182 # 計算每一代各種群的成績，
183 for i in range(maxgen):
184 scores = rankfunction(population)
185 # 打印歷代最好成績
186 print( ' the best score in genneration ',i, ' : ',scores[0][0])
187 # 成績完全吻合原函數的話，退出函數
188 if scores[0][0] == 0:
189 break
190 # 創建新一代各種群，成績前兩名的直接進入下一代
191 newpop = [scores[0][1], scores[1][1]]
192 while len(newpop) < popsize:
193 # pnew為引進隨機種群概率，未達此概率的，使用原種群的交配、變異生成新種群
194 if random() > pnew:
195 newpop.append(
196 self.mutate(self.crossover(scores[selectindex()][1], scores[selectindex()][1], probswap=breedingrate), pc,
197 probchange=mutationrate))
198 # 引入隨機種群
199 else:
200 newpop.append(self.makerandomtree(pc))
201 population = newpop
202 # 打印歷代最好種群
203 # scores[0][1].display()
204 return scores[0][1]
205
206
207
208 def main(argv):
209 e=opera()
210 def exampletree():
211 return node(e.ifw,[node(e.gtw,[paramnode(0),constnode(3)]),node(e.addw,[paramnode(1),constnode(5)]),node(e.subw,[paramnode(1),constnode(2)])])
212
213
214 # random1=e.makerandomtree(2)
215 # random1.evaluate([7,1])
216 # random1.evaluate([2,4])
217 # random2=e.makerandomtree(2)
218 # random2.evaluate([5,3])
219 # random2.evaluate([5,20])
220 # random1.display()
221 # random2.display()
222
223 # exampletree = e.exampletree()
224 # exampletree.display()
225 # print(exampletree.evaluate([6, 1]))
226 # print('after evaluate:')
227 # exampletree.display()
228 # exampletree.evaluate([2, 3])
229 # exampletree.evaluate([5, 3])
230 # exampletree.display()
231
232 a=opera()
233 row2=a.buildhiddenset()
234 # fig=plt.figure()
235 # ax=fig.add_subplot(1,1,1)
236 # plt.plot(np.random.randn(1000).cumsum())
237 # plt.show()
238
239
240
241 from mpl_toolkits.mplot3d import Axes3D
242 fig = plt.figure()
243 ax = fig.add_subplot(111, projection= ' 3d ')
244 X = [1, 1, 2, 2]
245 Y = [3, 4, 4, 3]
246 Z = [1, 2, 1, 1]
247 rx=[]
248 ry=[]
249 rz=[]
250 for i in row2:
251 rx.append(i[0])
252 ry.append(i[1])
253 rz.append(i[2])
254
255 ax.plot_trisurf(rx, ry, rz)
256 rz2=[]
257 rf = a.getrankfunction(row2)
258 final = a.evolve(2, 100, rf, mutationrate=0.2, breedingrate=0.1, pexp=0.7, pnew=0.1,maxgen=500)
259 print( ' __________________is it?_________________________ ')
260 final.display()
261 for j in range(0,len(rx)):
262 rz2.append(final.evaluate([rx[j],ry[j]]))
263 fig2 = plt.figure()
264 ax2 = fig2.add_subplot(111, projection= ' 3d ')
265 ax2.plot_trisurf(rx, ry, rz2)
266
267 plt.show()
268
269
270 # print(rf)
271 # final = a.evolve(2, 500, rf, mutationrate=0.2, breedingrate=0.1, pexp=0.7, pnew=0.1)
272 # print("final:",final)
273 # print(final.evaluate([1,8]))
274 # print(final.evaluate([2,9]))
275
276
277
278
279
280 if __name__== " __main__ ":
281 main(0)
看懂不一定寫的出來，這是這次寫這個程序最大的體會， 得定時拿出來復習復習。