神經網絡BP-GD算法和自動控制PID算法有類似之處，都是利用誤差回報對問題進行求解，不同的是自動控制調節的是系統的輸入，神經網絡調節的是系統本身。本文将引入誤差微分控制超參數kd_damp對NAGD算法進行優化，收斂的速度更快！波動更小！
自動控制PID算法與神經網絡Momentum算法比較：
http://www.sohu.com/a/242354509_297288
http://www.360doc.com/content/16/1010/08/36492363_597225745.shtml
https://blog.csdn.net/tsyccnh/article/details/76673073
NAGD算法在Momentum上考慮了誤差微分項，和PID算法更為接近。
Momentum法公式：
注意，上式中alpha移到梯度函數g()前，結果是一樣的，網上存在這兩種形式的公式！
NAGD公式一：
NAGD公式二：
注意，上式中alpha移到梯度函數g()前并去掉g()中的alpha，就是NAGD公式一的形式，結果是一樣的！
NAGD公式三：
all_loss = []
all_step = []
last_a = a
last_b = b
va = 0
vb = 0
gamma = 0.9
####增加誤差微分項控制,當kd_damp=0誤差微分項不起作用，相當于Momentum:：
kd_damp=0.0
####
for step in range(1,100):
loss = 0
all_da = 0
all_db = 0
a_ahead = a - gamma*va*kd_damp
b_ahead = b - gamma*vb*kd_damp
#-- 求loss
for i in range(0,len(x)):
y_p = a_ahead*x[i] + b_ahead
loss = loss + (y[i] - y_p)*(y[i] - y_p)/2
all_da = all_da + da(y[i],y_p,x[i])
all_db = all_db + db(y[i],y_p)
loss = loss/len(x)
……….
last_a = a
last_b = b
###
#-- 參數更新
# print('a = %.3f,b = %.3f' % (a,b))
va = gamma * va+ rate*all_da
vb = gamma * vb+ rate*all_db
a = a - va
b = b - vb
#--
對NAGD引入新的超參數kd_damp用于控制誤差微分項，當采用标準的NAGD算法時：kd_damp=1.0，結果如下圖：
當kd_damp=0誤差微分項不起作用，相當于Momentum：
當kd_damp=2.0增強誤差微分項的作用，優化收斂的速度更快，波動更小！
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# 本代碼是一個最簡單的線形回歸問題，優化函數為 momentum with Nesterov
rate = 0.01 # learning rate
def da(y,y_p,x):
return (y-y_p)*(-x)
def db(y,y_p):
return (y-y_p)*(-1)
def calc_loss(a,b,x,y):
tmp = y - (a * x + b)
tmp = tmp ** 2 # 對矩陣内的每一個元素平方
SSE = sum(tmp) / (2 * len(x))
return SSE
def draw_hill(x,y):
a = np.linspace(-20,20,100)
print(a)
b = np.linspace(-20,20,100)
x = np.array(x)
y = np.array(y)
allSSE = np.zeros(shape=(len(a),len(b)))
for ai in range(0,len(a)):
for bi in range(0,len(b)):
a0 = a[ai]
b0 = b[bi]
SSE = calc_loss(a=a0,b=b0,x=x,y=y)
allSSE[ai][bi] = SSE
a,b = np.meshgrid(a, b)
return [a,b,allSSE]
# 模拟資料
x = [30 ,35,37, 59, 70, 76, 88, 100]
y = [1100, 1423, 1377, 1800, 2304, 2588, 3495, 4839]
# 資料歸一化
x_max = max(x)
x_min = min(x)
y_max = max(y)
y_min = min(y)
for i in range(0,len(x)):
x[i] = (x[i] - x_min)/(x_max - x_min)
y[i] = (y[i] - y_min)/(y_max - y_min)
[ha,hb,hallSSE] = draw_hill(x,y)
hallSSE = hallSSE.T# 重要，将所有的losses做一個轉置。原因是矩陣是以左上角至右下角順序排列元素，而繪圖是以左下角為原點。
# 初始化a,b值
a = 10.0
b = -20.0
fig = plt.figure(1, figsize=(12, 8))
fig.suptitle('learning rate: %.2f method: Nesterov momentum'%(rate), fontsize=15)
# 繪制圖1的曲面
ax = fig.add_subplot(2, 2, 1, projection='3d')
ax.set_top_view()
ax.plot_surface(ha, hb, hallSSE, rstride=2, cstride=2, cmap='rainbow')
# 繪制圖2的等高線圖
plt.subplot(2,2,2)
ta = np.linspace(-20, 20, 100)
tb = np.linspace(-20, 20, 100)
plt.contourf(ha,hb,hallSSE,15,alpha=0.5,cmap=plt.cm.hot)
C = plt.contour(ha,hb,hallSSE,15,colors='black')
plt.clabel(C,inline=True)
plt.xlabel('a')
plt.ylabel('b')
plt.ion() # iteration on
all_loss = []
all_step = []
last_a = a
last_b = b
va = 0
vb = 0
gamma = 0.9
####增加誤差微分項控制,當kd_damp=0誤差微分項不起作用，相當于Momentum:：
kd_damp=2.0
####
for step in range(1,100):
loss = 0
all_da = 0
all_db = 0
a_ahead = a - gamma*va*kd_damp
b_ahead = b - gamma*vb*kd_damp
#-- 求loss
for i in range(0,len(x)):
y_p = a_ahead*x[i] + b_ahead
loss = loss + (y[i] - y_p)*(y[i] - y_p)/2
all_da = all_da + da(y[i],y_p,x[i])
all_db = all_db + db(y[i],y_p)
loss = loss/len(x)
### 繪圖區
# 繪制圖1中的loss點
ax.scatter(a, b, loss, color='black')
# 繪制圖2中的loss點
plt.subplot(2,2,2)
plt.scatter(a,b,s=5,color='blue')
plt.plot([last_a,a],[last_b,b],color='aqua')
# 繪制圖3中的回歸直線
plt.subplot(2, 2, 3)
plt.plot(x, y)
plt.plot(x, y, 'o')
x_ = np.linspace(0, 1, 2)
y_draw = a * x_ + b
plt.plot(x_, y_draw)
# 繪制圖4的loss更新曲線
all_loss.append(loss)
all_step.append(step)
plt.subplot(2,2,4)
plt.plot(all_step,all_loss,color='orange')
plt.xlabel("step")
plt.ylabel("loss")
last_a = a
last_b = b
###
#-- 參數更新
# print('a = %.3f,b = %.3f' % (a,b))
va = gamma * va+ rate*all_da
vb = gamma * vb+ rate*all_db
a = a - va
b = b - vb
#--
if step%1 == 0:
print("step: ", step, " loss: ", loss)
plt.show()
plt.pause(0.01)
plt.show()
plt.pause(99999999999)