最小二乘法求二次方程系数

例1：二次方程式计算 Y=a0+a1x+a2 x2 y=-6.3+2.4x+1.3x2

下表为自动计算系数，给出9组x 和y 的数值，自动计算出系数。

x y

x^2

x^3

x^4xy x^2y11-2.6111-2.6123.748167.41312.69278137.81424.1166425696.41538.[**************]4.[**************]74.[**************]96.[1**********]691

9120.[1**********]085求和

9

45

421.8

285

2025

[***********].[***********][1**********]97.4

系数系数值a0-6.30x y a12.40896.10a21.30

原理与多项式拟合说明附后。

-2.614.8113.4385.69551976.43635.86150.49768.622997.4

第一节 最小二乘法的基本原理和多项式拟合 一 最小二乘法的基本原理

从整体上考虑近似函数p (x ) 同所给数据点(x i , y i ) (i=0,1,…,m) 误差

r i =p (x i ) -y i (i=0,1,…,m) 的大小，常用的方法有以下三种：一是误差r i

r i =p (x i ) -y i (i=0,1,…,m) 绝对值的最大值max 0≤i ≤m ，即误差 向量r =(r 0, r 1, r m ) T 的∞—范数；二是误差绝对值的和∑i =0

m

r i

，即误差向量r 的1—

范数；三是误差平方和i =0的算术平方根，即误差向量r 的2—范数；前两种方法简单、自然，但不便于微分运算 ，后一种方法相当于考虑 2—范数的平方，因此在曲线拟合中常采用误差平方和i =0体大小。

∑r

m

2

i

∑r

m

2

i

来 度量误差r i (i=0，1，…，m) 的整

数据拟合的具体作法是：对给定数据 (x i , y i ) (i=0,1,…，m) ，在取定的函数类Φ中, 求p (x ) ∈Φ, 使误差r i =p (x i ) -y i (i=0,1,…,m) 的平方和最小，即

i =0

∑r i

m

2

=i =0

2

[]p (x ) -y ∑i i =min

m

从几何意义上讲，就是寻求与给定点(x i , y i ) (i=0,1,…,m) 的距离平方和为最

小的曲线 y =p (x ) （图6-1）。函数p (x ) 称为拟合 函数或最小二乘解，求拟合函数p (x ) 的方法称为曲线拟合的最小二乘法。 在曲线拟合中，函数类Φ可有不同的选取方法

.

6—1

二 多项式拟合

假设给定数据点(x i , y i ) (i=0,1,…,m) ，Φ为所有次数不超过n (n ≤m ) 的多项式构成的函数类，现求一

m

p n (x ) =∑a k x k ∈Φ

k =0

n

, 使得

2

⎛n ⎫2

I =∑[p n (x i ) -y i ]=∑ ∑a k x i k -y i ⎪=min

i =0i =0⎝k =0⎭ (1)

当拟合函数为多项式时，称为多项式拟合，满足式（1）的p n (x ) 称为最小二乘

m

拟合多项式。特别地，当n=1时，称为线性拟合或直线拟合。

显然

为a 0, a 1, a n 的多元函数，因此上述问题即为求I =I (a 0, a 1, a n ) 的极值 问题。由多元函数求极值的必要条件，得

m n

∂I

=2∑(∑a k x i k -y i ) x i j =0, j =0, 1, , n ∂a j i =0k =0

(2)

i =0

k =0

I =∑(∑a k x i k -y i ) 2

m n

即

i =0 k =0i =0

（3）是关于a 0, a 1, a n 的线性方程组，用矩阵表示为

∑(∑x

⎡

⎢m +1⎢m ⎢x

i

⎢∑i =0⎢ ⎢m

⎢∑x i n ⎢⎣i =0

n m

j +k i

) a k =∑x i j y i ,

m

j =0, 1, , n

(3)

∑x ∑x

i =0i =0

m

m

i

2i

∑x

i =0

m

n +1i

⎤⎡m ⎤

∑x ⎥y ⎢∑i ⎥a i =0⎡⎤i =0⎥0⎢m ⎥m

n +1⎥⎢a ⎥ ∑x i ⎢1⎥⎢∑x i y i ⎥

=⎢i =0⎥⎥i =0

⎢⎥

⎥⎢⎥⎢ ⎥

⎥a ⎢m ⎥m

n 2n ⎥⎣n ⎦⎢∑x i y i ⎥ ∑x i

⎢⎥⎥⎣i =0⎦ (4) i =0⎦

n

i

m

式（3）或式（4）称为正规方程组或法方程组。

