2023年8月3日发(作者:)
西北师大数学系
基于最小二乘法的多项式拟合
一 最小二乘法的基本原理设已知某物理过程
yf(x)的一组观测数据
(xi,f(xi)),
i1,2,,m. (1)
要求在某特定函数类(x)寻求一个函数(x)作为yf(x)的近似函数,使得二者在xi上的残差
i(xi)f(xi),
i1,2,,m (2)
按某种度量标准为最小,这就是拟合问题.
[0,1,,m]Tii要求残差按某种度量标准为最小,即要求由残差构成的残差向量的某种范数为最小,要求1,或m即
1i(xi)f(xi)i0i0m
maximax(xi)f(xi)ii
为最小,这本来都是很自然的,可是计算不太方便.通常要求:
或者
2(){[(xi)f(xi)]}2i2i0i0m12m12
22[(xi)f(xi)]22ii0i0mm (3)
为最小.这种要求误差平方和最小的拟合称为曲线拟合的最小二乘法.就是说,最小二乘法提供了一种数学方法,利用这种方法可以对实验数据实现在最小平方误差意义下的最好拟合.在曲线拟合中,函数类可有不同的选取方法.下面就常用的多项式拟合做介绍。
二 多项式拟合
(x,y)假设给定数据点ii(i=0,1,…,m),为所有次数不超过n(nm)的多项式构成的函数类,现求一pn(x)akxkk0n,使得
m2n2Ipn(xi)yiakxikyimini0i0k0 (4)
p(x)称为最小二乘拟合多项式。特别地,当n=1当拟合函数为多项式时,称为多项式拟合,满足式(4)的n时,称为线性拟合或直线拟合。显然
mnmI(akxikyi)2i0k0 西北师大数学系
n的多元函数,因此上述问题即为求为01必要条件,得
a,a,aII(a0,a1,an)的极值 问题。由多元函数求极值的mnI2(akxikyi)xij0,aji0k0j0,1,,n (5)
即
(xk0i0nmjki)akxijyi,i0mj0,1,,n (6)
(6)是关于a0,a1,an的线性方程组,用矩阵表示为
m1mxii0mxini0xxi0i0mmi2in1ixi0mmxyiai0i00mmn1axi1xiyii0i0ammn2nnxiyixii0 (7)
i0nim式(6)或式(7)称为正规方程组或法方程组。
可以证明,方程组(7)的系数矩阵是一个对称正定矩阵,故存在唯一解。从式(7)中解出n),从而可得多项式
ak(k=0,1,…,pn(x)akxkk0n (8)
可以证明,式(8)中的为最小二乘拟合多项式pn(x)满足式(4),即pn(x)为所求的拟合多项式。我们把i0pn(x)的平方误差,记作
222pmn(xi)yi2称r由式(5)可得
mpn(xi)yii0nmm
r22yak(xikyi)2ii0k0i0 (9)
多项式拟合的一般方法可归纳为以下几步:
(1) 由已知数据画出函数粗略的图形——散点图,确定拟合多项式的次数n;
(2) 列表计算xi0mji(j0,1,,2n)和xi0mjiyi(j0,1,,2n);
(3) 写出正规方程组,求出a0,a1,an;nk0
(4) 写出拟合多项式在实际应用中,npn(x)akxk。
m或nm;当nm时所得的拟合多项式就是拉格朗日或牛顿插值多项式。R()如表6-1,求电阻R与温度 T的近似函数关系。T例1 测得铜导线在温度i(℃)时的电阻ii 0
19.1
1
25.0
2
30.1
3
36.0
4
40.0
5
45.1
6
Ti(℃)
50.0 西北师大数学系
Ri()
76.30 77.80 79.25 80.80 82.35 83.90 85.10
解 画出散点图(图6-2),可见测得的数据接近一条直线,故取n=1,拟合函数为
Ra0a1T列表如下
i
0
1
2
3
4
5
6
Ti
19.1
25.0
30.1
36.0
40.0
45.1
50.0
245.3
Ri
76.30
77.80
79.25
80.80
82.35
83.90
85.10
565.5
Ti2
364.81
625.00
906.01
1296.00
1600.00
2034.01
2500.00
9325.83
TiRi
1457.330
1945.000
2385.425
2908.800
3294.000
3783.890
4255.000
20029.445
正规方程组为
解方程组得
245.3a0565.57245.39325.83a20029.4451
a070.572,故得R与T的拟合直线为
利用上述关系式,可以预测不同温度时铜导线的电阻值。例如,由R=0得T=-242.5,即预测温度时,铜导线无电阻。
a10.921
T=-242.5℃R70.5720.921T
6-2
例2 已知实验数据如下表
i 0
1
10
1
3
5
2
4
4
3
5
2
4
6
1
5
7
1
6
8
2
7
9
3
8
10
4
xi
yi
试用最小二乘法求它的二次拟合多项式。解 设拟合曲线方程为
ya0a1xa2x2
列表如下 西北师大数学系
I
0
1
2
3
4
5
6
7
8
xi
1
3
4
5
6
7
8
9
10
53
yi
10
5
4
2
1
1
2
3
4
32
xi2
1
9
16
25
36
49
64
81
100
381
xi3
1
27
64
125
216
343
512
729
1000
3017
xi4
1
81
256
625
1296
2401
4096
6561
10000
25317
xiyi
10
15
16
10
6
7
16
27
40
147
xi2yi
10
45
64
50
36
49
128
243
400
1025
得正规方程组
解得
52381a0329523813017a1471381301725317a21025
a013.4597,故拟合多项式为
a13.6053a20.2676
y13.45973.60530.2676x2
三 多项式拟合的MATLAB实现
用Polyfit函数P=polyfit(x,y,n)对数据进行拟合,返回n次多项式的系数,并用降序排列的向量表示,长度为n+1.
p(x)p1xnp2xn1pnxpn1
[p,s]=polyfit(x,y,n)返回多项式系数向量p和矩阵s。s与polyval函数一起用时,可以得到预测值的误差估计。
例如在MATLAB界面中输入一下命令
>> x=[0 0.0385 0.0963 0.1925 0.2888 0.385];
>> y=[0.042 0.104 0.186 0.338 0.479 0.612];
>> [p,s,mu]=polyfit(x,y,5)
输出结果为:
p =
Columns 1 through 5
0.0193 -0.0110 -0.0430 0.0073 0.2449
Column 6
0.2961
说明拟合的多项式为:0.0193x0.0110x0.043x0.0073xs =
R: [6x6 double]
df: 0
normr: 2.3684e-016
mu =
0.1669
0.1499
自由度为 0
标准偏差为 2.3684e-016
54320.2449x0.2961
2023年8月3日发(作者:)
西北师大数学系
基于最小二乘法的多项式拟合
一 最小二乘法的基本原理设已知某物理过程
yf(x)的一组观测数据
(xi,f(xi)),
i1,2,,m. (1)
要求在某特定函数类(x)寻求一个函数(x)作为yf(x)的近似函数,使得二者在xi上的残差
i(xi)f(xi),
i1,2,,m (2)
按某种度量标准为最小,这就是拟合问题.
[0,1,,m]Tii要求残差按某种度量标准为最小,即要求由残差构成的残差向量的某种范数为最小,要求1,或m即
1i(xi)f(xi)i0i0m
maximax(xi)f(xi)ii
为最小,这本来都是很自然的,可是计算不太方便.通常要求:
或者
2(){[(xi)f(xi)]}2i2i0i0m12m12
22[(xi)f(xi)]22ii0i0mm (3)
为最小.这种要求误差平方和最小的拟合称为曲线拟合的最小二乘法.就是说,最小二乘法提供了一种数学方法,利用这种方法可以对实验数据实现在最小平方误差意义下的最好拟合.在曲线拟合中,函数类可有不同的选取方法.下面就常用的多项式拟合做介绍。
二 多项式拟合
(x,y)假设给定数据点ii(i=0,1,…,m),为所有次数不超过n(nm)的多项式构成的函数类,现求一pn(x)akxkk0n,使得
m2n2Ipn(xi)yiakxikyimini0i0k0 (4)
p(x)称为最小二乘拟合多项式。特别地,当n=1当拟合函数为多项式时,称为多项式拟合,满足式(4)的n时,称为线性拟合或直线拟合。显然
mnmI(akxikyi)2i0k0 西北师大数学系
n的多元函数,因此上述问题即为求为01必要条件,得
a,a,aII(a0,a1,an)的极值 问题。由多元函数求极值的mnI2(akxikyi)xij0,aji0k0j0,1,,n (5)
即
(xk0i0nmjki)akxijyi,i0mj0,1,,n (6)
(6)是关于a0,a1,an的线性方程组,用矩阵表示为
m1mxii0mxini0xxi0i0mmi2in1ixi0mmxyiai0i00mmn1axi1xiyii0i0ammn2nnxiyixii0 (7)
i0nim式(6)或式(7)称为正规方程组或法方程组。
可以证明,方程组(7)的系数矩阵是一个对称正定矩阵,故存在唯一解。从式(7)中解出n),从而可得多项式
ak(k=0,1,…,pn(x)akxkk0n (8)
可以证明,式(8)中的为最小二乘拟合多项式pn(x)满足式(4),即pn(x)为所求的拟合多项式。我们把i0pn(x)的平方误差,记作
222pmn(xi)yi2称r由式(5)可得
mpn(xi)yii0nmm
r22yak(xikyi)2ii0k0i0 (9)
多项式拟合的一般方法可归纳为以下几步:
(1) 由已知数据画出函数粗略的图形——散点图,确定拟合多项式的次数n;
(2) 列表计算xi0mji(j0,1,,2n)和xi0mjiyi(j0,1,,2n);
(3) 写出正规方程组,求出a0,a1,an;nk0
(4) 写出拟合多项式在实际应用中,npn(x)akxk。
m或nm;当nm时所得的拟合多项式就是拉格朗日或牛顿插值多项式。R()如表6-1,求电阻R与温度 T的近似函数关系。T例1 测得铜导线在温度i(℃)时的电阻ii 0
19.1
1
25.0
2
30.1
3
36.0
4
40.0
5
45.1
6
Ti(℃)
50.0 西北师大数学系
Ri()
76.30 77.80 79.25 80.80 82.35 83.90 85.10
解 画出散点图(图6-2),可见测得的数据接近一条直线,故取n=1,拟合函数为
Ra0a1T列表如下
i
0
1
2
3
4
5
6
Ti
19.1
25.0
30.1
36.0
40.0
45.1
50.0
245.3
Ri
76.30
77.80
79.25
80.80
82.35
83.90
85.10
565.5
Ti2
364.81
625.00
906.01
1296.00
1600.00
2034.01
2500.00
9325.83
TiRi
1457.330
1945.000
2385.425
2908.800
3294.000
3783.890
4255.000
20029.445
正规方程组为
解方程组得
245.3a0565.57245.39325.83a20029.4451
a070.572,故得R与T的拟合直线为
利用上述关系式,可以预测不同温度时铜导线的电阻值。例如,由R=0得T=-242.5,即预测温度时,铜导线无电阻。
a10.921
T=-242.5℃R70.5720.921T
6-2
例2 已知实验数据如下表
i 0
1
10
1
3
5
2
4
4
3
5
2
4
6
1
5
7
1
6
8
2
7
9
3
8
10
4
xi
yi
试用最小二乘法求它的二次拟合多项式。解 设拟合曲线方程为
ya0a1xa2x2
列表如下 西北师大数学系
I
0
1
2
3
4
5
6
7
8
xi
1
3
4
5
6
7
8
9
10
53
yi
10
5
4
2
1
1
2
3
4
32
xi2
1
9
16
25
36
49
64
81
100
381
xi3
1
27
64
125
216
343
512
729
1000
3017
xi4
1
81
256
625
1296
2401
4096
6561
10000
25317
xiyi
10
15
16
10
6
7
16
27
40
147
xi2yi
10
45
64
50
36
49
128
243
400
1025
得正规方程组
解得
52381a0329523813017a1471381301725317a21025
a013.4597,故拟合多项式为
a13.6053a20.2676
y13.45973.60530.2676x2
三 多项式拟合的MATLAB实现
用Polyfit函数P=polyfit(x,y,n)对数据进行拟合,返回n次多项式的系数,并用降序排列的向量表示,长度为n+1.
p(x)p1xnp2xn1pnxpn1
[p,s]=polyfit(x,y,n)返回多项式系数向量p和矩阵s。s与polyval函数一起用时,可以得到预测值的误差估计。
例如在MATLAB界面中输入一下命令
>> x=[0 0.0385 0.0963 0.1925 0.2888 0.385];
>> y=[0.042 0.104 0.186 0.338 0.479 0.612];
>> [p,s,mu]=polyfit(x,y,5)
输出结果为:
p =
Columns 1 through 5
0.0193 -0.0110 -0.0430 0.0073 0.2449
Column 6
0.2961
说明拟合的多项式为:0.0193x0.0110x0.043x0.0073xs =
R: [6x6 double]
df: 0
normr: 2.3684e-016
mu =
0.1669
0.1499
自由度为 0
标准偏差为 2.3684e-016
54320.2449x0.2961
发布评论