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__
“:
281 main(0)
看懂不一定写的出来,这是这次写这个程序最大的体会, 得定时拿出来复习复习。