最小二乘法拟合正弦函数

y=sin2πxy=sin2\pi xy=sin2πx

y+ϵy+\epsilony+ϵ, ϵ\epsilonϵ~N(0,1)

#给y加上高斯分布的噪点

(X,Y)(X,Y)(X,Y)

attempt to fitting the function : y=a0+a1x+a2x2+..+anxny=a_0+a_1x+a_2x^2+..+a_nx^ny=a0​+a1​x+a2​x2+..+an​xn

#导入包import numpy as npimport scipy as spfrom scipy.optimize import leastsqimport matplotlib.pyplot as plt

#定义正弦函数def sine_func(x):return np.sin(2*np.pi*x)#定义多项式函数def fit_func(p,x): #p为数组f=np.poly1d(p) #poly1d函数能根据输入的数组p生成一个多项式return f(x)

np.poly1d([1,2,3])

1∗x2+2∗x1+3∗x01*x^2+2*x^1+3*x^01∗x2+2∗x1+3∗x0

policy: min⁡a(y−(ax+bx0))2\min_a (y-(ax+bx^0))^2mina​(y−(ax+bx0))2

#评价差 errordef error_func(p,x,y):error=fit_func(p,x)-y # error=将x值带入多项式后减去y值return error #返回真实值与我们拟合的曲线上对应的值的差x=np.linspace(0,1,10) #linspace函数产生0-1之间的步长为0.1的数组y_=sine_func(x) y=[np.random.normal(0,0.1)+y1 for y1 in y_] #得到加上噪点后的y值print(y)

[0.018863818084578765, 0.643352435292383, 0.9995270889071429, 0.8859076268379011, 0.3123492536633993, -0.2445361290110839, -0.857851035405584, -0.7992334756421094, -0.6318390256462963, -0.15169614809124077]

x_points=np.linspace(0,1,1000)#leastsqura..def fitting(M=0): #得到拟合M次的曲线p_init=np.random.rand(M+1) #生成M+1个服从[0,1)均匀分布的随机样本值p_lsq=leastsq(error_func,p_init,args=(x, y)) # 最小二乘法函数 三个参数：误差函数、函数参数列表、数据点print('parameters:', p_lsq[0]) #p_lsp函数返回值的第一个为拟合曲线的参数plt.plot(x_points,sine_func(x_points),label='real') # real曲线plt.plot(x_points,fit_func(p_lsq[0],x_points),label='fitting') #拟合曲线 plt.plot(x,y,'bo', label = 'noise') #训练值样本点 'b0'画出蓝色的圆圈plt.legend()return p_lsqfitting(M=0)

parameters: [0.01748444](array([0.01748444]), 3)

fitting(M=1)

parameters: [-1.37690968 0.70593928](array([-1.37690968, 0.70593928]), 1)

fitting(M=5)

parameters: [-5.67180721e+01 1.37260155e+02 -9.56060001e+01 9.15613970e+005.73126111e+00 1.54421206e-02](array([-5.67180721e+01, 1.37260155e+02, -9.56060001e+01, 9.15613970e+00,5.73126111e+00, 1.54421206e-02]),1)

fitting(M=9)

parameters: [ 2.44048348e+04 -1.07709860e+05 1.99353327e+05 -2.00963213e+051.19785217e+05 -4.27681566e+04 8.80453207e+03 -9.53138175e+024.62865992e+01 1.88638181e-02](array([ 2.44048348e+04, -1.07709860e+05, 1.99353327e+05, -2.00963213e+05,1.19785217e+05, -4.27681566e+04, 8.80453207e+03, -9.53138175e+02,4.62865992e+01, 1.88638181e-02]),2)

policy: min⁡a(y−(ax+bx0))2+λ∥a∥2\min_a (y-(ax+bx^0))^2+\lambda \Vert a \Vert_2mina​(y−(ax+bx0))2+λ∥a∥2​

penalty parameter: λ\lambdaλ

regularization: ∥a∥2\Vert a \Vert_2∥a∥2​

lam=0.001def error_plus_regu(p,x,y):error=y-fit_func(p,x)error=np.append(error,np.sqrt(0.5*lam*np.square(p))) #sqrt函数 求算数平方根 square函数 平方return error

M=3p_init=np.random.rand(M+1) p_lsq = leastsq(error_func, p_init, args=(x,y))p_lsq_regu=leastsq(error_plus_regu,p_init,args=(x,y)) #添加正则项后的参数学习结果：p_lsq_reguplt.plot(x_points,sine_func(x_points),label='real') #follow sin函数的plotplt.plot(x_points,fit_func(p_lsq_regu[0],x_points),label='regulation') #follow lasso的plotplt.plot(x_points,fit_func(p_lsq[0],x_points),label='fitting') # 原始fitting plotplt.plot(x,y,'bo', label = 'noise')#训练值样本点plt.legend()