可以证明，方程组（4）的系数矩阵是一个对称正定矩阵，故存在唯一解。从式（4）中解出a k (k=0,1,…，n) ，从而可得多项式

(5)

可以证明，式（5）中的p n (x ) 满足式（1），即p n (x ) 为所求的拟合多项式。我

k =0

p n (x ) =∑a k x k

n

们把i =0

∑[p

m

n

(x i ) -y i ]

2

称为最小二乘拟合多项式p n (x ) 的平方误差，记作

r

22

=∑[p n (x i ) -y i ]

i =0n

m

m

2

由式(2)可得

r

22

=∑y -∑a k (∑x i k y i )

2i i =0

k =0

i =0

m

(6)

多项式拟合的一般方法可归纳为以下几步：

(1) 由已知数据画出函数粗略的图形——散点图，确定拟合多项式的次数n ；

(2) 列表计算i =0

和i =0

(3) 写出正规方程组，求出a 0, a 1, a n ； (4) 写出拟合多项式

p n (x ) =∑a k x k

k =0n

∑x

m

j i

(j =0, 1, , 2n )

∑x

m

j i

y i

(j =0, 1, , 2n )

；

。

在实际应用中，n

例1 测得铜导线在温度T i (℃) 时的电阻R i (Ω) 如表6-1，求电阻R 与温度 T

数为

R =a 0+a 1T

245. 3⎤⎡a 0⎤⎡565. 5⎤⎡7

⎢245. 39325. 83⎥⎢a ⎥=⎢20029. 445⎥⎣⎦⎣1⎦⎣⎦

解方程组得

a 0=70. 572, a 1=0. 921

故得R 与T 的拟合直线为

R =70. 572+0. 921T

利用上述关系式，可以预测不同温度时铜导线的电阻值。例如，由R=0得T=-242.5，即预测温度 T=-242.5℃时，铜导线无电阻。

6-2

例2 例2 已知实验数据如下表

解 设拟合曲线方程为

2

y =a 0+a 1x +a 2x

52⎡9

⎢52381⎢⎢⎣3813017

381⎤⎡a 0⎤⎡32⎤

⎢a ⎥=⎢147⎥3017⎥⎥⎢1⎥⎢⎥

25317⎥⎦⎢⎣a 2⎥⎦⎢⎣1025⎥⎦

解得

a 0=13. 4597, a 1=-3. 6053a 2=0. 2676

故拟合多项式为

y =13. 4597-3. 6053+0. 2676x 2

*三 最小二乘拟合多项式的存在唯一性

定理1 设节点x 0, x 1, , x n 互异，则法方程组（4）的解存在唯一。

证 由克莱姆法则，只需证明方程组（4）的系数矩阵非奇异即可。 用反证法，设方程组（4）的系数矩阵奇异，则其所对应的齐次方程组

m

⎡

⎢m +1∑x i

i =0

⎢m m

2⎢x x ∑∑i

⎢i =0i

i =0

⎢ ⎢m m

⎢∑x i n ∑x i n +1 ⎢i =0⎣i =0

有非零解。式(7)可写为

⎤⎡m ⎤x y ∑⎥⎢∑i ⎥a i =0⎡⎤i =0⎥0⎢m ⎥m

n +1⎥⎢a ⎥x i ⎢1⎥⎢∑x i y i ⎥∑=⎢i =0⎥⎥i =0

⎢⎥

⎥⎢⎥⎢ ⎥

⎥a ⎢m ⎥m

n 2n ⎥⎣n ⎦⎢∑x i y i ⎥x i ∑⎢⎥⎥⎣i =0⎦ （7） i =0⎦

n

i

m

∑(∑x

k =0

i =0

n m

j +k i

) a k =0,

j =0, 1, , n

（8）

n

将式（8）中第j 个方程乘以a j (j=0,1,…，n) ，然后将新得到的n+1个方程左

⎡n m j +k ⎤a j ⎢∑(∑x i ) a k 0⎥=0∑⎦ 右两端分别 相加，得j =0⎣k =0i =0

因为

n m n n m

⎡n m j +k ⎤m n n 2j +k j k

[]a (x ) a =a a x =(a x )(a x ) =p (x ) ∑∑∑j i ∑∑j ⎢∑∑i k ⎥∑∑∑k j i k i n i

j =0i =0j =0k =0i =0⎣k =0i =0⎦i =0j =0k =0 其中

p n (x ) =∑a k x k

k =0n

所以

p n (x i ) =0 (i=0,1,…,m)

p n (x ) 是次数不超过n 的多项式，它有m+1＞n 个相异零点，由代数基本定理，必

