数据拟合方法

第二讲 数据拟合方法
在实验中，实验和戡测常常会产生大量的数据。

为了解释这些数据或者根据这些数据做出预测、判断，给决策者提供重要的依据。

需要对测量数据进行拟合，寻找一个反映数据变化规律的函数。

数据拟合方法与数据插值方法不同，它所处理的数据量大而且不能保证每一个数据没有误差，所以要求一个函数严格通过每一个数据点是不合理的。

数据拟合方法求拟合函数，插值方法求插值函数。

这两类函数最大的不同之处是，对拟合函数不要求它通过所给的数据点，而插值函数则必须通过每一个数据点。

例如，在某化学反应中，测
–33显然，连续函数关系是
客观存在的。

但是通过表中的数据不可能确切地得到这
种关系。

何况，由于仪器和
环境的影响，测量数据难免有误差。

因此只能寻求一个
近拟表达式
y = ϕ(t )
寻求合理的近拟表达式，以反映数据变化的规律，这种方法就是数据拟合方法。

数
据拟合需要解决两个问题：
第一，选择什么类型的函数
)(t ϕ作为拟合函数（数学模
型）；第二，对于选定的拟合
函数，如何确定拟合函数中的参数。

数学模型应建立在合理假设的基础上，假设的合理性首先体现在选择某种类型的拟合函数使之符合数据变化的趋势（总体的变化规律）。

拟合函数的选择比较灵活，可以选择线性函数、多项式函数、指数函数、三角函数或其它函数，这应根据数据分布的趋势作出选
假设拟合函数是线性函数，即拟合函数的图形是一条平面上的直线。

而表中的数据点未能精确地落在一条直线上的原因是实验数据的误差。

则下一步是确定函数
y= a + b x
中系数a 和b 各等于多少？从几何背景来考虑，就是要以a 和b 作为待定系数，确定一条平面直线使得表中数据所对应的10个点尽可能地靠近这条直线。

一般来讲，数据点将不会全部落在这条直线上，如果第k 个点的数据恰好落在这条直线上，则这个点的坐标满足直线的方程，即
a +
b x k = y k
如果这个点不在直线上，则它的坐标不满足直线方程，有一个绝对值为k k y bx a -+的差异（残差）。

于是全部点处的总误差是
∑=-+10
1
k k k
y bx
a
这是关于a 和b 的一个二元函数，合理的做法是选取a 和b ，使得这个函数取极小值。

但是在实际求解问题时为了操作上的方便，常常是求a 和b 使得函数
∑=-+=10
1
2)(),(k k k y bx a b a F
达到极小。

为了求该函数的极小值点，令
0=∂∂a F ，0=∂∂b
F
， 得
0)(210
1
=-+∑=k k k
y bx
a ，
∑==-+10
1
0)(2k k k k
x y bx
a
这是关于未知数a 和b 的线性方程组。

它们被称为法方程，又可以写成
⎪⎪⎩⎪⎪⎨⎧
=+=+∑∑∑∑∑=====10
11012101
10
1
10110k k k k k k k k k k k y x b x a x y b x a 求解这个二元线性方程组便得待定系数a 和b ，从而得线性拟合函数 y = a + b x 。

下图中
直线是数据的线性拟合的结果。

假设拟合函数不是线性函数，而是一个二次多项式函数。

即拟合函数的图形是一条平面上的抛物线，而表中的数据点未能精确地落在这条抛物线上的原因是实验数据的误差。

则下一步是确定函数
y = a 0 + a 1 x + a 2 x 2
中系数a 0、a 1和a 2各等于多少？从几何背景来考虑，就是要以a 0、a 1和a 2为待定系数，确定二次曲线使得表中数据所对应的10个点尽可能地靠近这条曲线。

一般来讲，数据点将不会全部落在这条曲线上，如果第k 个点的数据恰好落在曲线上，则这个点的坐标满足二次曲线的方程，即
a 0 + a 1 x k + a 2 x k 2 = y k
如果这个点不在曲线上，则它的坐标不满足曲线方程，有一个误差（残差）。

于是全部点处的总误差用残差平方和表示
∑=-++=10
1
22
210210])[(),,(k k k k y x a x a a a a a F
这是关于a 0、a 1和a 2的一个三元函数，合理的做法是选取a 0、a 1和a 2 ，使得这个函数取极小值。

为了求该函数的极小值点，令
00=∂∂a F ，01=∂∂a F ，02
=∂∂a F
得
⎪⎪⎪⎩⎪⎪⎪⎨⎧=-++=-++=-++∑∑∑===10
1
2
22101012
2101012
2100])[(20])[(20])[(2k k k k k k k k k k k k k k x y x a x a a x y x a x a a y x a x a a 这是关于待定系数a 0、a 1和a 2的线性方程组，写成等价的形式为
⎪⎪⎪⎩⎪⎪
⎪⎨⎧=++=++=++∑∑∑∑∑∑∑∑∑∑∑===========1012
1012410113101
02
10
110110123
12101
010
1
10122
1101010k k
k k k k k k k k k k k k k k k k k k
k k k k y x a x a x a x y x a x a x a x y a x a x a
这就是法方程，求解这一方程组可得二次拟合函数中的三个待定系数。

下图反映了例题所
三． 数据的n
下面的函数
y = a 0 + a 1 x + a 2 x 2 + …… + a n x n
这里要做一个假设，即多项式的阶数n 应小于题目所给数据的数目m （例题中m = 10）。

类似前面的推导，可得数据的n 次多项式拟合中拟合函数的系数应满足的正规方程组如下
⎥⎥
⎥
⎥⎥⎥⎥⎥⎥⎦
⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡∑∑∑∑∑∑∑∑∑∑∑=====+==+====m k k n k m k k k m k k n m k n k m k n k m k n k m k n k m k k m k k m k n k m k k y x y x y a a a x x x x x x x x m 11110121111112111 从这一方程组可以看出，线性拟合方法和二次拟合方法是多项式拟合的特殊情况。

从算法
上看，数据最小二乘拟合的多项式方法是解一个超定方程组
⎪⎪
⎩⎪⎪⎨
⎧=++++=++++=++++m n m n m m n
n n n y
x a x a x a a y x a x a x a a y x a x a x a a 2210
2
222221*********( m > n ) 的最小二乘解。

而多项式拟合所引出的正规方程组恰好是用超定方程组的系数矩阵的转置
矩阵去左乘超定方程组左、右两端所得。

正规方程组的系数矩阵是一个病态矩阵，这类方程组被称为病态方程组。

当系数矩阵或者是右端向量有微小的误差时，可能引起方程组准确解有很大的误差。

为了避免求解这样的线性方程组，在做多项式拟合时可以将多项式中的各次幂函数做正交化变换，使得所推出的正规方程的系数矩阵是对角矩阵。

四．点集{x 1，x 2，……，x m }上的正交多项式系
多项式q 0(x )，q 1(x )，q 2(x )，……，q n (x )在点集{x 1，x 2，……，x m }上的正交
∑==m
i i j i k j k x q x q q q 1
)()(),(
正交多项式系可以认为是幂函数系：1，x ，x 2，……，x n 通过正交变换而得到的一组函数。

正交多项式系构造的方法如下：
q 0(x )=1，q 0(x ) = x – a 1 ，(a 1 =
n x
m
i i
/1∑=)，
q k (x ) = (x - a k ) q k -1(x ) - b k q k-2(x ) ，( k = 2，3，……，n )
其中，
∑∑=-=-----==m
i i k m
i i k i k k k k k x q x q x q q q xq a 12
11
2
11111)(/)(),/(),(
∑∑=-=-----==m
i i k m i i k k k k k k x q x q q q q q b 1
2
21
2
12211)(/)(),/(),(
五．用正交多项式系组成拟合函数的多项式拟合
考虑拟合函数：)()()()(
x q a x q a x q a x +++= ϕ，将数据表
⎪⎪
⎩⎪⎪⎨
⎧=++++=++++=++++m
m n n m m m n n n n y x q a x q a x q a x q a y x q a x q a x q a x q a y x q a x q a x q a x q a )()()()()()()()()()()()(2211002
222221120011122111100 （m > n ） 其系数矩阵为
⎥⎥⎥⎥
⎦⎤
⎢⎢⎢⎢⎣⎡)()()()()()()()()()()()(21022221201121110m n m m m
n n x q x q x q x q x q x q x q x q x q x q x q x q
由于多项式q 0(x )，q 1(x )，q 2(x )，……，q n (x )在点集{x 1，x 2，……，x m }上的正交，所以
超定方程组的系数矩阵中不同列的列向量是相互正交的向量组。

于是用这一矩阵的转置矩阵去左乘超定方程组左、右两端得正规方程组
⎪⎪⎩⎪⎪⎨
⎧===)
,(),(),(),(),(),(1
1110000y q a q q y q a q q y q a q q n n n n => ⎪⎪⎩⎪⎪⎨
⎧===)
,/(),(),/(),()
,/(),(1
1110000n n n n q q y q a q q y q a q q y q a 其中，∑==
m
i i k
k k x q
q q 1
2)(),(，∑==m
i i i k k y x q y q 1
)(),(。

因为正规方程组中每一个方程都是
一元一次方程可以直接写出原超方程组的最小二乘解，所以拟合函数为
)()
,(),()(),()
,()(),(),()(11110000x q q q y q x q q q y q x q q q y q x n n n n
+++=
ϕ
这一结果与用次多项式拟合所得结果在理论是完全一样的，只是形式上不同、算法实现上
避免了解病态方程组。

六．指数函数的数据拟合 问题1：世界人中预测问题
下表给出了本世纪六十年代世界人口的统计数据(单位：亿)
有人根据表中数据，预测公元2000年世界人口会超过 60亿。

这一结论在六十年代末令人难以置信，但现在已成为事实。

试建立数学模型并根据表中数据推算出2000年世界人口的数量。

根据马尔萨斯人口理论，人口数量按指数递增的规律发展。

记人口数为 N (t )，则有指数函数N e a bt
=+。

现需要根据六十年代的人口数据确定函数表达式中两个常数a 、b 。

为了计算方便，对表达式两边取对数，得 ln N a bt =+，令 N y ln =。

于是bt a t y +=)(。

(1)计算出表中人口数据的对数值y k = ln N k ( k = 1，2， (9)
(2) 根据表中数据写出关于两个未知数a 、b 的9个方程的超定方程组（方程数多于未知数个数的方程组）
a + b t k = y k ( k = 1，2， (9)
其中，t 1 =1960，t 2 =1961，t 3 =1962，……，t 9 =1968； y 1= ln29.72，y 2 = ln 30.61，……，y 9 = ln34.83。

(3) 利用MA TLAB 解线性方程组Ax=c 的命令A \c 计算出a 、b 的值，并写出人口增长函数。

利用人口增长函数计算出2000年世界人口数据：N (2000) 七．多元线性函数的数据拟合 问题2 人的耗氧能力的数据拟合。

人的耗氧能力y (ml/min ·kg)与下列变量有关
x 1 年龄 x 2 体重
x 3 1.5英里跑步所用时间 x 4 静止时心速 x 5 跑步时最大心速
某健身中心对31个自愿者进行测试，得到31组数据（每一组数据有6个数）
y k x 1k x 2k x 3k x 4k x 5k (k = 1，…31)
令耗氧能力为因变量，其它的指标为自变量，建立线性模型
y = a 0 + a 1 x 1+ a 2 x 2+ a 3 x 3+ a 4 x 4+ a 5 x 5
为了确定6个系数，利用已记录的数据得超定方程组
a 0 + a 1 x 1k + a 2 x 2k + a 3 x 3k + a 4 x 4k + a 5 x 5k = y k ( k = 1，2， (31)
这一方程组包含6个未知数 a 0 ，a 1， a 2 ，a 3 ，a 4 ，a 5 ，但却有31个方程。

写出超定方程组的系数矩阵和右端向量如下
⎥
⎥
⎥⎥⎦⎤⎢
⎢⎢
⎢⎣⎡=31,531
,431
,331
,231
,152********
5141312111111x x x x x x x x x x
x x x x x A ，⎥⎥⎥⎥⎦
⎤
⎢⎢⎢⎢⎣⎡=3121y y y y 由最小二乘法可得正规方程组
y A AX A T T =
其中，X = [a 0 a 1 a 2 a 3 a 4 a 5]T。