模拟退火算法SAA求解常见问题
使用SAA进行数组排序、求解旅行商问题(TSP)、求解背包问题(Knapsack)、多项式回归(Polynomial Regression)。 本文仅仅是为了思考这些问题在SAA算法下如何求解,不考虑求解效率问题。
目录
SAA算法
此算法最小化目标值,如果需要最大化目标值,则需要修改一下公式(1)
Step1: 生成随机解 X, 带入目标函数求解目标值 V. 生成初始温度 T
Step2: 按照一个随机策略,由X生成新的解 X ′ X' X′, 并计算新的目标值 V ′ V' V′
Step3: 按照公式, 计算出接受新解的概率, P= { 1 , i f V ′ < V e − ( V ′ − V ) / T e l s e \begin{cases}1, \qquad if \quad V'< V \\ e^{-(V'-V) / T} \qquad else\end{cases} {1,ifV′<Ve−(V′−V)/Telse 式(1)
Step4: 生成随机数 pr,0 < pr < 1, 且如果 pr < P, 则接收新的解,令 X= X ′ X' X′, V= V ′ V' V′
Step5: 重复 Step2 - Step4 步骤,N次
Step6: 更新温度, T=T * a, 0< a < 1 或者其他合理的温度下降公式
Step7: 重复 Step2 到 Step6, 直到符合停止条件,例如: 温度小于一个给定的值,或者 V小于一个给定的值, 停止条件可以根据实际情况改变
SAA算法抽象类封装
import abc
import time
import numpy as np
import matplotlib.pyplot as plt
# plt 显示中文设置
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus']=False
class SAA(metaclass=abc.ABCMeta):
def __init__(self, x=None, max_epoch=float('inf'), t0=0, loop_times=30, rate=0.999):
self.x = x # np.ndarray
if x is not None:
self.value = self.objective_func(x)
else:
self.value = None
self.best_x = self.copy_x(x)
self.best_value = self.value
self.run_data = [] # 存储运行数据
self.max_epoch = max_epoch # 最大优化轮次
self.t = t0 # 初始温度
self.t0 = t0
self.epoch = 0
self.loop_times = loop_times
self.rate = rate
def probability(self, delta, t):
if t <= 0 or delta < 0:
return 1
return np.exp(-(delta / t))
def temperature_down(self, T, rate=0.999, num=0):
a = self.rate or rate
# 指数下降温度
return T * a
def temperature_init(self):
if self.t0:
self.t = self.t0
return self.t0
self.t = np.sum(np.array(self.x, 'float64')) ** 2 + 10000
self.t0 = self.t
# print("===init_t", self.t)
@staticmethod
def copy_x(x):
return x.copy()
@abc.abstractmethod
def update_x(self, x):
return x
@abc.abstractmethod
def objective_func(self, x):
return 0
def epoch_callback(self):
t = self.t
v = self.value
bv = self.best_value
self.run_data.append({"epoch": self.epoch, "temperature": t, "value": v, "best_value": bv})
return
def stop(self, min_temprature=0.000001):
n = self.epoch
if self.t <= min_temprature or n > self.max_epoch:
return True
return False
def show_run_curve(self):
x = [i['epoch'] for i in self.run_data]
y = [i['value'] for i in self.run_data]
y2 = [i['best_value'] for i in self.run_data]
plt.plot(x, y, label='value')
plt.plot(x, y2, label='best_value')
plt.title('优化曲线')
plt.legend()
plt.show()
def run(self, clear_run_data=True):
self.temperature_init()
self.epoch = 0
loop_times = self.loop_times
print(f"初始温度:{self.t}, loop_times: {loop_times}")
if clear_run_data:
self.run_data.clear()
if self.value is None:
self.value = self.objective_func(self.x)
self.best_value = self.value
best_value = self.best_value
while True:
for _ in range(loop_times):
x = self.x
# old_value = self.objective_func(x)
old_value = self.value
new_x = self.update_x(x)
new_value = self.objective_func(new_x)
delta = new_value - old_value
p = self.probability(delta, self.t)
pr = np.random.random()
if best_value > new_value: # 沿着value 下降的方向优化
self.best_x = self.copy_x(new_x)
self.best_value = new_value
best_value = new_value
if pr < p: # 按概率更新解
self.x = new_x
self.value = new_value
# print(f"==={n}概率接受解p: {p} T: {self.t} δ:{delta} ({pr},{p})")
print(f"==epoch: {self.epoch} t:{self.t} {old_value} -> {new_value}", end='\r')
# 回调
self.epoch_callback()
# 温度下降
self.t = self.temperature_down(self.t, num=self.epoch)
# 停机判断
if self.stop():
print(f"===epoch: {self.epoch} 退出:T:{self.t}, P:{p}, {old_value}")
break
self.epoch += 1
self.x = x
def run_with_t0_update(self, times=5, show_curve=True):
for i in range(times):
self.t0 *= 1 / (10 ** i)
self.run()
self.show_run_curve()
# show_path(tsp102.best_x, points_xy2)
排序
定义从小到大的排序目标函数为:左边大于右边的个数. 优化目标函数达到最小值.
import numpy as np
class SORT(SAA):
def update_x(self, x):
x = x.copy()
i0 = np.random.randint(0, x.shape[-1])
i1 = np.random.randint(0, x.shape[-1])
while x.shape[-1] > 2 and i0 == i1:
i1 = np.random.randint(0, x.shape[-1])
x[i0], x[i1] = x[i1], x[i0]
return x
def objective_func(self, x):
n = 0
for i in range(x.shape[-1] - 1):
n += int(x[i+1] < x[i])
time.sleep(0.00000001)
return n
def stop(self, min_temprature=0.000001):
# 排序问题的最优值是可知的, 可以提前退出计算
if self.best_value <= 0:
return True
n = self.epoch
if self.t <= min_temprature or n > self.max_epoch:
return True
return False
s = SORT(x=np.array([0, 1, 2, 3, 1,3233, 2, 29, 10, 889, 9]), t0=100)
# s = SORT(x=np.array([600000000, 5, 1, 2, 3, -4]))
# s = SORT(x=np.array([60000, 1, 2]))
s.run()
初始温度:100, loop_times: 30
===epoch: 5672 退出:T:0.3427793934809477, P:0.0029415406738396072, 0
s.show_run_curve()
print(list(s.x))
[0, 1, 1, 2, 2, 3, 9, 10, 29, 889, 3233]
print(list(s.best_x))
[0, 1, 1, 2, 2, 3, 9, 10, 29, 889, 3233]
TSP问题求解
TSP问题,即旅行商问题(Traveling Salesman Problem),是数学领域中的一个著名问题。
它是一个组合优化问题,旨在找到访问给定城市集合的最短可能路径,同时满足每个城市只能访问一次,并且最终要回到起点。
TSP问题可以用以下方式定义:
问题定义:给定一个包含n个城市的集合,找到访问每个城市一次并最终回到起点的最短路径。
目标:使得所走的路程最短。
限制:每个城市只能访问一次。
TSP问题是一个NP-完全问题(NPC),这意味着随着城市数量的增加,问题的计算复杂性会显著增加,使得在合理时间内找到精确解变得非常困难。
此外,TSP问题还有多个变体,包括对称TSP和非对称TSP,以及带有时间窗口和资源限制的TSP问题。
import numpy as np
class TSP(SAA):
def __init__(self, matrix, points=None, *args, **kwargs):
self.matrix = matrix
self.points = points
super().__init__(*args, **kwargs)
def update_x(self, x):
x = x.copy()
i0 = np.random.randint(1, x.shape[-1])
i1 = np.random.randint(1, x.shape[-1])
while x.shape[-1] > 3 and i0 == i1:
i1 = np.random.randint(1, x.shape[-1])
x[i0], x[i1] = x[i1], x[i0]
return x
def objective_func(self, x):
n = self.matrix[x[0], x[-1]]
for i in range(x.shape[-1] - 1):
n += self.matrix[x[i], x[i+1]]
return n
@staticmethod
def x_to_chars(x, points=''):
# points = self.points or points
return [points[i] for i in x]
@staticmethod
def generate_matrix(points='', distance_matrix=True):
d = len(points)
mat = np.random.random((d, d))
if distance_matrix:
# 使用numpy.triu生成上三角矩阵
mat += np.triu(mat, k=1)
mat += mat.T
# 填充对角线上的值
np.fill_diagonal(mat, 0)
return mat
return mat
@staticmethod
def generate_matrix_by_points(n=6, width=10):
x = np.random.random(n) * width
y = np.random.random(n) * width
plt.scatter(x, y)
plt.show()
points = list(zip(x, y))
mat = np.zeros((n, n), dtype='float32')
for i in range(n):
for j in range(i + 1, n):
d = ((points[i][0] - points[j][0]) ** 2 + (points[i][1] - points[j][1]) ** 2) ** 0.5
mat[i][j] = d
mat[j][i] = d
return mat, points
def show_run_curve(self):
plt.rcParams['font.sans-serif'] = ['SimHei']
x = [i['epoch'] for i in self.run_data]
y = [i['value'] for i in self.run_data]
y2 = [i['best_value'] for i in self.run_data]
plt.plot(x, y, label='value')
plt.plot(x, y2, label='best_value')
plt.title('TSP优化曲线')
plt.legend()
plt.show()
def run_with_t0_update(self, times=5, show_curve=True, points_xy=None):
for i in range(times):
self.t0 *= 1 / (10 ** i)
self.run()
self.show_run_curve()
if points_xy:
show_path(self.best_x, points_xy)
points = list('abcdefghjklmn')
# m = {(i,j): np.random.random() * 100 for i, _ in enumerate(points) for j, _ in enumerate(points)}
# print(m)
np.random.seed(200)
m = TSP.generate_matrix(points=points)
# print(m)
m, points_xy = TSP.generate_matrix_by_points(n = len(points))
tsp = TSP(m, points=points, x=np.array(range(len(points))))
tsp.run()
初始温度:16084.0, loop_times: 30
===epoch: 23489 退出:T:9.99338209184425e-07, P:0.0, 34.03287887573242417968750164
tsp.x_to_chars(tsp.x, points=tsp.points), tsp.x_to_chars(tsp.best_x, points=tsp.points)
(['a', 'b', 'e', 'j', 'g', 'm', 'n', 'h', 'c', 'k', 'l', 'd', 'f'],
['a', 'd', 'l', 'k', 'c', 'h', 'n', 'm', 'g', 'j', 'e', 'b', 'f'])
优化曲线
tsp.show_run_curve()
def show_path(path, points):
x = [i[0] for i in points]
y = [i[1] for i in points]
plt.scatter(x, y)
path_x = [points[i][0] for i in path]
path_y = [points[i][1] for i in path]
# 绘制箭头
s = 0.18
for i in range(0, len(path) - 1):
x_start, y_start = path_x[i], path_y[i]
x_end, y_end = path_x[i + 1], path_y[i + 1]
a = x_end - x_start
b = y_end - y_start
if a == 0 and b == 0:
continue
l = 1 - s / ((a ** 2 + b ** 2) ** 0.5)
a *= l
b *= l
plt.arrow(x_start, y_start,
a, b,
head_width=0.2, head_length=s, fc='red', ec='red')
# plt.plot(path_x, path_y, 'r--')
plt.show()
show_path(tsp.best_x, points_xy)
2-opt方法搜索
class TSPByKopt(TSP):
def update_x(self, x):
return self.opt2(x)
@staticmethod
def opt2(x):
x = x.copy()
i0 = np.random.randint(1, x.shape[-1])
i1 = np.random.randint(1, x.shape[-1])
while x.shape[-1] > 3 and i0 == i1:
i1 = np.random.randint(1, x.shape[-1])
while i0 < i1:
x[i0], x[i1] = x[i1], x[i0]
i0 += 1
i1 -= 1
return x
tsp2 = TSPByKopt(m, points=points,x=np.array(range(len(points))))
tsp2.run()
初始温度:16084.0
===epoch: 23489 退出:T:9.99338209184425e-07, P:1.0, 33.72549819946289498199462896
tsp2.x_to_chars(tsp2.x, points=tsp2.points), tsp2.x_to_chars(tsp2.best_x, points=tsp2.points)
(['a', 'd', 'l', 'k', 'c', 'h', 'n', 'm', 'g', 'j', 'e', 'b', 'f'],
['a', 'd', 'l', 'k', 'c', 'h', 'n', 'm', 'g', 'j', 'e', 'b', 'f'])
tsp2.show_run_curve()
show_path(tsp2.best_x, points_xy)
points = list(range(100))
m2, points_xy2 = TSP.generate_matrix_by_points(n = len(points))
tsp101 = TSP(m2, points=points,x=np.array(range(len(points))))
两点对换与2-opt方法比较
tsp101.run()
初始温度:24512500.0
===epoch: 30814 退出:T:9.997865508578958e-07, P:0.0, 145.0625610351562540283203125
tsp101.show_run_curve()
show_path(tsp101.best_x, points_xy2)
points = list(range(100))
# m2, points_xy2 = TSPByKopt.generate_matrix_by_points(n = len(points))
tsp102 = TSPByKopt(m2, points=points,x=np.array(range(len(points))))
tsp102.run()
初始温度:24512500.0
===epoch: 30814 退出:T:9.997865508578958e-07, P:1.0, 105.5768737792968873779296888
tsp102.show_run_curve()
show_path(tsp102.best_x, points_xy2)
tsp102.best_value, tsp102.best_value
(105.57687, 105.57687)
tsp102.run_with_t0_update()
初始温度:0.00245125
===epoch: 7800 退出:T:9.99450426296383e-07, P:0.0, 93.37339019775395860290527343
初始温度:0.000245125
===epoch: 5499 退出:T:9.990169333340238e-07, P:1.0, 92.5606842041015684204101566
初始温度:2.45125e-06
===epoch: 896 退出:T:9.991496610472398e-07, P:1.0, 92.5606842041015642041015675
求解背包问题
我们有n种物品,物品j的重量为wj,价格为pj
goods = 10 # 总物品数
np.random.seed(10)
# 需要求解的变量, 初始变量
selects = np.random.randint(0, 2, size=10)
# 每个物品可以被选择的次数
select_times = np.random.randint(1, 5, size=10)
weights = np.random.random(goods) # 物品重量
values = np.random.randint(1, 20, size=10) # 物品价值
total = 10 # 背包总容量
print(weights)
print(values)
print(selects)
print(f"select_times {select_times}")
[0.68535982 0.95339335 0.00394827 0.51219226 0.81262096 0.61252607
0.72175532 0.29187607 0.91777412 0.71457578]
[14 14 13 2 5 19 14 12 11 10]
[1 1 0 1 0 1 1 0 1 1]
select_times [1 2 2 3 1 2 1 3 1 3]
class Knapsack(SAA):
def __init__(self, weights, values, select_times, total, x=None, *args, **kwargs):
self.weights = weights
self.values=values
self.select_times = select_times
self.total = total
if x is not None:
self.constraint(x)
super().__init__(x=x, *args, **kwargs)
# def probability(self, delta, t):
# if t <= 0 or delta < 0:
# return 1
# return np.exp(-(delta / t))
def update_x(self, x):
x = x.copy()
self.constraint(x)
i = np.random.randint(0, len(x))
if np.random.random() > 0.5:
times = np.random.randint(0, self.select_times[i] + 1)
x[i] = 0
x[i] = max(min(int(self.total - x.sum()), times), 0)
else:
d = (np.random.random() > 0.5) * 2 - 1
a = min(self.select_times[i], x[i] + d)
x[i] = 0
x[i] = max(min(int(self.total - x.sum()), a), 0)
if x[i] > self.select_times[i]:
print("解异常")
self.constraint(x)
return x
def constraint(self, x):
"""
背包问题约束条件
"""
for i in range(len(x)):
while (x[i] > 0 and x.sum() > self.total) or (x[i] > self.select_times[i]):
x[i] = max(0, x[i] - 1)
if x[i] > self.select_times[i]:
raise Exception("约束异常")
def objective_func(self, x):
v = 0
for i in range(len(x)):
v += x[i] * self.values[i]
# 转换成求最小值
m = self.values.sum() * self.select_times.sum()
return m - v
# return -v
knapsack = Knapsack(weights, values, select_times, total, x=selects)
knapsack.run()
knapsack.show_run_curve()
初始温度:10049.0
===epoch: 23019 退出:T:9.992215880401829e-07, P:1.0, 2032
print(knapsack.best_x,knapsack.best_value)
print(knapsack.select_times)
[1 2 2 0 0 2 1 2 0 0] 2022
[1 2 2 3 1 2 1 3 1 3]
最短路径
s_points = list(range(50))
# shortest path problem
class SPP(TSP):
"""
邻接矩阵为无穷大表示不连通
"""
def update_x(self, x):
# assert x.shape[-1] > 3, "点数量过少"
x = x.copy()
self.constraint(x)
if np.random.random() > 0.6:
i0 = np.random.randint(1, x.shape[-1] - 1) # 起止点为0 和 -1, 不交换
i1 = np.random.randint(1, x.shape[-1] - 1)
while x.shape[-1] > 4 and i0 == i1:
i1 = np.random.randint(1, x.shape[-1] - 1)
if np.random.random() > 0.5:
x[i0], x[i1] = x[i1], x[i0]
else:
# 直接取消中间的路径, 如果不加此方法,则很难找到最优解
while i0 < i1:
x[i0 + 1] = x[i0]
i0 += 1
else:
i0 = np.random.randint(1, x.shape[-1] - 1)
s = np.random.randint(0, self.matrix.shape[-1] - 1)
x[i0] = s
self.constraint(x)
return x
def objective_func(self, x):
n = 0
inf = float('inf')
# print(self.matrix)
for i in range(x.shape[-1] - 1):
if self.matrix[x[i], x[i+1]] == inf:
n += np.abs(self.matrix).sum()
n += self.matrix[x[i], x[i+1]]
return n
def constraint(self, x):
"""
约束条件
"""
for i in range(x.shape[-1]- 2):
if self.matrix[x[i], x[i + 1]] == float('inf'):
x[i + 1] = x[i].copy()
np.random.seed(6)
m, p_list = SPP.generate_matrix_by_points(n=len(s_points))
x = np.array(range(len(s_points)))
np.random.shuffle(x)
max_xy = np.argwhere(m == m.max())[0]
print(max_xy)
x[0], x[-1] = max_xy
print(x)
# 起止点设置不连通
m[x[0], x[-1]] = float('inf')
shortest = SPP(matrix=m, points=s_points, x=x)
[31 32]
[31 10 24 38 28 33 11 4 22 13 46 17 12 37 40 34 21 36 7 27 9 26 1 18
48 6 47 0 25 30 29 20 14 2 43 41 16 31 45 23 32 44 49 42 35 39 15 8
3 32]
# shortest.t0 = 10 ** 3
shortest.run_with_t0_update(times=2, points_xy=p_list)
# shortest.show_run_curve()
初始温度:1607696.0
===epoch: 28091 退出:T:9.99767870003539e-07, P:0.0, 11.777696013450623181731224065
初始温度:160769.6
===epoch: 25790 退出:T:9.993342393559143e-07, P:1.0, 11.77769601345062396013450623
退火算法查找最短路径与Dijkstra算法比较
def dijkstra(matrix, i0=0):
n = matrix.shape[0]
find = [-1 for i in range(n)]
distance = [float("inf") for i in range(n)]
pre = [-1 for i in range(n)]
distance[i0] = 0
find[i0] = 1
for j in range(n):
if j != i0:
distance[j] = matrix[i0, j]
pre[j] = i0
# print(distance)
i = i0
while True:
# 遍历每一个未确认的点,尝试更新未确认点的距离,即:通过未确认的点,则距离变短
for j in range(n):
if find[j] > 0 or i == j:
continue
if matrix[i,j] >= float("inf"):
continue
new_d = matrix[i,j] + distance[i]
if new_d < distance[j]: # 当从i经过j点时,距离可以更短,则更新j的最小距离
distance[j] = new_d
pre[j] = i
# 在未确认的点中找到最小距离的点
m = float("inf")
i = None
for j in range(n):
if find[j] > 0:
continue
if distance[j] < m:
i = j
m = distance[j]
if i is None:
break
# print(f"next: {i}")
find[i] = 1
print('=' * 10)
print(find)
print(distance)
print(pre)
return distance, pre
def to_path(s, last=0):
p = [last]
while s[last] != -1:
last = s[last]
p.insert(0, last)
return p
distance, pre_points = dijkstra(m, i0=x[0])
dijkstra_path = to_path(pre_points, last=x[-1])
==========
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[7.304405, 7.0522614, 2.2776055, 10.450214, 9.338489, 8.000453, 4.8732204, 6.442265, 10.385517, 6.835383, 10.106222, 9.702004, 9.271525, 4.3813324, 6.7429028, 4.347516, 4.0250196, 9.10099, 7.93373, 6.641005, 9.448183, 7.1655736, 1.9282913, 2.9891455, 9.282758, 5.8116393, 8.759328, 6.4195356, 6.034359, 9.861004, 3.476944, 0, 11.777696, 3.2385242, 6.152703, 4.0818205, 6.212417, 10.283808, 9.299705, 8.050401, 2.5899367, 7.2155046, 6.103101, 4.4901824, 5.4372907, 6.782229, 8.695412, 10.481743, 7.812679, 3.9552972]
[31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, -1, 37, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31, 31]
print(dijkstra_path)
[31, 37, 32]
show_path(dijkstra_path, p_list)
多项式回归(Polynomial Regression)
def polynomial_func(x, params):
"""
params type: [(float, float)]
"""
v = 0
for a, b in params:
v += a * np.power(x, b)
return v
## 生成一个随机的多项式函数,测试代码
pa = [(np.random.random(), i) for i in range(10)]
pa = np.array(pa)
pa, pa.shape
(array([[0.17746591, 0. ],
[0.89390912, 1. ],
[0.32112076, 2. ],
[0.12548015, 3. ],
[0.34522942, 4. ],
[0.69496884, 5. ],
[0.02675415, 6. ],
[0.5419919 , 7. ],
[0.29978091, 8. ],
[0.66899116, 9. ]]),
(10, 2))
x = np.arange(-2, 2, 0.001)
y = np.apply_along_axis(lambda a: polynomial_func(a, pa), 0, x)
plt.plot(x, y)
plt.show()
拟合二次函数
np.random.seed(10)
train_x = np.arange(-2, 2, 0.01)
train_y = -np.power(train_x, 2) + np.random.random(len(train_x)) * 0.3
plt.scatter(train_x, train_y)
plt.show()
class PolynomialRegression(SAA):
def __init__(self, train_x, train_y, *args, **kwargs):
self.train_x = train_x
self.train_y = train_y
super().__init__(*args, **kwargs)
def update_x(self, x):
i = np.random.randint(0, x.shape[0])
if np.random.random() < 0.01:
x[i, 0] = np.random.random() * 2 - 1
else:
a = (np.random.random() * 2 - 1) * 0.1
x[i, 0] *= a
return x
def objective_func(self, x):
pre_y = np.apply_along_axis(lambda a: polynomial_func(a, x), 0, self.train_x)
d = pre_y - self.train_y
return (d ** 2).sum() / self.train_y.shape[0]
def predict(self, x, params=None):
if params is None:
params = self.best_x
return np.apply_along_axis(lambda a: polynomial_func(a, params), 0, x)
pol = PolynomialRegression(train_x, train_y, x=pa)
pol.run_with_t0_update(times=5)
初始温度:12410.3870039805, loop_times: 30
===epoch: 23230 退出:T:9.991767101361876e-07, P:0.0, 0.0171797724867069288576094586636
初始温度:1241.0387003980502, loop_times: 30
===epoch: 20928 退出:T:9.99743078972076e-07, P:0.0, 0.0169463164575696669857609458566
初始温度:12.410387003980503, loop_times: 30
===epoch: 16325 退出:T:9.99875903159779e-07, P:0.0, 0.0170044256091905779839097938762
初始温度:0.012410387003980504, loop_times: 30
===epoch: 9421 退出:T:9.99575009871613e-07, P:0.0, 0.01694526861586353778185395333466
初始温度:1.2410387003980505e-06, loop_times: 30
===epoch: 215 退出:T:9.998406312320972e-07, P:0.0, 0.01694526861586353757609458666
pre_y = pol.predict(train_x)
np.round(pol.best_x, 2)
array([[ 0. , 0. ],
[-0. , 1. ],
[-0.94, 2. ],
[ 0. , 3. ],
[ 0. , 4. ],
[-0. , 5. ],
[-0. , 6. ],
[-0. , 7. ],
[-0. , 8. ],
[ 0. , 9. ]])
plt.scatter(train_x, train_y, c='g')
plt.scatter(train_x, pre_y, c='r')
plt.show()
拟合三次函数
np.random.seed(10)
train_x = np.arange(-2, 2, 0.01)
train_y = np.power(train_x, 3) + np.random.random(len(train_x))
plt.scatter(train_x, train_y)
plt.show()
pol = PolynomialRegression(train_x, train_y, x=pa)
pol.run()
初始温度:12024.999999982396, loop_times: 30
===epoch: 23198 退出:T:9.996464766364758e-07, P:0.0, 0.20693877934235541598321377945
np.round(pol.best_x, 2)
array([[ 0.46, 0. ],
[-0. , 1. ],
[-0. , 2. ],
[ 0.88, 3. ],
[-0. , 4. ],
[ 0. , 5. ],
[-0. , 6. ],
[-0. , 7. ],
[ 0. , 8. ],
[ 0. , 9. ]])
pol.show_run_curve()
pre_y = pol.predict(train_x)
plt.scatter(train_x, train_y, c='g')
plt.scatter(train_x, pre_y, c='r')
plt.show()