须有a 0=a 1= a n =0，与齐次方程组有非零解的假设矛盾。因此正规方程组（4）

a a , , a n 必有唯一解 。定理2 设0, 1是正规方程组（4）的解，则是满足式（1）的最小二乘拟合多项式。

p n (x ) =∑a k x k

k =0

n

n

证 只需证明，对任意一组数

m i =0

b 0, b 1, , b n

m

组成的多项式

2

Q n (x ) =∑b k x k

k =0

，恒有

∑[Q n (x i ) -y i ]

即可。

2

≥∑[p n (x i ) -y i ]

i =0

∑[Q

i =0

m i =0

m

n

(x i ) -y i ]-∑[p n (x i ) -y i ]

2

i =0

2

m

m

2

=∑[Q n (x i ) -p n (x i ) ]+2∑[Q n (x i ) -p n (x i ) ]⋅[p n (x i ) -y i ]

i =0

≥0+2∑∑

i =0j =0

m n

[

n ⎧m

⎡⎛n ⎡n ⎤⎫j ⎤⎫⎪⎪k k

(b j -a j ) x i ⋅⎢∑a k x i -y i ⎥=2∑⎨(b j -a j )∑⎢ ∑a k x i -y i ⎪x i ⎥⎬

j =0⎪i =0⎣⎝k =0⎣k =0⎦⎭⎦⎪⎩⎭

j

]

因为a k (k=0,1,…，n) 是正规方程组（4）的解，所以满足式（2），因此有

2

[][]Q (x ) -y -p (x ) -y ∑n i i ∑n i i ≥0

2

i =0

i =0

m

m

故p n (x ) 为最小二乘拟合多项式。

*四 多项式拟合中克服正规方程组的病态

在多项式拟合中，当拟合多项式的次数较高时，其正规方程组往往是病态的。而且

①正规方程组系数矩阵的阶数越高，病态越严重；

②拟合节点分布的区间[x 0, x m ]偏离原点越远，病态越严重； ③x i (i=0,1,…，m) 的数量级相差越大，病态越严重。 为了克服以上缺点，一般采用以下措施：

①尽量少作高次拟合多项式，而作不同的分段低次拟合；

②不使用原始节点作拟合，将节点分布区间作平移，使新的节点x i 关于原 点对称，可大大降低正规方程组的条件数，从而减低病态程度。 平移公式为：

x +x m

x i =x i -0, i =0, 1, , m

2 (9) ③对平移后的节点x i (i=0,1,…，m), 再作压缩或扩张处理：

x i *=p x i ,

p =(m +1)

i =0, 1, , m （10） , （r 是拟合次数） （11）

其中

∑(x )

i

i =0

m

2r

*

x i 经过这样调整可以使的数量级不太大也不太小，特别对于等距节点

x i =x 0+ih (i =0, 1, , m ) ，作式（10）和式（11）两项变换后，其正规方程

组的系数矩阵设 为A ，则对1～4次多项式拟合，条件数都不太大，都可以得到满意的结果。

④在实际应用中还可以利用正交多项式求拟合多项式。一种方法是构造离散正交多项式；另一种方法是利用切比雪夫节点求出函数值后再使用正交多项式。这两种方法都使正规方程 组的系数矩阵为对角矩阵，从而避免了正规方程组的病态。我们只介绍第一种，见第三节。 例如 m=19,x 0=328,h=1, x 1=x 0+ih，i=0,1,…，19，即节点 分布在［328,347］，作二次多项式拟合时

① 直接用x i 构造正规方程组系数矩阵A 0，计算可得

cond 2(A 0) =2. 25⨯1016

严重病态，拟合结果完全不能用。 ② 作平移变换

328+347

x i =x i -, i =0, 1, , 19

2

用x i 构造正规方程组系数矩阵A 1，计算可得

cond 2(A 1) =4. 483868⨯1016

比cond 2(A 0) 降低了13个数量级，病态显著改善，拟合效果较好。 ③ 取压缩因子

p =

20

∑(x )

i

i =0

19

≈0. 1498

4

*x i 作压缩变换 =p x i ,

i =0, 1, , 19

*

x 用i 构造正规方程组系数矩阵A 2，计算可得 cond 2(A 2) =6. 839

又比cond 2(A 1) 降低了3个数量级，是良态的方程组，拟合效果十分理想。

*

p (x ) 中使用原来节点所对应的变量x ，可写为n 如有必要，在得到的拟合多项式

x 0+x m

)) 2

仍为一个关于x 的n 次多项式，正是我们要求的拟合多项式。

Q n (x ) =p n (p ⋅(x -