模拟退火算法(Simulated Annealing)

模拟退火算法

说说什么是算法？有人理解是输入数据经过一些步骤的处理，最后输出。在我理解来看就是一种优化问题或者说是数值分析，无论机器学习中的算法还是深度学习的算法，寻找最优模型都是解最小损失的函数。最后计算得到的模型都是对常用的模型的套用，得到具体模型前，关键还是求解损失函数，具体化模型参数。模拟退火不是模型，是一种算法，用于求解最优值。

背景

摘自维基百科模拟退火： 模拟退火来自冶金学的专有名词退火。退火是将材料加热后再经特定速率冷却，目的是增大晶粒的体积，并且减少晶格中的缺陷。材料中的原子原来会停留在使内能有局部最小值的位置，加热使能量变大，原子会离开原来位置，而随机在其他位置中移动。退火冷却时速度较慢，使得原子有较多可能可以找到内能比原先更低的位置。 模拟退火的原理也和金属退火的原理近似：我们将热力学的理论套用到统计学上，将搜寻空间内每一点想像成空气内的分子；分子的能量，就是它本身的动能；而搜寻空间内的每一点，也像空气分子一样带有“能量”，以表示该点对命题的合适程度。算法先以搜寻空间内一个任意点作起始：每一步先选择一个“邻居”，然后再计算从现有位置到达“邻居”的概率。 可以证明，模拟退火算法所得解依概率收敛到全局最优解。
不知到你什么感受，反正我开始看这段内容，然后理解算法有十万八千里的感觉。一句话，没有毛用。

公式

先介绍这个式子，$T$表示温度，$k$玻尔兹曼常量，$\Delta E$表示变化前后的能量的差值

上面这个公式$T$是温度，退火表示温度的变低的过程，所以温度应用随时间变化。于是，$T$从某一高温状态T向低温状态冷却时的降温函数

$T$₀表示初始温度，$t$表示时间。言归正传，将上面的公式利用起来总结为以下概率公式，$dE$表示新值减旧值得到的状态差，假定温度$T$短时间不变化（不重要）

摆脱掉温度的干扰，找最优解的公式中，上面公式理解为：若函数上的y值变化为新的y值后变小了，则以概率为1的可能性，朝这个方向走，反之，可以一种随机概率其他方向走，就好比周志华《机器学习》中摇臂机问题要以$ϵ$的概率探索一样向其他方向走。
特别说明：重要的是思想，公式并不重要，上面的概率函数完全可以改用其他函数，温度衰减过程也可以用其他方式代替。

代码流程

代码由于是采用的是2维函数，方向单一所以比较简单：

import numpy as np
import matplotlib.pyplot as plt
#1、假设一个多项式为误差函数，需要找最优点
def Y(x):
    return (x - 2) * (x - 4) *(x+7)* (x + 11) * (x - 9)
    
#2、[-12,10]选随机选一个点做起始点
X=np.random.uniform(-12,10)
#-------下面注释由于最后跟温度没有关系，只需要最终的点，所以不需要两个循环
T=1000#表示初始温度
#Tend=1#表示最后温度，温度越低降温速率越低
# while T>Tend:
#     for i in range(20):
#         print(X)
#         oldY=Y(X)
#         preX=X+np.random.uniform(-1,1)
#         if preX>-12 and preX<10:
#           if Y(preX)<oldY:
#                X=preX
#           else:
#                if np.random.rand()>np.exp(-(Y(preX)-oldY)/T):
#                   X=preX
    # t+=1
    # T=1000/t
#--------------------------------------------------------------------------------
for i in range(200):
    oldY=Y(X)#记录旧点
    preX=X+np.random.uniform(-1,1)#随机选择邻域
    if preX>-12 and preX<10:#验证可能的新点在区间
      if Y(preX)<oldY:#按Metropolis准则
           X=preX
      else:
           if np.random.rand()>np.exp(-(Y(preX)-oldY)/T):#这个过程表示以一定的概率选择其他方向
              X=preX
              T+=0.95*T#衰减不一定按上面的函数来
print(X)
print(Y(X))
#!!!完
#下面是绘图,看效果，明显可以看到最优值,不是每次都能选到最优值
x=np.linspace(-12,10,100)
y=[Y(x) for x in x]
#图的操作
plt.plot(x,y,color='b')
gc=plt.gca()
gc.spines['top'].set_color(None)
gc.spines['right'].set_color(None)
gc.spines['left'].set_position(('data',0))
gc.spines['bottom'].set_position(('data',0))
plt.plot(X,Y(X),'or')
plt.show()

注意：不是每次都能选到最优值，所以需要取多次结果，取最小值，找最大值，把条件反过来写即可

[1]:https://www.cnblogs.com/ECJTUACM-873284962/p/8468831.html
[2]:https://blog.csdn.net/weixin_40562999/article/details/80853418
[3]:https://www.imooc.com/article/details/id/30160