使用 MSE 将二次函数拟合到数据
Fit quadratic function to data using MSE
所以我的想法是(从神经网络人员那里借来的)如果我有数据集 D,我可以通过首先计算误差相对于参数 (a, b,和 c),然后做一个最小化错误的小更新。我的问题是,下面的代码实际上并不符合曲线。对于线性的东西,类似的方法可行,但二次方似乎由于某种原因失败了。你能看出我做错了什么吗(假设或只是实现错误)
编辑:问题不够具体:以下代码不能很好地处理数据偏差。出于某种原因,它会以某种方式更新 a 和 b 参数,从而使 c 落在后面。这种方法类似于机器人学(使用雅可比矩阵寻找路径)和神经网络(基于错误寻找参数)所以它不是不合理的算法,现在的问题是,为什么这个具体实现没有产生我期望的结果。
在下面的Python代码中,我使用math作为m,MSE是一个计算两个数组之间的均方误差的函数。除此之外,代码是独立的
代码:
def quadraticRegression(data, dErr):
a = 1 #Starting values
b = 1
c = 1
a_momentum = 0 #Momentum to counter steady state error
b_momentum = 0
c_momentum = 0
estimate = [a*x**2 + b*x + c for x in range(len(data))] #Estimate curve
error = MSE(data, estimate) #Get errors 'n stuff
errorOld = 0
lr = 0.0000000001 #learning rate
while abs(error - errorOld) > dErr:
#Fit a (dE/da)
deda = sum([ 2*x**2 * (a*x**2 + b*x + c - data[x]) for x in range(len(data)) ])/len(data)
correction = deda*lr
a_momentum = (a_momentum)*0.99 + correction*0.1 #0.99 is to slow down momentum when correction speed changes
a = a - correction - a_momentum
#fit b (dE/db)
dedb = sum([ 2*x*(a*x**2 + b*x + c - data[x]) for x in range(len(data))])/len(data)
correction = dedb*lr
b_momentum = (b_momentum)*0.99 + correction*0.1
b = b - correction - b_momentum
#fit c (dE/dc)
dedc = sum([ 2*(a*x**2 + b*x + c - data[x]) for x in range(len(data))])/len(data)
correction = dedc*lr
c_momentum = (c_momentum)*0.99 + correction*0.1
c = c - correction - c_momentum
#Update model and find errors
estimate = [a*x**2 +b*x + c for x in range(len(data))]
errorOld = error
print(error)
error = MSE(data, estimate)
return a, b, c, error
对我来说,您的代码看起来完全正确!至少算法是正确的。我已将您的代码更改为使用 numpy
而不是纯 Python 进行快速计算。我还配置了一些参数,比如改变动量和学习率,也实现了 MSE
.
然后我用matplotlib
画了剧情动画。最后在动画上看起来你的回归实际上试图使曲线适合数据。虽然在最后一次迭代中拟合 sin(x)
for x
in [0; 2 * pi]
看起来像线性近似,但仍然尽可能接近数据点以求二次曲线。但是对于 sin(x)
for x
in [0; pi]
它看起来像是理想的近似值(它从大约 12-th
迭代开始拟合)。
i-th
帧动画只是用dErr = 0.7 ** (i + 15)
做回归。
我的动画对于 运行 脚本来说有点慢,但是如果你像这样添加参数 save
python script.py save
它会 render/save 到 line.gif
情节绘图动画。如果您 运行 没有参数的脚本,它将在您的 PC 屏幕上实时显示动画 plotting/fitting。
完整代码在图形之后,代码需要通过 运行ning 安装一些 python 模块一次 python -m pip install numpy matplotlib
。
接下来是 sin(x)
for x
in (0, pi)
:
接下来是 sin(x)
for x
in (0, 2 * pi)
:
接下来是 abs(x)
for x
in (-1, 1)
:
# Needs: python -m pip install numpy matplotlib
import math, sys
import numpy as np, matplotlib.pyplot as plt, matplotlib.animation as animation
from matplotlib.animation import FuncAnimation
x_range = (0., math.pi, 0.1) # (xmin, xmax, xstep)
y_range = (-0.2, 1.2) # (ymin, ymax)
num_iterations = 50
def f(x):
return np.sin(x)
def derr(iteration):
return 0.7 ** (iteration + 15)
def MSE(a, b):
return (np.abs(np.array(a) - np.array(b)) ** 2).mean()
def quadraticRegression(*, x, data, dErr):
x, data = np.array(x), np.array(data)
assert x.size == data.size, (x.size, data.size)
a = 1 #Starting values
b = 1
c = 1
a_momentum = 0.1 #Momentum to counter steady state error
b_momentum = 0.1
c_momentum = 0.1
estimate = a*x**2 + b*x + c #Estimate curve
error = MSE(data, estimate) #Get errors 'n stuff
errorOld = 0.
lr = 10. ** -4 #learning rate
while abs(error - errorOld) > dErr:
#Fit a (dE/da)
deda = np.sum(2*x**2 * (a*x**2 + b*x + c - data))/len(data)
correction = deda*lr
a_momentum = (a_momentum)*0.99 + correction*0.1 #0.99 is to slow down momentum when correction speed changes
a = a - correction - a_momentum
#fit b (dE/db)
dedb = np.sum(2*x*(a*x**2 + b*x + c - data))/len(data)
correction = dedb*lr
b_momentum = (b_momentum)*0.99 + correction*0.1
b = b - correction - b_momentum
#fit c (dE/dc)
dedc = np.sum(2*(a*x**2 + b*x + c - data))/len(data)
correction = dedc*lr
c_momentum = (c_momentum)*0.99 + correction*0.1
c = c - correction - c_momentum
#Update model and find errors
estimate = a*x**2 +b*x + c
errorOld = error
#print(error)
error = MSE(data, estimate)
return a, b, c, error
fig, ax = plt.subplots()
fig.set_tight_layout(True)
x = np.arange(x_range[0], x_range[1], x_range[2])
#ax.scatter(x, x + np.random.normal(0, 3.0, len(x)))
line0, line1 = None, None
do_save = len(sys.argv) > 1 and sys.argv[1] == 'save'
def g(x, derr):
a, b, c, error = quadraticRegression(x = x, data = f(x), dErr = derr)
return a * x ** 2 + b * x + c
def dummy(x):
return np.ones_like(x, dtype = np.float64) * 100.
def update(i):
global line0, line1
de = derr(i)
if line0 is None:
assert line1 is None
line0, = ax.plot(x, f(x), 'r-', linewidth=2)
line1, = ax.plot(x, g(x, de), 'r-', linewidth=2, color = 'blue')
ax.set_ylim(y_range[0], y_range[1])
if do_save:
sys.stdout.write(str(i) + ' ')
sys.stdout.flush()
label = 'iter {0} derr {1}'.format(i, round(de, math.ceil(-math.log(de) / math.log(10)) + 2))
line1.set_ydata(g(x, de))
ax.set_xlabel(label)
return line1, ax
if __name__ == '__main__':
anim = FuncAnimation(fig, update, frames = np.arange(0, num_iterations), interval = 200)
if do_save:
anim.save('line.gif', dpi = 200, writer = 'imagemagick')
else:
plt.show()
所以我的想法是(从神经网络人员那里借来的)如果我有数据集 D,我可以通过首先计算误差相对于参数 (a, b,和 c),然后做一个最小化错误的小更新。我的问题是,下面的代码实际上并不符合曲线。对于线性的东西,类似的方法可行,但二次方似乎由于某种原因失败了。你能看出我做错了什么吗(假设或只是实现错误)
编辑:问题不够具体:以下代码不能很好地处理数据偏差。出于某种原因,它会以某种方式更新 a 和 b 参数,从而使 c 落在后面。这种方法类似于机器人学(使用雅可比矩阵寻找路径)和神经网络(基于错误寻找参数)所以它不是不合理的算法,现在的问题是,为什么这个具体实现没有产生我期望的结果。
在下面的Python代码中,我使用math作为m,MSE是一个计算两个数组之间的均方误差的函数。除此之外,代码是独立的
代码:
def quadraticRegression(data, dErr):
a = 1 #Starting values
b = 1
c = 1
a_momentum = 0 #Momentum to counter steady state error
b_momentum = 0
c_momentum = 0
estimate = [a*x**2 + b*x + c for x in range(len(data))] #Estimate curve
error = MSE(data, estimate) #Get errors 'n stuff
errorOld = 0
lr = 0.0000000001 #learning rate
while abs(error - errorOld) > dErr:
#Fit a (dE/da)
deda = sum([ 2*x**2 * (a*x**2 + b*x + c - data[x]) for x in range(len(data)) ])/len(data)
correction = deda*lr
a_momentum = (a_momentum)*0.99 + correction*0.1 #0.99 is to slow down momentum when correction speed changes
a = a - correction - a_momentum
#fit b (dE/db)
dedb = sum([ 2*x*(a*x**2 + b*x + c - data[x]) for x in range(len(data))])/len(data)
correction = dedb*lr
b_momentum = (b_momentum)*0.99 + correction*0.1
b = b - correction - b_momentum
#fit c (dE/dc)
dedc = sum([ 2*(a*x**2 + b*x + c - data[x]) for x in range(len(data))])/len(data)
correction = dedc*lr
c_momentum = (c_momentum)*0.99 + correction*0.1
c = c - correction - c_momentum
#Update model and find errors
estimate = [a*x**2 +b*x + c for x in range(len(data))]
errorOld = error
print(error)
error = MSE(data, estimate)
return a, b, c, error
对我来说,您的代码看起来完全正确!至少算法是正确的。我已将您的代码更改为使用 numpy
而不是纯 Python 进行快速计算。我还配置了一些参数,比如改变动量和学习率,也实现了 MSE
.
然后我用matplotlib
画了剧情动画。最后在动画上看起来你的回归实际上试图使曲线适合数据。虽然在最后一次迭代中拟合 sin(x)
for x
in [0; 2 * pi]
看起来像线性近似,但仍然尽可能接近数据点以求二次曲线。但是对于 sin(x)
for x
in [0; pi]
它看起来像是理想的近似值(它从大约 12-th
迭代开始拟合)。
i-th
帧动画只是用dErr = 0.7 ** (i + 15)
做回归。
我的动画对于 运行 脚本来说有点慢,但是如果你像这样添加参数 save
python script.py save
它会 render/save 到 line.gif
情节绘图动画。如果您 运行 没有参数的脚本,它将在您的 PC 屏幕上实时显示动画 plotting/fitting。
完整代码在图形之后,代码需要通过 运行ning 安装一些 python 模块一次 python -m pip install numpy matplotlib
。
接下来是 sin(x)
for x
in (0, pi)
:
接下来是 sin(x)
for x
in (0, 2 * pi)
:
接下来是 abs(x)
for x
in (-1, 1)
:
# Needs: python -m pip install numpy matplotlib
import math, sys
import numpy as np, matplotlib.pyplot as plt, matplotlib.animation as animation
from matplotlib.animation import FuncAnimation
x_range = (0., math.pi, 0.1) # (xmin, xmax, xstep)
y_range = (-0.2, 1.2) # (ymin, ymax)
num_iterations = 50
def f(x):
return np.sin(x)
def derr(iteration):
return 0.7 ** (iteration + 15)
def MSE(a, b):
return (np.abs(np.array(a) - np.array(b)) ** 2).mean()
def quadraticRegression(*, x, data, dErr):
x, data = np.array(x), np.array(data)
assert x.size == data.size, (x.size, data.size)
a = 1 #Starting values
b = 1
c = 1
a_momentum = 0.1 #Momentum to counter steady state error
b_momentum = 0.1
c_momentum = 0.1
estimate = a*x**2 + b*x + c #Estimate curve
error = MSE(data, estimate) #Get errors 'n stuff
errorOld = 0.
lr = 10. ** -4 #learning rate
while abs(error - errorOld) > dErr:
#Fit a (dE/da)
deda = np.sum(2*x**2 * (a*x**2 + b*x + c - data))/len(data)
correction = deda*lr
a_momentum = (a_momentum)*0.99 + correction*0.1 #0.99 is to slow down momentum when correction speed changes
a = a - correction - a_momentum
#fit b (dE/db)
dedb = np.sum(2*x*(a*x**2 + b*x + c - data))/len(data)
correction = dedb*lr
b_momentum = (b_momentum)*0.99 + correction*0.1
b = b - correction - b_momentum
#fit c (dE/dc)
dedc = np.sum(2*(a*x**2 + b*x + c - data))/len(data)
correction = dedc*lr
c_momentum = (c_momentum)*0.99 + correction*0.1
c = c - correction - c_momentum
#Update model and find errors
estimate = a*x**2 +b*x + c
errorOld = error
#print(error)
error = MSE(data, estimate)
return a, b, c, error
fig, ax = plt.subplots()
fig.set_tight_layout(True)
x = np.arange(x_range[0], x_range[1], x_range[2])
#ax.scatter(x, x + np.random.normal(0, 3.0, len(x)))
line0, line1 = None, None
do_save = len(sys.argv) > 1 and sys.argv[1] == 'save'
def g(x, derr):
a, b, c, error = quadraticRegression(x = x, data = f(x), dErr = derr)
return a * x ** 2 + b * x + c
def dummy(x):
return np.ones_like(x, dtype = np.float64) * 100.
def update(i):
global line0, line1
de = derr(i)
if line0 is None:
assert line1 is None
line0, = ax.plot(x, f(x), 'r-', linewidth=2)
line1, = ax.plot(x, g(x, de), 'r-', linewidth=2, color = 'blue')
ax.set_ylim(y_range[0], y_range[1])
if do_save:
sys.stdout.write(str(i) + ' ')
sys.stdout.flush()
label = 'iter {0} derr {1}'.format(i, round(de, math.ceil(-math.log(de) / math.log(10)) + 2))
line1.set_ydata(g(x, de))
ax.set_xlabel(label)
return line1, ax
if __name__ == '__main__':
anim = FuncAnimation(fig, update, frames = np.arange(0, num_iterations), interval = 200)
if do_save:
anim.save('line.gif', dpi = 200, writer = 'imagemagick')
else:
plt.show()