Download presentation
Presentation is loading. Please wait.
1
第 五 章 插 值 法 与曲线拟合 插值法
2
上一讲的简单回顾 ● 插值多项式的存在惟一性: 满足插值条件 Pn(xi)=f(xi), ( i=0,1,2,…,n)
n次插值多项式Pn(x)=a0+a1x+a2x2+……+anxn 存在而且惟一. ● 插值余项: Rn(x)= f(x)- Pn(x)= , ● Lagrange插值多项式 ,i =0,1,…,n 其中 称为Lagrange插值基函数
3
§ 牛顿途径 对于n+1个不同的节点 ,考虑n次多项式 (6) 如果满足: ,那么它就是n+1个点上的n次插值多项式,对于这样的 ,有
4
本讲主要内容: 缺点: 不具有承袭性,即每当增加一个节点时,不仅要增加求 和的项数,而且以前的各项也必须重新计算.
优点: 具有严格的规律性,便于记忆. 缺点: 不具有承袭性,即每当增加一个节点时,不仅要增加求 和的项数,而且以前的各项也必须重新计算. 为了克服这一缺点,本讲将建立具有承袭性的插值公式— Newton插值公式. 本讲主要内容: ● Newton插值多项式的构造 ● 差商的定义及性质 ● 差分的定义及性质 ● 等距节点Newton插值公式
5
一 基函数 问题1 求作n次多项式 使满足 (1) (2) 为了使 的形式得到简化,引入如下记号 (3)
6
称为Newton插值的以x0,,x1,…,xn为节点的基函数,即
可以证明,这样选取的基函数是线性无关的,由此得 出的问题4.1的解便于求值,而且新增加一个节点 xn+1时 只需加一个新项 即 而 (4)
7
依据条件(2),可以依次确定系数c0,c1,…,cn..例如, 取x=x0,,得 取x=x1 ,得
8
. 为了得到计算系数ci的一般方法,下面引进差商的概念.
9
二 差商的定义 给定[a,b]中互不相同的点x0,x1,x2,…,以及f(x)在这些点处相
二 差商的定义 给定[a,b]中互不相同的点x0,x1,x2,…,以及f(x)在这些点处相 应的函数值f(x0),f(x1),f(x2),…,用记号 表示f(x)在x0及x1两点的一阶差商. 用记号 表示f(x)在x0,x1,x2三点的二阶差商. 一般地,有了k-1阶差商之后, 可以定义f(x)在x0,x1,..,xk的k阶差商
10
三 Newton插值公式 由差商定义,有 f(x)= f[x0]+(x-x0)f[x,x0]
f[x,x0]= f[x0,x1]+(x-x1)f[x,x0,x1] f[x,x0,x1]= f[x0,x1,x2]+(x-x2)f[x,x0,x1,x2] ……….. f[x,x0,…xn-1]= f[x0,…,xn]+(x-xn)f[x,x0,….,xn] 将以上各式,由下而上逐步代入,得到 f(x)= f(x0)+(x-x0) f[x0,x1]+(x-x0)(x-x1) f[x0,x1,x2] +…+(x-x0)…(x-xn-1) f[x0,…,xn] (5) +(x-x0)…(x-xn-1)(x-xn)f[x,x0,…xn]
11
记 Nn(x)= f(x0)+(x-x0) f[x0,x1]+(x-x0)(x-x1) f[x0,x1,x2] …+(x-x0)…(x-xn-1) f[x0,…,xn] (6) Rn(x)= (x-x0)…(x-xn) f[x,x0,…,xn]= f[x,x0,…,xn]wn+1(x) (7) 则(5)可表示为 f(x)= Nn(x)+ Rn(x) (8) 显然, Nn(x)是次数不超过n的多项式,且有 Rn(xi)= f[x,x0,…,xn]wn+1(xi)=0, i=0,1,…,n 即 Nn(xi)= f(xi) , i=0,1,…,n 由此可知,如此构造的函数Nn(x)就是问题1的解,且 ci =f[x0,…,xi] , i=0,1,…,n (9)
12
公式(6)称为函数f(x)在节点x0,…,xn上的n阶Newton 插值公式, (7)式称为Newton插值公式余项,即截断误差.注
意到,余项表达式(7)只要求被插值函数f(x)在插值区间[a,b] 上连续. 由函数f(x)的插值多项式的惟一性,函数f(x) 的Newton 插值多项式与Lagrange插值多项式实为同一个多项式,即 Nn(x) ≡ Ln(x) 两者不过是表现形式不同而已.由此有: 若f(x)∈Cn+1[a,b], 则有 Rn(x)= f[x,x0,…,xn]wn+1(x)= (10) ,
13
四 差商的性质 性质1(差商与函数值的关系) 证明: 性质2(对称性) 其中i0,…,ik是0,1,…k的任排列. 证明: 由性质1可知.
14
性质3(差商与导数的关系) f(x)∈Ck[a,b],
证明: 由式(10)即得. 性质4(多项式的差商)设f(x)为n次多项式,则其一阶差商 是x的n-1次多项式 证明: 推论 n次多项式pn(x)的k阶差商pn[x0,…xk],当k≤n时 是n-k次多项式,当k>n时恒为0
15
运用Newton插值公式(6)进行计算时,先计算f(x)关于 节点x0,…,xn的各阶差商.计算过程如下表所示 .
xk f(xk) 一阶差商 二阶差商 三阶差商 n 阶差商
16
计算Nn(x)时,常采用秦九韶算法,即 下面给出Newton插值法的计算机算法.开始时,f(k)存放函数值f(xk),运算完毕f(k)存放k阶差商f[x0,…,xk]
17
Newton插值算法 (1) 输入: xi, fi; di=fi (i=0,1,…,n); 计算差商 对i=1,2,…,n 做 (3.1) 对j=i,i+1,…,n 做 fj=(dj-dj-1)/(xj-xj-i); (3.2) 对j=i,i+1,…,n 做 dj=fj; (4) 计算插值N(u) (4.1)输入插值点u; (4.2) v=0; (4.3) 对i=n,n-1,…,1,0 做 v=v(u-xi)+fi; (5) 输出u,v.
18
五 等距节点的Newton插值公式与差分 (6)可以简化.
当插值节点x0,…,xn为等距分布时,Newton插值公式 (6)可以简化. 设插值节点xj=x0+jh, j=0,1,…,n; h=(b-a)/n称为步长,且 x0=a, xn=b. 令x=x0+th,则当x0≤x≤xn时,0≤t≤n. 基函数 此时差商也可进一步化简,为此引进差分的概念. 定义 称 △f(xi)=f(xi+h)-f(xi) 和 ▽f(xi)= f(xi)-f(xi-h) 分别为函数f(x)在点xi处的一阶向前差分和一阶向后差分.
19
一般地,称k阶差分的差分为k+1阶差分,如二阶向前和向
后差分分别为 计算各阶差分可按如下差分表进行.
20
向前差分表
21
性质3(多项式的差分)若f(x)∈Pn(n次多项式类), 则
差分具有如下性质: 性质1(差分与函数值的关系) 各阶差分均可表示为函值 fj=f(xj)的线性组合: 性质2(前差与后差的关系): 性质3(多项式的差分)若f(x)∈Pn(n次多项式类), 则 其中
22
性质4(差分与差商的关系): 性质5(差分与导数的关系 利用这些性质,可将Newton公式 Nn(x)= f(x0)+(x-x0) f[x0,x1]+(x-x0)(x-x1) f[x0,x1,x2] …+(x-x0)…(x-xn-1) f[x0,…,xn] 进一步简化
23
(11) 称公式(11)为Newton向前差分插值公式,其余项为 (12) 如果将式(6)改为按节点xn,xn-1,…,x0的次序排列的Newton插 值公式,即
24
令x=xn-th, 则当x0≤x≤xn时,0≤t≤n.利用差商与向后
差分的关系, 式(13)可简化为 (13) (14) 称式(14)为Newton向后差分插值公式,其余项为
25
若要计算的插值点x较靠近点x0,则用向前插值公式(4. 8),这时t=(x-x0)/n的值较小,数值稳定性较好
若要计算的插值点x较靠近点x0,则用向前插值公式(4.8),这时t=(x-x0)/n的值较小,数值稳定性较好.反之,若x靠近xn,,,则用向后插值公式(14). 利用向前与向后差分的关系(差分性质2): 式(14)可表示成 (15) 这样,计算靠近x0或xn的点的值时,都只需构造向前差分表
26
例 给定f(x)在等距节点上的函数值表如下: xi 0.4 0.6 0.8 1.0
f(xi) 分别用Newton向前和向后差分公式,求f(0.5)及f(0.9)的近似值. 解 先构造向前差分表如下: xi fi △fi △2fi △3fi x0=0.4, h=0.2, x3=1.0. 由(4.8)和(4.12),分别用差分表中对角线上的值和最后一行的值,得Newton向前和向后插值公式如 下:
27
(1) (2) 当 x=0.5时,用公式(1),这时t=(x-x0)/h=0.5. 将t=0.5代入(1),得 f (0.5)≈N3(0.5)= 当x=0.9时, 用公式(2), 这时t=(x3-x)/h=0.5. 将t=0.5代入(2), 得 f (0.9)≈N3(0.9)=
28
课后练习题 P217: 第2题, 第4题, P219: 第11题.
29
§6 曲线拟合 6.1.2 曲线拟合问题 仍然是已知 x1 … xm ; y1 … ym, 求一个简单易算的近似函数 f(x) 来拟合这些数据。 但是① m 很大; ② yi 本身是测量值,不准确,即 yi f (xi) 这时没必要取 f(xi) = yi , 而要使 i=f(xi) yi 总体上 尽可能地小。 这种构造近似函数 的方法称为曲线拟合,f(x) 称为拟合函数 称为“残差”
30
“使 i=P(xi) yi 尽可能地小”有不同的准则
常见做法: 使 最小 不可导,求解困难,P283 使 最小 使 最小
31
6.2 线性拟合问题 6.2.1 ||.||2 意义下的线性拟合(线性最小二乘问题) 确定拟合函数
,对于一组数据(xi, yi) (i = 1, 2, …, m) 使得 达到极小,这里 n <= m。 Denote:
32
称方程组Ax=b为超定方程组
33
记 E 实际上是 c0, c1, …, cn 的多元函数,在 E 的极值点应有
34
得到关于c1,c2,…,cn的方程组 法方程组(或正规方程组)
35
例1 数据 ti fi
36
6.3 线性最小二乘问题 设A是m×n阶矩阵(m>n), 称线性方程组 Ax=b (1) 为超定方程组; 这里x∈Rn,b∈Rm.
6.3 线性最小二乘问题 设A是m×n阶矩阵(m>n), 称线性方程组 Ax=b (1) 为超定方程组; 这里x∈Rn,b∈Rm. 如果A的秩r(A)=n, 称A为列满秩矩阵. 记残向量r=b-Ax,考虑确定一个向量x, 使‖r‖2 2=‖b-Ax‖2 2, 达到最小的问题称为线性最小二乘问题, 这样的x称为方程组(1)的最小二乘解.
37
6.3.4 最小二乘解的存在惟一性 结论1 :设A是m×n阶矩阵,x∈Rn, b∈Rm. 由线性方程组理论可知,线性方程组
最小二乘解的存在惟一性 结论1 :设A是m×n阶矩阵,x∈Rn, b∈Rm. 由线性方程组理论可知,线性方程组 Ax=b (24) 有解的充分必要条件是 r (A)= r (A|b) (25)
38
定理6. 3. 7假设方程组(24)有解,令x是其一个解. 那么,方程组(24)的所有解的集合为{x}+N(A)
定理6.3.7假设方程组(24)有解,令x是其一个解. 那么,方程组(24)的所有解的集合为{x}+N(A). 方程组(24)有 惟一解的充分必要条件是null(A)=0. 这里, null (A)表示A的核子空间的维数.
39
证明: 首先证明任意的向量y∈{x}+N(A)都是方程组(24)的解.
事实上,将y记为y=x+z, 其中z∈N(A), 即Az=0,x∈{x}. 因此, Ay=Ax+Az=b,即y满足方程组(24). 反过来, 若y满足方程组(24), 有 Ay-Ax =A(y-x)= 0, 即y-x∈N(A). 记y=x+(y-x),从而有y∈{x}+N(A). 惟一性. 因为齐次方程组Ax=0有惟一零解的 充 分必要条件是A为满秩矩阵,即null (A)=0.
40
定理6.3.8当m>n时, 超定方程组(1)的最小二乘解总是存在的. 最小二乘解惟一的充分必要条件是null (A)=0.
41
证: 记b=b1+b2, 其中b1∈R(A),b2∈N(AT).
对任意x∈Rn, Ax∈R(A), b1-Ax∈R(A). 因此, ‖r‖22=‖b-Ax‖22=‖(b1-Ax)+b2‖22. 由定理6.3.3的推论1和定理6.3.2, ‖r‖22=‖b1-Ax‖22+‖b2‖22. 要使‖r‖22达到最小等价于确定x,使‖b1-Ax‖22 为0, 即求方程组Ax=b1的解x. 因为b1,Ax, b1-Ax都是R(A)中的向量,因此,可以 把b1看成由A的列向量线性表示, 即b1=Ax. 换句话说,方程组Ax=b1的解总是存在的,从而方程 组(1)的最小二乘解也总是存在的. 惟一性的证明可直接由定理6.3.7得到.
42
6.3.1 正交性的有关性质 在线性代数欧氏空间理论中, 将R3中两个向量x,y之间的夹角φ满足的关系式
正交性的有关性质 在线性代数欧氏空间理论中, 将R3中两个向量x,y之间的夹角φ满足的关系式 xTy=‖x‖2‖y‖2cosφ (2) 推广到Rn. 设x,y∈Rn, 由Cauchy不等式 -1≤ ≤1 从而得到Rn中两个向量之间的夹角为 φ=arccos (3)
43
定理6.3.1 设x, y是Rn中的向量, x与y正交 的充分必要条件为xTy=0.
44
定理6.3.2:设x, y∈Rn, 且x⊥y,那 么: ‖x+y‖22=‖x‖22+‖y‖22.
= xTx+2yTx+yTy 而xTy=yTx=0, 因此 ‖x+y‖22=‖x‖22+‖y‖22 注:推广到Rn中的向量组α1,α2,…,αk, 如果αiTαj=0 (i≠j), 称α1,α2,…,αk是 正交向量组. 特别地: 如果‖αi‖2=1(i=1,2,…,k), 即 αiTαj=δij,称α1,α2,…,αk 为标准正交向量组.
45
设U是Rn中的子空间, x∈Rn. 如果x与U中任意向量正交, 称向量x与子空间U正交, 记为x⊥U.
设U,V是Rn中两个子空间, 如果任意x∈U和任意y∈V是正交的, 称子空间U与子空间V正交, 记为U⊥V. 设U,V是Rn中互补的子空间. 如果U⊥V, 那么称U,V互为正交补子空间, 记U=V⊥或V=U⊥. 可以证明, 一个子空间的正交补子空间是惟一的.
46
定理6.3.3 设A是n×k阶矩阵,x∈Rn, 那么下列三种情况是
等价的: ①x⊥R(A); ②ATx=0; ③x∈N(AT). 这里,N(AT)={ATx=0, x∈Rn}称为AT的核子空间. 证:由N(AT)的定义, ②与③显然等价. 下面证明①与②等价. 记A=(α1,α2,…,αk), 那么,αi∈R(A) (i=1,2,…,k). 假设x⊥R(A), 即αiTx=0 (i=1,2,…,k). 从而ATx=0. 另一方面,如果ATx=0, 那么有z∈Rk, 使Az=y∈R(A). 这时,yTx=zTATx=0,即x⊥y. 由z的任意性, 得Az是任意的, 因此x⊥R(A). 由这个定理, 容易得到: 推论1设A是n×k阶矩阵, 那么R(A)有惟一的正交补子 空间N(AT).
47
6.3 线性最小二乘问题 设A是m×n阶矩阵(m>n), 称线性方程组 Ax=b (1)
6.3 线性最小二乘问题 设A是m×n阶矩阵(m>n), 称线性方程组 Ax=b (1) 为超定方程组; 这里x∈Rn, b∈Rm. 如果A的秩r (A) =n, 称A为列满秩矩阵. 记残向量r=b-Ax, 考虑确定一个向量x, 使‖r‖2 2=‖b-Ax‖2 2, 达到最小的问题称为线性最小二乘问题, 这样的x称为方程组(1)的最小二乘解.
48
6.3.1 正交性的有关性质 在线性代数欧氏空间理论中, 将R3中两个向量x,y之间的夹角φ满足的关系式
正交性的有关性质 在线性代数欧氏空间理论中, 将R3中两个向量x,y之间的夹角φ满足的关系式 xTy=‖x‖2‖y‖2cosφ (2) 推广到Rn. 设x,y∈Rn, 由Cauchy不等式 -1≤ ≤1 从而得到Rn中两个向量之间的夹角为 φ=arccos (3)
49
定理6.3.1 设x, y是Rn中的向量, x与y正交 的充分必要条件为xTy=0.
50
定理6.3.2:设x, y∈Rn, 且x⊥y,那 么: ‖x+y‖22=‖x‖22+‖y‖22.
= xTx+2yTx+yTy 而xTy=yTx=0, 因此 ‖x+y‖22=‖x‖22+‖y‖22 注:推广到Rn中的向量组α1,α2,…,αk, 如果αiTαj=0 (i≠j), 称α1,α2,…,αk是 正交向量组. 特别地: 如果‖αi‖2=1(i=1,2,…,k), 即 αiTαj=δij,称α1,α2,…,αk 为标准正交向量组.
51
设U是Rn中的子空间, x∈Rn. 如果x与U中任意向量正交, 称向量x与子空间U正交, 记为x⊥U.
设U,V是Rn中两个子空间, 如果任意x∈U和任意y∈V是正交的, 称子空间U与子空间V正交, 记为U⊥V. 设U,V是Rn中互补的子空间. 如果U⊥V, 那么称U,V互为正交补子空间, 记U=V⊥或V=U⊥. 可以证明, 一个子空间的正交补子空间是惟一的.
52
定理6.3.3 设A是n×k阶矩阵,x∈Rn, 那么下列三种情况是
等价的: ①x⊥R(A); ②ATx=0; ③x∈N(AT). 这里,N(AT)={ATx=0, x∈Rn}称为AT的核子空间. 证:由N(AT)的定义, ②与③显然等价. 下面证明①与②等价. 记A=(α1,α2,…,αk), 那么,αi∈R(A) (i=1,2,…,k). 假设x⊥R(A), 即αiTx=0 (i=1,2,…,k). 从而ATx=0. 另一方面,如果ATx=0, 那么有z∈Rk, 使Az=y∈R(A). 这时,yTx=zTATx=0,即x⊥y. 由z的任意性, 得Az是任意的, 因此x⊥R(A). 由这个定理, 容易得到: 推论1设A是n×k阶矩阵, 那么R(A)有惟一的正交补子 空间N(AT).
53
矩阵的QR分解 定理6.3.4 设A=(α1,α2,…,αn)是列满秩矩阵,αi∈Rm (i=1,2,…,n)且m≥n. 那么, A有惟一的QR分解,记为 A=QR, (4) 这里,Q是有n个标准正交列的m×n阶矩阵,R是有正对角元的n阶上三角矩阵. 证:由A是列满秩矩阵可知, ATA是n阶正定矩阵, 因此有 惟一的Cholesky分解: ATA=RTR, (5) 这里R是有正对角元的上三角矩阵, R-1存在. 令Q=AR-1, (6) 那么, QTQ=R-TATAR-1=In, 即Q是有标准正交列的m×n阶矩阵. 由(6)式, (4)式成立, 且由(5)式的惟一性, 分解式(4)也是惟一的.
54
6.3.3 Householder矩阵与矩阵的正交三角化
定义1设w是欧氏空间Rn中的单位向量, 形如H=I-2wwT (10) 的n阶矩阵称为Householder矩阵, 也称 为反射(镜像)矩阵或称为Householder变 换(反射变换、镜像变换).
55
定理6.3.5 设H是Rn中的Householder
矩阵, 那么, ①HT=H (对称性); ②H-1=HT (正交性); ③H2=I (对合性). 证:直接验证, 性质①成立. 由H=I-2wwT, wTw=1和性质①, H2=HTH=(I-2wwT)T (I-2wwT) =I-4wwTwwT+4wwT=I, 即HT=H-1,性质②和③成立.
56
定理6.3.6设Rn中有非零向量x≠y 且‖x‖2=‖y‖2,那么,存在 Householder矩阵H,使Hx=y 证:不妨令w=(x-y)/ ‖x-y‖2 ,有wTw=1.设H=I- 2wwT ,那么 Hx=(I-2wwT) 由已知, xTx=yTy. 所以 (x-y)T(x-y)=xTx- 2yTx+yTy=2(xTx-yTx). 从而, Hx=x-(x-y)=y. 证明中如果令u=x-y, 那么Householder矩阵形为(11)式,即H=I-β-1uuT,β=uTu. 同样可验证Hx=y成立.
57
例1 设R3中向量x=(0,4,3)T和向量σe1,其中σ2=‖x‖22,即σ= =±5
例1 设R3中向量x=(0,4,3)T和向量σe1,其中σ2=‖x‖22,即σ= =±5. 按(11)式构造一个3阶Householder矩阵H,使Hx=σe1. 解:令u=x-σe1.设σ=-5, 于是 u=(-5, 4, 3)T 又由uTu= =50, 所以β=25. 由H=I-β-1uuT, 得
58
容易验证Hx=σe1. 在Rn中,如果已知非零向量x=(x1,x2,…,xn)T和向量σe1≠x, 其中,σ=±‖x‖2, 即‖σe1‖2=‖x‖2 由定理6.3.6, 存在一个形如(11)式或(10)式的 Householder矩阵H,使Hx=σe1. 由定理6.3.6, 令u=x-σe (13) 又由σ2=‖x‖22=xTx和(12)式, β=1/2× (x-σe1)T(x-σe1)= 1/2× (xTx-2σx1+σ2) =σ(σ-x1) (14) 记H=I-β-1uuT,容易验证Hx=σe1.
Similar presentations