第七章函数逼近与曲线(面)拟合

合集下载

Ch3函数逼近与曲线拟合省名师优质课赛课获奖课件市赛课一等奖课件

定理若 f ( x) C[a, b], 则总存在Pn*( x) H , 使得 || f ( x) Pn* ( x) || En
定义（偏差点）
若 f ( x) C[a, b], P( x) Hn , 若在x x0上有
|
P( x0 )
f ( x0 ) |
max |
a xb
P(x)
内积与内积空间
定义设X 是数域K上的线性空间，对u, v X , 存在K中的一个数(u, v)与之对应，满足（1）(u, v) (v, u)
(2) ( u, v) (u, v), K
(3) (u v, w) (u, w) (v, w) (4) (u, u) 0,(u, u) 0 u 0 称( , ) 0为X 上的内积，定义了内积的线性空间称为内积空间.
n ( x)
1 2n1
Tn ( x)与零的偏差最小，其偏差为
1 2n1
.
证明
n ( x)
1 2n1
Tn ( x)
xn
P* n1
(
x
)
max
1 x1
|
n
(
x
)
|
1 2n1
max
1 x1
|
Tn
(
x
)
|
1 2n1
, 又知道
k
xk cos n (k 0,1, 2,
, n)是切比雪夫交错点组，
由此知道Pn*1 ( x)是xn的区间[- 1,1]上的最佳逼近多项式，
1 2n n!
dn dx n
{( x2
1)n }
勒让德多项式旳性质
1. 正交性
0,
1 －1
Pn

常用函数的逼近和曲线拟合

常用函数的逼近和曲线拟合在数学中，函数逼近和曲线拟合都是常见的问题。

函数逼近是指找到一个已知函数，尽可能地接近另一个函数。

而曲线拟合则是给定一组数据点，找到一条曲线来描述这些数据点的分布。

本文将讨论常用的函数逼近和曲线拟合方法。

一、函数逼近1. 插值法插值法是最简单的函数逼近方法之一。

它的基本思想是：给定一组已知点，通过构造一个多项式，使得该多项式在这些点处的函数值与已知函数值相等。

插值法的优点是精度高，缺点是易产生龙格现象。

常用的插值多项式有拉格朗日插值多项式和牛顿插值多项式。

拉格朗日插值多项式的形式为：$f(x)=\sum_{i=0}^{n}y_{i}\prod_{j=i,j\neq i}^{n}\frac{x-x_{j}}{x_{i}-x_{j}}$其中，$x_{i}$是已知点的横坐标，$y_{i}$是已知点的纵坐标，$n$是已知点的数量。

牛顿插值多项式的形式为：$f(x)=\sum_{i=0}^{n}f[x_{0},x_{1},...,x_{i}]\prod_{j=0}^{i-1}(x-x_{j})$其中，$f[x_{0},x_{1},...,x_{i}]$是已知点$(x_{0},y_{0}),(x_{1},y_{1}),...,(x_{i},y_{i})$的差商。

2. 最小二乘法最小二乘法是一种常用的函数逼近方法。

它的基本思想是：给定一组数据点，找到一个函数，在这些数据点上的误差平方和最小。

通常采用线性模型，例如多项式模型、指数模型等。

最小二乘法的优点是适用性广泛，缺点是对于非线性模型要求比较高。

最小二乘法的一般形式为：$F(x)=\sum_{i=0}^{n}a_{i}\varphi_{i}(x)$其中，$a_{i}$是待求的系数，$\varphi_{i}(x)$是一组已知的基函数，$n$是基函数的数量。

最小二乘法的目标是使得$\sum_{i=1}^{m}[f(x_{i})-F(x_{i})]^{2}$最小，其中$m$是数据点的数量。

逼近和拟合专题教育课件

则在[a,b]上g( x) 0；
就称( x)为[a,b]上的权函数.
例2 设f ( x), g( x) C[a,b], ( x)为[a,b]上的权函数 ,则可
定义内积
( f , g) ab( x) f ( x)g( x)dx. 1,( f , g) ab f (x)g(x)dx.
容易验证内积定义中的四个性质，并导出范数
项式正交.
（4）有递推关系
pn1( x) ( x n ) pn( x) n pn1( x), n 0,1,, (2.4)
其中 p0( x) 1，p1( x) 0，
n ( xpn , pn ) /( pn , pn ), n ( pn , pn ) /( pn1, pn1)，n 1,2,,
称为Gram矩阵，则G非奇异的充要条件是 u1, u2,, un线性
无关.
证明：1) G非奇异以G为系数矩阵的齐次线性方程组
n
n
( ju j , uk ) (u j , uk ) j 0,
j1
j1
只有零解。
k 1,,n.
n
n
n
2) juj 0 ( juj , juj ) 0
有限维空间 vs 无限维空间.
Rn, C[a,b],
定理 1(维尔斯特拉斯 ) 如果f ( x) C[a,b], 那么 0,
多项式p( x),使得
| f ( x) p( x) | , 对于一切a x b.
伯恩斯坦(1912)给出一种构造性证明：伯恩斯坦多项式
Bn (
f
,
x)
n
k0
f
(n ,n )
根据定理 3，0,,n线性无关 det(G) 0.
§2 正交多项式

工科研究生“数值分析”课程教学大纲及教学日历

工科研究生“数值分析”课程的教学大纲序号：课程编号：课程名称：数值分析/ Numerical Analysis学时：40 学分： 2.5责任教师：王开荣，何光辉，董海云，李东，温罗生适用专业：工科研究生各专业先修课程：高等数学、线性代数课程教材：《应用数值分析》，王开荣，杨大地，高等教育出版社，2010年7月参考教材：1. 关治, 陆金甫，《数值方法》清华大学出版社，2006.2.2. Numerical Analysis Using MATLAB,Fourth Edition,电子工业出版社(影印版)，2005年7月。

一、课程的性质、目的和任务学习数值分析课程能培养学生运用数学的方法和借助计算机解决工程计算问题的能力。

其任务是通过近似计算，使得许多难以求解的数学问题得以简化、可行。

并得到满足误差要求的近似解。

本课程的目的和任务是使工科研究生掌握工程应用中的数值计算方法，为具有不同工程背景的学生能运用这些近似计算方法处理在工程技术及其科学研究中出现的计算问题奠定坚实的基础。

通过学习要求学生能正确理解数值分析的所有的概念和算法，掌握算法的构造思想及其基本算法的步骤。

能应用工具软件Matlab独立完成常用的算法的编程及数值计算。

通过典型的数值算例验证所编程序的正确性，并且应用到实际问题中。

二、课程的教学内容和基本要求1．误差(4学时)(1)了解误差的来源和误差的概念;(2)理解误差的传播和算法中应避免的问题;2．线性方程组的直接解法(6学时)(1)掌握Guass消去法，理解范数的概念;(2)熟练运用Gauss列主元素法，三角分解法，追赶法;3．线性方程组的迭代法(4学时)(1)理解迭代法的收敛条件，掌握Jacobi迭代法;(2)熟练运用Seidel，SOR迭代法;4．方阵的特征值与特征向量的计算(2学时)(1)了解QR方法;(2)熟练运用乘幂法和反幂法，Jacobi方法;5．非线性方程求根(4学时)(1)掌握二分法;(2)熟练使用Newton法;6．插值法(6时)(1)掌握Lagrange插值，Newton插值，Hermite插值；(2)熟练运用分段插值，样条插值；7．函数逼近与数据拟合(2时)(1)掌握多项式逼近，拟合；(2)熟练运用正交多项式逼近，拟合;8．数值积分(6时)(1)掌握Newton-Cotes公式，Gauss求积公式；(2) 熟练运用Romberg积分公式，复化Gauss型公式;9．常微分方程初值问题的数值解法(4时)(1)掌握Euler方法，Runge-kutta方法，Admas预测-校正法；(2)了解稳定性、收敛性和计算误差估计，高阶方程及方程组.10．总复习(2时)四、考试方式考试以笔试、闭卷的方式进行。

小波分析之函数逼近与曲线拟合

(ϕ 0 , ϕ 0 ) (ϕ 0 , ϕ1 ) ⋯ (ϕ 0 , ϕ n ) (ϕ1 , ϕ 0 ) (ϕ1 , ϕ1 ) ⋯ (ϕ1 , ϕ n ) G = ⋮ ⋮ ⋮ (ϕ , ϕ ) (ϕ , ϕ ) ⋯ (ϕ , ϕ ) n 1 n n n 0
∑α
j =1
n
j
u j = 0.
Th3的证明
G非奇异当且仅当齐次方程组
∑ (u
j =1
n
j
, u k ) α j = 0, k = 1, 2, ⋯ , n.
只有零解,即 ( ∑ α j u j , u k ) = 0, k = 1, 2, ⋯ , n.
j =1
n
只有零解,即
∑α
j =1
n
j
u j = 0 只有零解,即
C[a.b]带权内积的定义
• [a,b]上的非负函数 ρ (x) 称为[a,b]上的权函数，若满足： b ① ∫ x k ρ ( x)dx 存在且为有限值(k=0,1,…). a ② 对[a,b]上的非负连续函数g (x) ，若
∫
b
a
ρ ( x ) g ( x ) dx = 0 ,
则 g ( x) ≡ 0
2
= =
∑
n
x
2 i
i = 1
∞
m a x
1 ≤ i ≤ n
x
i
n与C[a,b]上范数的扩充关系 R
= = =
• 向量范数：范数：范数
x x x
1
∑
n
x
i = 1
i
2
∑
n∞
max
x
1 ≤ i ≤ n

函数逼近与曲线(面)拟合(matlab)

7.1 曲线拟合、误差及其MATLAB 程序例7.1.1 已知函数)2(sin 7145)(23x x x x f y π+-==和一组数据),(i i y x 列入表7–1中，比较最大误差，平均误差，均方根误差和误差平方和.表7–1 例7.1.1的一组数据),(y x解由给定的函数和数据，在MATLAB 工作窗口输入>> x=[-2.5,-1.7,-1.1,-0.8,0,0.1,0.5,3.6]; n=length(x);y=[-43.50 5.69 11.34 14.16 0 1.02 -6.37 185.84];f=5.*x.^3-14.*x+7.*(sin(2*pi*x)).^2; fy=abs(f-y);fy2=fy.^2; [x',y',f',fy',fy2'], Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n), E=sum(fy2)运行后屏幕显示如下x y f fy fy2-2.5000 -43.5000 -43.1250 0.3750 0.1406-1.7000 5.6900 5.5666 0.1234 0.0152-1.1000 11.3400 11.1634 0.1766 0.0312-0.8000 14.1600 14.9716 0.8116 0.65860 0 0 0 00.1000 1.0200 1.0234 0.0034 0.00000.5000 -6.3700 -6.3750 0.0050 0.00003.6000 185.8400 185.2984 0.5416 0.2933Ew = E1 = E2 = E =0.8116 0.2546 0.3773 1.13907.2 曲线拟合的线性最小二乘法及其MATLAB 程序例7.2.1 给出一组数据点),(i i y x 列入表7–2中，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.表7–2 例7.2.1的一组数据),(y x解（1）在MATLAB 工作窗口输入程序>> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('例7.2.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序>> syms a1 a2 a3 a4x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];fi=a1.*x.^3+ a2.*x.^2+ a3.*x+ a4运行后屏幕显示关于a 1,a 2, a 3和a 4的线性方程组fi =[ -125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4,a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4]编写构造误差平方和的MATLAB 程序>> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];fi=[-125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4,5832/125*a1+324/25*a2+18/5*a3+a4];fy=fi-y; fy2=fy.^2; J=sum(fy.^2)运行后屏幕显示误差平方和如下J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a 2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2为求4321,,,a a a a 使J 达到最小，只需利用极值的必要条件0=∂∂ka J )4,3,2,1(=k ，得到关于4321,,,a a a a 的线性方程组，这可以由下面的MA TLAB 程序完成，即输入程序>> syms a1 a2 a3 a4J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4...+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a 4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2;Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4);Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3), Ja41=simple(Ja4),运行后屏幕显示J 分别对a 1, a 2 ,a 3 ,a 4的偏导数如下Ja11=56918107/10000*a1+32097579/25000*a2+1377283/2500*a3+23667/250*a4-8442429/625Ja21 =32097579/25000*a1+1377283/2500*a2+23667/250*a3+67*a4+767319/625Ja31 =1377283/2500*a1+23667/250*a2+67*a3+18/5*a4-232638/125Ja41 =23667/250*a1+67*a2+18/5*a3+18*a4+14859/25解线性方程组Ja 11 =0，Ja 21 =0，Ja 31 =0，Ja 41 =0，输入下列程序>>A=[56918107/10000, 32097579/25000, 1377283/2500, 23667/250; 32097579/25000, 1377283/2500, 23667/250, 67; 1377283/2500, 23667/250, 67, 18/5; 23667/250, 67, 18/5, 18];B=[8442429/625, -767319/625, 232638/125, -14859/25];C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f 及其系数C 如下C = 5.0911 -14.1905 6.4102 -8.2574f=716503695845759/140737488355328*x^3-7988544102557579/562949953421312*x^2+1804307491277693/281474976710656*x-4648521160813215/562949953421312故所求的拟合曲线为8.25746.410214.19055.0911)(23-+-=x x x x f .（4）编写下面的MATLAB 程序估计其误差，并作出拟合曲线和数据的图形.输入程序>> xi=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];n=length(xi);f=5.0911.*xi.^3-14.1905.*xi.^2+6.4102.*xi -8.2574;x=-2.5:0.01: 3.6;F=5.0911.*x.^3-14.1905.*x.^2+6.4102.*x -8.2574;fy=abs(f-y); fy2=fy.^2; Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n)plot(xi,y,'r*'), hold on, plot(x,F,'b-'), hold offlegend('数据点(xi,yi)','拟合曲线y=f(x)'),xlabel('x'), ylabel('y'),title('例7.2.1的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w ，平均误差E 1和均方根误差E 2及其数据点),(i i y x 和拟合曲线y =f (x )的图形（略）.Ew = E1 = E2 =3.105 4 0.903 4 1.240 97.3 函数)(x r k 的选取及其MATLAB 程序例7.3.1 给出一组实验数据点),(i i y x 的横坐标向量为x =(-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5, -2.1,-1.5, -2.7,-3.6)，纵横坐标向量为y =(459.26,52.81,198.27,165.60,59.17,41.66,25.92, 22.37,13.47, 12.87, 11.87,6.69,14.87,24.22)，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.解（1）在MATLAB 工作窗口输入程序>>x=[-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5,-2.1,-1.5, -2.7,-3.6];y=[459.26,52.81,198.27,165.60,59.17,41.66,25.92,22.37,13.47, 12.87, 11.87,6.69,14.87,24.22];plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('例7.3.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序>> syms a bx=[-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5,-2.1,-1.5,-2.7,-3.6]; fi=a.*exp(-b.*x)运行后屏幕显示关于a 和b 的线性方程组fi =[ a*exp(17/2*b), a*exp(87/10*b), a*exp(71/10*b),a*exp(34/5*b), a*exp(51/10*b), a*exp(9/2*b), a*exp(18/5*b), a*exp(17/5*b), a*exp(13/5*b), a*exp(5/2*b), a*exp(21/10*b), a*exp(3/2*b), a*exp(27/10*b), a*exp(18/5*b)]编写构造误差平方和的MATLAB 程序如下>>y=[459.26,52.81,198.27,165.60,59.17,41.66,25.92,22.37,13.47,12.87, 11.87, 6.69,14.87,24.22];fi =[ a*exp(17/2*b), a*exp(87/10*b), a*exp(71/10*b), a*exp(34/5*b), a*exp(51/10*b), a*exp(9/2*b), a*exp(18/5*b), a*exp(17/5*b), a*exp(13/5*b), a*exp(5/2*b), a*exp(21/10*b), a*exp(3/2*b), a*exp(27/10*b), a*exp(18/5*b)];fy=fi-y;fy2=fy.^2;J=sum(fy.^2)运行后屏幕显示误差平方和如下J =(a*exp(17/2*b)-22963/50)^2+(a*exp(87/10*b)-5281/100)^2+(a*exp(71/10*b)-19827/100)^2+(a*exp(34/5*b)-828/5)^2+(a*exp(51/10*b)-5917/100)^2+(a*exp(9/2*b)-2083/50)^2+(a*exp(18/5*b)-648/25)^2+(a*exp(17/5*b)-2237/100)^2+(a*exp(13/5*b)-1347/100)^2+(a*ex p(5/2*b)-1287/100)^2+(a*exp(21/10*b)-1187/100)^2+(a*exp(3/2*b)-669/100)^2+(a*exp(27/10*b)-1487/100)^2+(a*exp(18/5*b)-1211/50)^2为求b a ,使J 达到最小，只需利用极值的必要条件，得到关于b a ,的线性方程组，这可以由下面的MA TLAB 程序完成，即输入程序>> syms a bJ=(a*exp(17/2*b)-22963/50)^2+(a*exp(87/10*b)-5281/100)^2+(a*exp(71/10*b)-19827/100)^2+(a*exp(34/5*b)-828/5)^2+(a*exp(51/10*b)-5917/100)^2+(a*exp(9/2*b)-2083/50)^2+(a*exp(18/5*b)-648/25)^2+(a*exp(17/5*b)-2237/100)^2+(a*exp(13/5*b)-1347/100)^2+(a*exp(5/2*b)-1287/100)^2+(a*exp(21/10*b)-1187/100)^2+(a*exp(3/2*b )-669/100)^2+(a*exp(27/10*b)-1487/100)^2+(a*exp(18/5*b)-1211/50)^2;Ja=diff(J,a); Jb=diff(J,b);Ja1=simple(Ja), Jb1=simple(Jb),运行后屏幕显示J 分别对b a ,的偏导数如下Ja1 =2*a*exp(3*b)+2*a*exp(17*b)+2*a*exp(87/5*b)+2*exp(68/5*b)*a+2*exp(9*b)*a+2*a*exp(34/5*b)-669/50*exp(3/2*b)-1487/50*exp(27/10*b)-2507/25*exp(18/5*b)-22963/25*exp(17/2*b)-5281/50*exp(87/10*b)-19827/50*exp(71/10*b)-2237/50*exp(17/5*b)-1656/5*exp(34/5*b)-1347/50*exp(13/5*b)-5917/50*exp(51/10*b)-1287/50*exp(5/2*b )-2083/25*exp(9/2*b)-1187/50*exp(21/10*b)+4*a*exp(36/5*b)+2*a*e xp(26/5*b)+2*a*exp(71/5*b)+2*a*exp(51/5*b)+2*a*exp(5*b)+2*a*exp (21/5*b)+2*a*exp(27/5*b)Jb1 =1/500*a*(2100*a*exp(21/10*b)^2+8500*a*exp(17/2*b)^2+6800*a*exp(34/5*b)^2-10035*exp(3/2*b)-40149*exp(27/10*b)-180504*exp (18/5*b)-3903710*exp(17/2*b)-459447*exp(87/10*b)-1407717*exp(71/10*b)-76058*exp(17/5*b)-1126080*exp(34/5*b)-35022*exp(13/5*b)-301767*exp(51/10*b)-32175*exp(5/2*b)-187470*exp(9/2*b)-24927*ex p(21/10*b)+7100*a*exp(71/10*b)^2+5100*a*exp(51/10*b)^2+4500*a*exp(9/2*b)^2+7200*a*exp(18/5*b)^2+3400*a*exp(17/5*b)^2+2600*a*ex p(13/5*b)^2+2500*a*exp(5/2*b)^2+1500*a*exp(3/2*b)^2+2700*a*exp(27/10*b)^2+8700*a*exp(87/10*b)^2)用解二元非线性方程组的牛顿法的MATLAB 程序求解线性方程组J a1 =0，J b1 =0，得a = b=2.811 0 0.581 6故所求的拟合曲线（7.13）为0811.2)(=x f e x 5816.0-. （7.14）（4）根据（7.2），（7.3），（7.4）和（7.14）式编写下面的MA TLAB 程序估计其误差，并做出拟合曲线和数据的图形.输入程序>> xi=[-8.5 -8.7 -7.1 -6.8 -5.10 -4.5 -3.6 -3.4 -2.6 -2.5-2.1 -1.5 -2.7 -3.6];y=[459.26 52.81 198.27 165.60 59.17 41.66 25.92 22.3713.47 12.87 11.87 6.69 14.87 24.22];n=length(xi); f=2.8110.*exp(-0.5816.*xi); x=-9:0.01: -1;F=2.8110.*exp(-0.5816.*x); fy=abs(f-y); fy2=fy.^2;Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n), plot(xi,y,'r*'), hold on plot(x,F,'b-'), hold off,legend('数据点(xi,yi)','拟合曲线y=f(x)')xlabel('x'), ylabel('y'),title('例7.3.1的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w = 390.141 5，平均误差E 1=36.942 2和均方根误差E 2=106.031 7及其数据点),(i i y x 和拟合曲线y =f (x )的图形（略）.7.4 多项式拟合及其MATLAB 程序例7.4.1 给出一组数据点),(i i y x 列入表7–3中，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.表7–3 例7.4.1的一组数据),(y x解（1）首先根据表7–3给出的数据点i i ，用下列MATLAB 程序画出散点图. 在MATLAB 工作窗口输入程序>> x=[-2.9 -1.9 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[53.94 33.68 20.88 16.92 8.79 8.98 4.17 9.1219.88];plot(x,y,'r*'), legend('数据点(xi,yi)')xlabel('x'), ylabel('y'),title('例7.4.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）用作线性最小二乘拟合的多项式拟合的MATLAB 程序求待定系数k a )3,2,1(=k .输入程序>> a=polyfit(x,y,2)运行后输出（7.16）式的系数a =2.8302 -7.3721 9.1382故拟合多项式为2138.91372.72830.2)(2+-=x x x f .（4）编写下面的MATLAB 程序估计其误差，并做出拟合曲线和数据的图形.输入程序>> xi=[-2.9 -1.9 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[53.94 33.68 20.88 16.92 8.79 8.98 4.17 9.12 19.88];n=length(xi); f=2.8302.*xi.^2-7.3721.*xi+9.1382x=-2.9:0.001:3.6;F=2.8302.*x.^2-7.3721.*x+8.79;fy=abs(f-y); fy2=fy.^2; Ew=max(fy), E1=sum(fy)/n,E2=sqrt((sum(fy2))/n), plot(xi,y,'r*', x,F,'b-'),legend('数据点(xi,yi)','拟合曲线y=f(x)')xlabel('x'), ylabel('y'),title('例7.4.1 的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w ，平均误差E1和均方根误差E 2及其数据点(x i ,y i )和拟合曲线y =f (x )的图形（略）.Ew = E1 = E2 =0.745 7, 0.389 2, 0.436 37.5 拟合曲线的线性变换及其MATLAB 程序例7.5.1 给出一组实验数据点),(i i y x 的横坐标向量为x =（7.5 6.8 5.10 4.53.6 3.4 2.6 2.5 2.1 1.5 2.7 3.6），纵横坐标向量为y =（359.26 165.60 59.17 41.66 25.92 22.37 13.47 12.87 11.87 6.69 14.87 24.22），试用线性变换和线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.解（1）首先根据给出的数据点),(i i y x ，用下列MATLAB 程序画出散点图.在MATLAB 工作窗口输入程序>> x=[7.5 6.8 5.10 4.5 3.6 3.4 2.6 2.5 2.1 1.5 2.73.6];y=[359.26 165.60 59.17 41.66 25.92 22.37 13.47 12.87 11.87 6.69 14.87 24.22];plot(x,y,'r*'), legend('数据点(xi,yi)')xlabel('x'), ylabel('y'),title('例7.5.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（2）根据数据散点图，取拟合曲线为a y =e bx )0,0(≠>b a , （7.19）其中b a ,是待定系数.令b B a A y Y ===,ln ,ln ，则（7.19）化为Bx A Y +=.在MATLAB 工作窗口输入程序>> x=[7.5 6.8 5.10 4.5 3.6 3.4 2.6 2.5 2.1 1.5 2.73.6];y=[359.26 165.60 59.17 41.66 25.92 22.37 13.47 12.87 11.87 6.69 14.87 24.22];Y=log(y); a=polyfit(x,Y,1); B=a(1);A=a(2); b=B,a=exp(A)n=length(x); X=8:-0.01:1; Y=a*exp(b.*X); f=a*exp(b.*x);plot(x,y,'r*',X,Y,'b-'), xlabel('x'),ylabel('y')legend('数据点(xi,yi)','拟合曲线y=f(x)')title('例7.5.1 的数据点(xi,yi)和拟合曲线y=f(x)的图形')fy=abs(f-y); fy2=fy.^2; Ew=max(fy), E1=sum(fy)/n,E2=sqrt((sum(fy2))/n)运行后屏幕显示a y =e bx 的系数b =0.624 1，a =2.703 9,数据),(i i y x 与拟合函数f的最大误差Ew =67.641 9，平均误差E 1=8.677 6和均方根误差E 2=20.711 3及其数据点),(i i y x 和拟合曲线9703.2)(=x f e x 1624.0的图形（略）.7.6 函数逼近及其MATLAB 程序最佳均方逼近的MATLAB 主程序function [yy1,a,WE]=zjjfbj(f,X,Y,xx)m=size(f);n=length(X);m=m(1);b=zeros(m,m); c=zeros(m,1);if n~=length(Y)error('X 和Y 的维数应该相同')endfor j=1:mfor k=1:mb(j,k)=0;for i=1:nb(j,k)=b(j,k)+feval(f(j,:),X(i))*feval(f(k,:),X(i));endendc(j)=0;for i=1:nc(j)=c(j)+feval(f(j,:),X(i))*Y(i);endenda=b\c;WE=0;for i=1:nff=0;for j=1:mff=ff+a(j)*feval(f(j,:),X(i));endWE=WE+(Y(i)-ff)*(Y(i)-ff);endif nargin==3return ;endyy=[];for i=1:ml=[];for j=1:length(xx)l=[l,feval(f(i,:),xx(j))];endyy=[yy l'];endyy=yy*a; yy1=yy'; a=a';WE;例7.6.1 对数据X 和Y , 用函数2,,1x y x y y ===进行逼近，用所得到的逼近函数计算在 6.5=x 处的函数值，并估计误差.其中X =(1 3 4 5 6 7 8 9); Y =(-11 -13 -11 -7 -1 7 17 29).解在MATLAB 工作窗口输入程序>> X=[ 1 3 4 5 6 7 8 9]; Y=[-11 -13 -11 -7 -1 7 17 29];f=['fun0';'fun1';'fun2']; [yy,a,WE]=zjjfbj(f,X,Y,6.5)运行后屏幕显示如下yy =2.75000000000003a =-7.00000000000010 -4.99999999999995 1.00000000000000WE =7.172323350269439e-027例7.6.2 对数据X 和Y ，用函数2,,1x y x y y ===，x y cos =，=y e x，x y sin =进行逼近，其中X =（0 0.50 1.00 1.50 2.00 2.50 3.00），Y =（0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645）.解在MATLAB 工作窗口输入程序>> X=[ 0 0.50 1.00 1.50 2.00 2.50 3.00];Y=[0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645];f=['fun0';'fun1';'fun2';'fun3';'fun4';'fun5'];xx=0:0.2:3;[yy,a,WE]=zjjfbj(f,X,Y, xx), plot(X,Y,'ro',xx,yy,'b-') 运行后屏幕显示如下（图略）yy = Columns 1 through 7-0.0005 0.2037 0.3939 0.5656 0.7141 0.8348 0.9236Columns 8 through 140.9771 0.9926 0.9691 0.9069 0.8080 0.6766 0.5191Columns 15 through 160.3444 0.1642a = 0.3828 0.4070 -0.3901 0.0765 -0.4598 0.5653 WE = 1.5769e-004即，最佳逼近函数为y=0.3828+0.4070*x -0.3901*x^2+0.0765*exp(x) -0.4598*cos(x) +0.5653*sin(x).7.7 三角多项式逼近及其MATLAB 程序计算三角多项式的MATLAB 主程序function [A,B,Y1,Rm]=sanjiao(X,Y,X1,m)n= length(X)-1;max1=fix((n-1)/2);if m > max1m=max1;endA=zeros(1,m+1);B=zeros(1,m+1);Ym=(Y(1)+Y(n+1))/2; Y(1)=Ym; Y(n+1)=Ym; A(1)=2*sum(Y)/n; for i=1:mB(i+1)=sin(i*X)*Y'; A(i+1)=cos(i*X)*Y';endA=2*A/n; B=2*B/n; A(1)=A(1)/2;Y1=A(1);for k=1:mY1=Y1+A(k+1)*cos(k*X1)+ B(k+1)*sin(k*X1);Tm=A(1)+A(k+1).*cos(k*X)+ B(k+1).*sin(k*X); k=k+1;end Y;Tm; Rm=(sum(Y-Tm).^2)/n;例7.7.1 根据],[ππ-上的350,60,13=n 个等距横坐标点n i x i π+π-=2 ),,2,1,0(n i =和函数3sin 2)(x x f =. （1）求)(x f 的6阶三角多项式逼近，计算均方误差;（2）将这三个三角多项式分别与)(x f 的傅里叶级数nx n n x f n n sin 19)1(318)(121∑∞=+--π= 的前6项进行比较;（3）利用三角多项式分别计算X i = -2, 2.5的值;（4）在同一坐标系中，画出函数)(x f ，350,60,13=n 的三角多项式和数据点的图形.解（1）输入程序>> X1=-pi:2*pi/13:pi;Y1=2*sin(X1/3);X1i=[-2,2.5];[A1,B1,Y11,Rm1]=sanjiao(X1,Y1,X1i,6),X2=-pi:2*pi/60:pi;Y2=2*sin(X2/3);[A2,B2,Y12,Rm2]=sanjiao(X2,Y2,X1i,6)X3=-pi:2*pi/350:pi;Y3=2*sin(X3/3);[A3,B3,Y13,Rm3]=sanjiao(X3,Y3,X1i,6)X1i=[-2,2.5];Y1=2*sin(X1i/3)for n=1:6bi=(-1)^(n+1)*18*sqrt(3)*n/(pi*(9*n^2-1))end（2）画图，输入程序>>X1=-pi:2*pi/13:pi;Y1=2*sin(X1/3);Xi=-pi:0.001:pi; f=2*sin(Xi/3);[A1,B1,Y1i,R1m]=sanjiao(X1,Y1,Xi,6);X2=-pi:2*pi/60:pi; Y2=2*sin(X2/3); X3=-pi:2*pi/350:pi;Y3=2*sin(X3/3);[A2,B2,Y2i,R2m]=sanjiao(X2,Y2,Xi,6);[A3,B3,Y3i,R3m]=sanjiao(X3,Y3,Xi,6);plot(X1,Y1,'r*', Xi, Y1i,'b-',Xi, Y2i,'g--', Xi, Y3i, 'm:',Xi, f, 'k-.')xlabel('x'),ylabel('y')legend('数据点(xi,yi)','n=13的三角多项式','n=60的三角多项式','n=350的三角多项式','函数f(x)')title('例7.7.1 的数据点(xi,yi)、n=13，60，350的三角多项式T3和函数f(x)的图形')运行后图形（略）.7.8 随机数据点上的二元拟合及其MATLAB 程序例7.8.1 设节点(X,Y,Z )中的X 和Y 分别是在区间]3,3[-和]5.3,5.2[-上的50个随机数，Z 是函数Z =7-3x 3e 22y - -x 在(X,Y )的值，拟合点(X I ,Y I )中的X I =-3:0.2:3, Y I =-2.5:0.2:3.5.分别用二元拟合方法中最近邻内插法、三角基线性内插法、三角基三次内插法和MATLAB 4网格化坐标方法计算在(X I ,Y I )处的值，作出它们的图形，并与被拟和曲面进行比较.解（1）最近邻内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'nearest') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用最近邻内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用最近邻内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面及其插值Z I （略）.（2）三角基线性内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'linear') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基线性内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用三角基线性内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I （略）.（3）三角基三次内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'cubic') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基三次内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用三角基三次内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I （略）.（4）MATLAB 4网格化坐标方法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z 的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'v4') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用MATLAB 4网格化坐标方法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用MATLAB 4网格化坐标方法拟合函数Z =7-3x 3e 22y - -x在两组不同节点处的曲面和节点的图形及其插值ZI （略）.（5）作被拟合曲面Z =7-3x 3e 22y - -x 和节点的图形.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z 的值.X1=-3.:0.1:3.;Y1=-2.5:0.1:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=7-3* XI.^3 .* exp(-XI.^2 - YI.^2);mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('被拟合函数曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示被拟合函数Z =7-3x 3e 22y - -x 的曲面和节点的图形及其函数值ZI （略）.7.9 随机数据点上的n 元拟合及其MATLAB 程序例7.9.1 首先利用MATLAB 函数rand 产生随机数据X 1,Y 1,Z 1,然后用线性变换b at u += (其中5,5-==b a )将随机数据X 1,Y 1 ,Z 1 变换为节点坐标(X,Y,Z ),再用函数)1(373+-=z y x w e222z y x ---生成数据W , 用三元最近邻内插法方法计算函数w 在插值点x i ,10:5.0:3-= y i ,13:5.0:2-=z i =y i 处拟合数据的值，并作其图形.解输入程序>> X1=-5+5*rand(10,1);Y1=-5+5*rand(10,1);Z1=Y1;[X,Y,Z] = meshgrid(X1,Y1,Z1);W=7-3* X.^3 .* Y.*(Z+1).* exp(-X.^2 - Y.^2- Z.^2);xi=-3:0.5:10;yi=-2:0.5:13;zi=yi;[XI,YI,ZI] = meshgrid(xi,yi,zi);W1=griddata3(X, Y, Z, W, XI, YI, ZI, 'nearest');slice(XI,YI,ZI,W1,[-2 4 9.5],9,[-2 2 9]),%shading flat%lighting flatxlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数W=7-3X^3Y(Z+1)exp(-X^2 - Y^2- Z^2) '); hold oncolorbar('horiz')view([-30 45])运行后屏幕显示三元线性拟合值及其图形（略）.例7.9.2 设节点(X,Y,Z,W )中的X,Y 和Z 分别是在区间]3,3[-和]5.3,5.2[-,Y =Z上的15个随机数，W 是函数x w +=2e 222z y x ---在(X,Y,Z )的值，拟合点(x i ,y i ,z i )中的x i =-3:0.2:3， y i =-2.5:0.2:3.5，z i =y i ，用'linear '方法计算拟合数据的值，并作其图形.解输入程序>> x=rand(15,1); y=rand(15,1);X1=-3+(3-(-3))*x;Y1=-2.5+(3.5-(-2.5))*y;Z1=Y1;[X,Y,Z] = meshgrid(X1,Y1,Z1);W=2+X.* exp(-X.^2 - Y.^2- Z.^2);xi=-3:0.2:3; yi=-2.5:0.2:3.5; zi=yi;[X2,Y2,Z2]=meshgrid(xi,yi,zi);W1=griddata3(X, Y, Z, W, X2,Y2,Z2,'linear');slice(X2,Y2,Z2,W1,[-1 0 1.5],2,[-2 3]),shading flat,lighting flat,xlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数W=2+X exp(-X^2 - Y^2- Z^2)');hold on,colorbar('horiz'), view([-3 5])运行后屏幕显示三元线性拟合值及其图形（略）.。

函数逼近泛函

函数逼近1. 函数的定义在数学中，函数逼近是一种通过使用一组已知函数来近似描述一个未知函数的方法。

函数逼近的目标是找到一个或多个已知函数，使其在某个范围内与未知函数的值尽可能接近。

2. 函数逼近的用途函数逼近在许多领域中都有广泛的应用，特别是在数值计算和数据分析中。

以下是几个常见的应用场景：2.1 插值插值是一种通过已知数据点之间构建一个连续函数来估计未知数据点的方法。

函数逼近可以用于选择合适的插值函数，并通过最小化插值误差来提高插值精度。

2.2 曲线拟合曲线拟合是一种通过找到一个或多个已知函数，使其与给定数据点之间的差异最小化来估计未知曲线的方法。

函数逼近可以用于选择最佳拟合函数，并通过调整参数来优化拟合结果。

2.3 数据压缩对于大规模数据集，使用较少数量的已知函数来表示整个数据集可以有效地进行数据压缩。

函数逼近可以将原始数据转换为更紧凑且易于存储的表示形式，同时保持数据的关键特征。

2.4 数据平滑函数逼近可以用于平滑噪声数据，通过选择适当的平滑函数来减小数据中的不规则性和噪声。

这对于信号处理和图像处理等领域尤为重要。

3. 函数逼近的工作方式函数逼近的工作方式可以分为两个主要步骤：选择适当的已知函数和确定最佳拟合参数。

3.1 选择已知函数在函数逼近中，首先需要选择一组已知函数。

这些已知函数可以是多项式、三角函数、指数函数或其他常见的数学函数。

选择已知函数时，需要考虑未知函数的特点和拟合需求。

3.2 确定最佳拟合参数确定最佳拟合参数是通过最小化误差来实现的。

常见的误差度量方法包括均方误差（Mean Squared Error）和最大误差（Maximum Error）。

通过调整已知函数中的参数，使其与未知函数在给定数据点上的值之间的误差最小化。

优化算法通常用于找到最佳拟合参数。

常见的优化算法包括梯度下降法、牛顿法和Levenberg-Marquardt算法等。

这些算法通过迭代的方式调整参数，直到达到最小化误差的目标。

计算方法-函数拟合

7
7
拟合曲线的误差平方和：Q ei2 ((xi ) yi )2 3.09524
i 1
i 1
7
S yy ( yi y)2 64 i 1
R 1 Q 0.9755 Syy
15
➢ 多元线性拟合
除了一个变量的拟合外，还有多个变量的线性拟合函数，
(x1, x2,, xm ) a0 a1x1 a2x2 amxm
常见做法
使 max | (xi ) yi | 最小 1im
太复杂
······ ·· ···· ( x) ？ (x) ?
m
使 | (xi ) yi | 最小
不可导，求解困难
i1
m
使 | (xi ) yi |2 最小
i1
最小二乘法：目前最好的多项式曲线拟合算法
按误差平方和达到极小构造拟合函数的方法称为最佳平方逼近。这一最小准则称为最小二乘原理。用最小二乘原理估计拟合函数待定参数由高斯(1975)提出。
6
写成方程组形式
n
n
n
n
n
( xi0 )？a0 ( xi1)？a1 ( xi2 )？a2 ( xim )？am xi0 yi
i1
i 1
i 1
i 1
i 1
n
n
n
n
n
( xi1)？a0 ( xi2 )？a1 ( xi3 )？a2 ( xim1)a？m xi1 yi
i 1
①
②
③
12
n
2 ( yi i (xi ))(i (xi ) y)
i (xi ) a0 a1xi1 amxim
i 1 n
2 ( yi a0 a1xi1 amxim )(a0 a1xi1 amxim y)

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

7.1 曲线拟合、误差及其MATLAB 程序例7.1.1 已知函数)2(sin 7145)(23x x x x f y π+-==和一组数据),(i i y x 列入表7–1中，比较最大误差，平均误差，均方根误差和误差平方和.表7–1 例7.1.1的一组数据),(y x解由给定的函数和数据，在MATLAB 工作窗口输入>> x=[-2.5,-1.7,-1.1,-0.8,0,0.1,0.5,3.6]; n=length(x);y=[-43.50 5.69 11.34 14.16 0 1.02 -6.37 185.84];f=5.*x.^3-14.*x+7.*(sin(2*pi*x)).^2; fy=abs(f-y);fy2=fy.^2; [x',y',f',fy',fy2'], Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n), E=sum(fy2)运行后屏幕显示如下x y f fy fy2-2.5000 -43.5000 -43.1250 0.3750 0.1406-1.7000 5.6900 5.5666 0.1234 0.0152-1.1000 11.3400 11.1634 0.1766 0.0312-0.8000 14.1600 14.9716 0.8116 0.65860 0 0 0 00.1000 1.0200 1.0234 0.0034 0.00000.5000 -6.3700 -6.3750 0.0050 0.00003.6000 185.8400 185.2984 0.5416 0.2933Ew = E1 = E2 = E =0.8116 0.2546 0.3773 1.13907.2 曲线拟合的线性最小二乘法及其MATLAB 程序例7.2.1 给出一组数据点),(i i y x 列入表7–2中，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.表7–2 例7.2.1的一组数据),(y x解（1）在MATLAB 工作窗口输入程序>> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('例7.2.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序>> syms a1 a2 a3 a4x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];fi=a1.*x.^3+ a2.*x.^2+ a3.*x+ a4运行后屏幕显示关于a 1,a 2, a 3和a 4的线性方程组fi =[ -125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4,a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4]编写构造误差平方和的MA TLAB 程序>> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.5068.04];fi=[-125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4,5832/125*a1+324/25*a2+18/5*a3+a4];fy=fi-y; fy2=fy.^2; J=sum(fy.^2)运行后屏幕显示误差平方和如下J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a 2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2为求4321,,,a a a a 使J 达到最小，只需利用极值的必要条件0=∂∂k a J )4,3,2,1(=k ，得到关于4321,,,a a a a 的线性方程组，这可以由下面的MA TLAB 程序完成，即输入程序>> syms a1 a2 a3 a4J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4...+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a 4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2;Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4);Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3), Ja41=simple(Ja4),运行后屏幕显示J 分别对a 1, a 2 ,a 3 ,a 4的偏导数如下Ja11=56918107/10000*a1+32097579/25000*a2+1377283/2500*a3+23667/250*a4-8442429/625Ja21 =32097579/25000*a1+1377283/2500*a2+23667/250*a3+67*a4+767319/625Ja31 =1377283/2500*a1+23667/250*a2+67*a3+18/5*a4-232638/125Ja41 =23667/250*a1+67*a2+18/5*a3+18*a4+14859/25解线性方程组Ja 11 =0，Ja 21 =0，Ja 31 =0，Ja 41 =0，输入下列程序>>A=[56918107/10000, 32097579/25000, 1377283/2500, 23667/250; 32097579/25000, 1377283/2500, 23667/250, 67; 1377283/2500, 23667/250, 67, 18/5; 23667/250, 67, 18/5, 18];B=[8442429/625, -767319/625, 232638/125, -14859/25];C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f 及其系数C 如下C = 5.0911 -14.1905 6.4102 -8.2574f=716503695845759/140737488355328*x^3-7988544102557579/562949953421312*x^2+1804307491277693/281474976710656*x-4648521160813215/562949953421312故所求的拟合曲线为8.25746.410214.19055.0911)(23-+-=x x x x f . （4）编写下面的MA TLAB 程序估计其误差，并作出拟合曲线和数据的图形.输入程序>> xi=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];n=length(xi);f=5.0911.*xi.^3-14.1905.*xi.^2+6.4102.*xi -8.2574;x=-2.5:0.01: 3.6;F=5.0911.*x.^3-14.1905.*x.^2+6.4102.*x -8.2574;fy=abs(f-y); fy2=fy.^2; Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n)plot(xi,y,'r*'), hold on, plot(x,F,'b-'), hold offlegend('数据点(xi,yi)','拟合曲线y=f(x)'),xlabel('x'), ylabel('y'),title('例7.2.1的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w ，平均误差E 1和均方根误差E 2及其数据点),(i i y x 和拟合曲线y =f (x )的图形（略）.Ew = E1 = E2 =3.105 4 0.903 4 1.240 97.3 函数)(x r k 的选取及其MATLAB 程序例7.3.1 给出一组实验数据点),(i i y x 的横坐标向量为x =(-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5, -2.1,-1.5, -2.7,-3.6)，纵横坐标向量为y =(459.26,52.81,198.27,165.60,59.17,41.66,25.92, 22.37,13.47, 12.87, 11.87,6.69,14.87,24.22)，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.解（1）在MATLAB 工作窗口输入程序>>x=[-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5,-2.1,-1.5, -2.7,-3.6];y=[459.26,52.81,198.27,165.60,59.17,41.66,25.92,22.37,13.47, 12.87, 11.87,6.69,14.87,24.22];plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('例7.3.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序>> syms a bx=[-8.5,-8.7,-7.1,-6.8,-5.10,-4.5,-3.6,-3.4,-2.6,-2.5,-2.1,-1.5,-2.7,-3.6]; fi=a.*exp(-b.*x)运行后屏幕显示关于a和b的线性方程组fi =[ a*exp(17/2*b), a*exp(87/10*b), a*exp(71/10*b),a*exp(34/5*b), a*exp(51/10*b), a*exp(9/2*b), a*exp(18/5*b),a*exp(17/5*b), a*exp(13/5*b), a*exp(5/2*b), a*exp(21/10*b),a*exp(3/2*b), a*exp(27/10*b), a*exp(18/5*b)]编写构造误差平方和的MA TLAB程序如下>>y=[459.26,52.81,198.27,165.60,59.17,41.66,25.92,22.37,13.47,12.87, 11.87, 6.69,14.87,24.22];fi =[ a*exp(17/2*b), a*exp(87/10*b), a*exp(71/10*b),a*exp(34/5*b), a*exp(51/10*b), a*exp(9/2*b),a*exp(18/5*b), a*exp(17/5*b), a*exp(13/5*b),a*exp(5/2*b), a*exp(21/10*b), a*exp(3/2*b),a*exp(27/10*b), a*exp(18/5*b)];fy=fi-y;fy2=fy.^2;J=sum(fy.^2)运行后屏幕显示误差平方和如下J =(a*exp(17/2*b)-22963/50)^2+(a*exp(87/10*b)-5281/100)^2+(a*exp(71/10*b)-19827/100)^2+(a*exp(34/5*b)-828/5)^2+(a*exp(51/10*b)-5917/100)^2+(a*exp(9/2*b)-2083/50)^2+(a*exp(18/5*b)-648/25)^2+(a*exp(17/5*b)-2237/100)^2+(a*exp(13/5*b)-1347/100)^2+(a*exp(5/2*b)-1287/100)^2+(a*exp(21/10*b)-1187/100)^2+(a*exp(3/2*b)-669/100)^2+(a*exp(27/10*b)-1487/100)^2+(a*exp(18/5*b)-1211/50)^2为求ba,的线性方程组，a,使J达到最小，只需利用极值的必要条件，得到关于b这可以由下面的MA TLAB程序完成，即输入程序>> syms a bJ=(a*exp(17/2*b)-22963/50)^2+(a*exp(87/10*b)-5281/100)^2+(a*exp(71/10*b)-19827/100)^2+(a*exp(34/5*b)-828/5)^2+(a*exp(51/10*b)-5917/100)^2+(a*exp(9/2*b)-2083/50)^2+(a*exp(18/5*b)-648/ 25)^2+(a*exp(17/5*b)-2237/100)^2+(a*exp(13/5*b)-1347/100)^2+(a*exp(5/2*b)-1287/100)^2+(a*exp(21/10*b)-1187/100)^2+(a*exp(3/2*b)-669/100)^2+(a*exp(27/10*b)-1487/100)^2+(a*exp(18/5*b)-1211/50)^2;Ja=diff(J,a); Jb=diff(J,b);Ja1=simple(Ja), Jb1=simple(Jb),运行后屏幕显示J分别对ba,的偏导数如下Ja1 =2*a*exp(3*b)+2*a*exp(17*b)+2*a*exp(87/5*b)+2*exp(68/5*b)*a+2*exp(9*b)*a+2*a*exp(34/5*b)-669/50*exp(3/2*b)-1487/50*exp(27/10*b)-2507/25*exp(18/5*b)-22963/25*exp(17/2*b)-5281/50*exp(87/10*b)-19827/50*exp(71/10*b)-2237/50*exp(17/5*b)-1656/5*exp(34/5*b)-1347/50*exp(13/5*b)-5917/50*exp(51/10*b)-1287/50*exp(5/2*b)-2083/25*exp(9/2*b)-1187/50*exp(21/10*b)+4*a*exp(36/5*b)+2*a*exp(26/5*b)+2*a*exp(71/5*b)+2*a*exp(51/5*b)+2*a*exp(5*b)+2*a*exp (21/5*b)+2*a*exp(27/5*b)Jb1 =1/500*a*(2100*a*exp(21/10*b)^2+8500*a*exp(17/2*b)^2+6800*a*exp(34/5*b)^2-10035*exp(3/2*b)-40149*exp(27/10*b)-180504*exp (18/5*b)-3903710*exp(17/2*b)-459447*exp(87/10*b)-1407717*exp(71/10*b)-76058*exp(17/5*b)-1126080*exp(34/5*b)-35022*exp(13/5*b)-301767*exp(51/10*b)-32175*exp(5/2*b)-187470*exp(9/2*b)-24927*exp(21/10*b)+7100*a*exp(71/10*b)^2+5100*a*exp(51/10*b)^2+4500*a*exp(9/2*b)^2+7200*a*exp(18/5*b)^2+3400*a*exp(17/5*b)^2+2600*a*ex p(13/5*b)^2+2500*a*exp(5/2*b)^2+1500*a*exp(3/2*b)^2+2700*a*exp(27/10*b)^2+8700*a*exp(87/10*b)^2)用解二元非线性方程组的牛顿法的MATLAB 程序求解线性方程组J a1 =0，J b1 =0，得a = b=2.811 0 0.581 6故所求的拟合曲线（7.13）为0811.2)(=x f e x5816.0-. （7.14）（4）根据（7.2），（7.3），（7.4）和（7.14）式编写下面的MA TLAB 程序估计其误差，并做出拟合曲线和数据的图形.输入程序>> xi=[-8.5 -8.7 -7.1 -6.8 -5.10 -4.5 -3.6 -3.4 -2.6 -2.5-2.1 -1.5 -2.7 -3.6];y=[459.26 52.81 198.27 165.60 59.17 41.66 25.92 22.3713.47 12.87 11.87 6.69 14.87 24.22];n=length(xi); f=2.8110.*exp(-0.5816.*xi); x=-9:0.01: -1;F=2.8110.*exp(-0.5816.*x); fy=abs(f-y); fy2=fy.^2;Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n), plot(xi,y,'r*'), hold on plot(x,F,'b-'), hold off,legend('数据点(xi,yi)','拟合曲线y=f(x)')xlabel('x'), ylabel('y'),title('例7.3.1的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w = 390.141 5，平均误差E 1=36.942 2和均方根误差E 2=106.031 7及其数据点),(i i y x 和拟合曲线y =f (x )的图形（略）.7.4 多项式拟合及其MATLAB 程序例7.4.1 给出一组数据点),(i i y x 列入表7–3中，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.表7–3 例7.4.1的一组数据),(y x解（1）首先根据表7–3给出的数据点),(i i y x ，用下列MATLAB 程序画出散点图.在MATLAB 工作窗口输入程序>> x=[-2.9 -1.9 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[53.94 33.68 20.88 16.92 8.79 8.98 4.17 9.1219.88];plot(x,y,'r*'), legend('数据点(xi,yi)')xlabel('x'), ylabel('y'),title('例7.4.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（3）用作线性最小二乘拟合的多项式拟合的MATLAB 程序求待定系数k a )3,2,1(=k .输入程序>> a=polyfit(x,y,2)运行后输出（7.16）式的系数a =2.8302 -7.3721 9.1382故拟合多项式为2138.91372.72830.2)(2+-=x x x f . （4）编写下面的MATLAB 程序估计其误差，并做出拟合曲线和数据的图形.输入程序>> xi=[-2.9 -1.9 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[53.94 33.68 20.88 16.92 8.79 8.98 4.17 9.12 19.88];n=length(xi); f=2.8302.*xi.^2-7.3721.*xi+9.1382x=-2.9:0.001:3.6;F=2.8302.*x.^2-7.3721.*x+8.79;fy=abs(f-y); fy2=fy.^2; Ew=max(fy), E1=sum(fy)/n,E2=sqrt((sum(fy2))/n), plot(xi,y,'r*', x,F,'b-'),legend('数据点(xi,yi)','拟合曲线y=f(x)')xlabel('x'), ylabel('y'),title('例7.4.1 的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w ，平均误差E1和均方根误差E 2及其数据点(x i ,y i )和拟合曲线y =f (x )的图形（略）.Ew = E1 = E2 =0.745 7, 0.389 2, 0.436 37.5 拟合曲线的线性变换及其MATLAB 程序例7.5.1 给出一组实验数据点),(i i y x 的横坐标向量为x =（7.5 6.8 5.10 4.53.6 3.4 2.6 2.5 2.1 1.5 2.7 3.6），纵横坐标向量为y =（359.26 165.60 59.17 41.66 25.92 22.37 13.47 12.87 11.87 6.69 14.87 24.22），试用线性变换和线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线.解（1）首先根据给出的数据点),(i i y x ，用下列MATLAB 程序画出散点图.在MATLAB 工作窗口输入程序>> x=[7.5 6.8 5.10 4.5 3.6 3.4 2.6 2.5 2.1 1.5 2.73.6];y=[359.26 165.60 59.17 41.66 25.92 22.37 13.47 12.8711.87 6.69 14.87 24.22];plot(x,y,'r*'), legend('数据点(xi,yi)')xlabel('x'), ylabel('y'),title('例7.5.1的数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图（略）.（2）根据数据散点图，取拟合曲线为a y =e bx )0,0(≠>b a , （7.19）其中b a ,是待定系数.令b B a A y Y ===,ln ,ln ，则（7.19）化为Bx A Y +=.在MATLAB 工作窗口输入程序>> x=[7.5 6.8 5.10 4.5 3.6 3.4 2.6 2.5 2.1 1.5 2.73.6];y=[359.26 165.60 59.17 41.66 25.92 22.37 13.4712.87 11.87 6.69 14.87 24.22];Y=log(y); a=polyfit(x,Y,1); B=a(1);A=a(2); b=B,a=exp(A)n=length(x); X=8:-0.01:1; Y=a*exp(b.*X); f=a*exp(b.*x);plot(x,y,'r*',X,Y,'b-'), xlabel('x'),ylabel('y')legend('数据点(xi,yi)','拟合曲线y=f(x)')title('例7.5.1 的数据点(xi,yi)和拟合曲线y=f(x)的图形')fy=abs(f-y); fy2=fy.^2; Ew=max(fy), E1=sum(fy)/n,E2=sqrt((sum(fy2))/n)运行后屏幕显示a y =e bx 的系数b =0.624 1，a =2.703 9,数据),(i i y x 与拟合函数f的最大误差Ew =67.641 9，平均误差E 1=8.677 6和均方根误差E 2=20.711 3及其数据点),(i i y x 和拟合曲线9703.2)(=x f e x 1624.0的图形（略）.7.6 函数逼近及其MATLAB 程序最佳均方逼近的MATLAB 主程序function [yy1,a,WE]=zjjfbj(f,X,Y,xx)m=size(f);n=length(X);m=m(1);b=zeros(m,m); c=zeros(m,1);if n~=length(Y)error('X 和Y 的维数应该相同')endfor j=1:mfor k=1:mb(j,k)=0;for i=1:nb(j,k)=b(j,k)+feval(f(j,:),X(i))*feval(f(k,:),X(i));endendc(j)=0;for i=1:nc(j)=c(j)+feval(f(j,:),X(i))*Y(i);endenda=b\c;WE=0;for i=1:nff=0;for j=1:mff=ff+a(j)*feval(f(j,:),X(i));endWE=WE+(Y(i)-ff)*(Y(i)-ff);endif nargin==3return ;endyy=[];for i=1:ml=[];for j=1:length(xx)l=[l,feval(f(i,:),xx(j))];endyy=[yy l'];endyy=yy*a; yy1=yy'; a=a';WE;例7.6.1 对数据X 和Y , 用函数2,,1x y x y y ===进行逼近，用所得到的逼近函数计算在 6.5=x 处的函数值，并估计误差.其中X =(1 3 4 5 6 7 8 9); Y =(-11 -13 -11 -7 -1 7 17 29).解在MATLAB 工作窗口输入程序>> X=[ 1 3 4 5 6 7 8 9]; Y=[-11 -13 -11 -7 -1 7 17 29];f=['fun0';'fun1';'fun2']; [yy,a,WE]=zjjfbj(f,X,Y,6.5)运行后屏幕显示如下yy =2.75000000000003a =-7.00000000000010 -4.99999999999995 1.00000000000000WE =7.172323350269439e-027例7.6.2 对数据X 和Y ，用函数2,,1x y x y y ===，x y cos =，=y e x ，xy sin =进行逼近，其中X =（0 0.50 1.00 1.50 2.00 2.50 3.00），Y =（0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645）.解在MATLAB 工作窗口输入程序>> X=[ 0 0.50 1.00 1.50 2.00 2.50 3.00];Y=[0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645];f=['fun0';'fun1';'fun2';'fun3';'fun4';'fun5'];xx=0:0.2:3;[yy,a,WE]=zjjfbj(f,X,Y, xx), plot(X,Y,'ro',xx,yy,'b-') 运行后屏幕显示如下（图略）yy = Columns 1 through 7-0.0005 0.2037 0.3939 0.5656 0.7141 0.83480.9236Columns 8 through 140.9771 0.9926 0.9691 0.9069 0.8080 0.67660.5191Columns 15 through 160.3444 0.1642a = 0.3828 0.4070 -0.3901 0.0765 -0.4598 0.5653 WE = 1.5769e-004即，最佳逼近函数为y=0.3828+0.4070*x-0.3901*x^2+0.0765*exp(x) -0.4598*cos(x) +0.5653*sin(x).7.7 三角多项式逼近及其MATLAB 程序计算三角多项式的MATLAB 主程序function [A,B,Y1,Rm]=sanjiao(X,Y,X1,m)n= length(X)-1;max1=fix((n-1)/2);if m > max1m=max1;endA=zeros(1,m+1);B=zeros(1,m+1);Ym=(Y(1)+Y(n+1))/2; Y(1)=Ym; Y(n+1)=Ym; A(1)=2*sum(Y)/n; for i=1:mB(i+1)=sin(i*X)*Y'; A(i+1)=cos(i*X)*Y';endA=2*A/n; B=2*B/n; A(1)=A(1)/2;Y1=A(1);for k=1:m Y1=Y1+A(k+1)*cos(k*X1)+ B(k+1)*sin(k*X1);Tm=A(1)+A(k+1).*cos(k*X)+ B(k+1).*sin(k*X); k=k+1;endY;Tm; Rm=(sum(Y-Tm).^2)/n;例7.7.1 根据],[ππ-上的350,60,13=n 个等距横坐标点n i x i π+π-=2 ),,2,1,0(n i =和函数3sin 2)(xx f =. （1）求)(x f 的6阶三角多项式逼近，计算均方误差;（2）将这三个三角多项式分别与)(x f 的傅里叶级数nx n n x f n n sin 19)1(318)(121∑∞=+--π=的前6项进行比较;（3）利用三角多项式分别计算X i = -2, 2.5的值;（4）在同一坐标系中，画出函数)(x f ，350,60,13=n 的三角多项式和数据点的图形.解（1）输入程序>> X1=-pi:2*pi/13:pi;Y1=2*sin(X1/3);X1i=[-2,2.5];[A1,B1,Y11,Rm1]=sanjiao(X1,Y1,X1i,6),X2=-pi:2*pi/60:pi;Y2=2*sin(X2/3);[A2,B2,Y12,Rm2]=sanjiao(X2,Y2,X1i,6)X3=-pi:2*pi/350:pi;Y3=2*sin(X3/3);[A3,B3,Y13,Rm3]=sanjiao(X3,Y3,X1i,6)X1i=[-2,2.5];Y1=2*sin(X1i/3)for n=1:6bi=(-1)^(n+1)*18*sqrt(3)*n/(pi*(9*n^2-1))end（2）画图，输入程序>>X1=-pi:2*pi/13:pi;Y1=2*sin(X1/3);Xi=-pi:0.001:pi; f=2*sin(Xi/3);[A1,B1,Y1i,R1m]=sanjiao(X1,Y1,Xi,6);X2=-pi:2*pi/60:pi; Y2=2*sin(X2/3); X3=-pi:2*pi/350:pi;Y3=2*sin(X3/3);[A2,B2,Y2i,R2m]=sanjiao(X2,Y2,Xi,6);[A3,B3,Y3i,R3m]=sanjiao(X3,Y3,Xi,6);plot(X1,Y1,'r*', Xi, Y1i,'b-',Xi, Y2i,'g--', Xi, Y3i, 'm:',Xi, f, 'k-.')xlabel('x'),ylabel('y')legend('数据点(xi,yi)','n=13的三角多项式','n=60的三角多项式','n=350的三角多项式','函数f(x)')title('例7.7.1 的数据点(xi,yi)、n=13，60，350的三角多项式T3和函数f(x)的图形')运行后图形（略）.7.8 随机数据点上的二元拟合及其MATLAB 程序例7.8.1 设节点(X,Y,Z )中的X 和Y 分别是在区间]3,3[-和]5.3,5.2[-上的50个随机数，Z 是函数Z =7-3x 3e 22y - -x 在(X,Y )的值，拟合点(X I ,Y I )中的X I =-3:0.2:3, Y I =-2.5:0.2:3.5.分别用二元拟合方法中最近邻内插法、三角基线性内插法、三角基三次内插法和MATLAB 4网格化坐标方法计算在(X I ,Y I )处的值，作出它们的图形，并与被拟和曲面进行比较.解（1）最近邻内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'nearest') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用最近邻内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用最近邻内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面及其插值Z I （略）.（2）三角基线性内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'linear') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基线性内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用三角基线性内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I （略）.（3）三角基三次内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'cubic') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基三次内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用三角基三次内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I （略）.（4）MATLAB 4网格化坐标方法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=griddata(X,Y,Z,XI,YI, 'v4') %计算在每个插值点（XI,YI ）处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用MATLAB 4网格化坐标方法拟合函数z =7-3 x^3 exp(-x^2 - y^2)的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示用MATLAB 4网格化坐标方法拟合函数Z =7-3x 3e 22y - -x在两组不同节点处的曲面和节点的图形及其插值ZI （略）.（5）作被拟合曲面Z =7-3x 3e 22y - -x 和节点的图形.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y ， x ，y .X=-3+(3-(-3))*x;%利用x 生成随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点（X,Y ）处计算Z的值.X1=-3.:0.1:3.;Y1=-2.5:0.1:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标（XI,YI ）网格化.ZI=7-3* XI.^3 .* exp(-XI.^2 - YI.^2);mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('被拟合函数曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z). hold of %结束在当前图形上添加新图形.运行后屏幕显示被拟合函数Z =7-3x 3e 22y - -x 的曲面和节点的图形及其函数值ZI （略）.7.9 随机数据点上的n 元拟合及其MATLAB 程序例7.9.1 首先利用MATLAB 函数rand 产生随机数据X 1,Y 1,Z 1,然后用线性变换b at u += (其中5,5-==b a )将随机数据X 1,Y 1 ,Z 1 变换为节点坐标(X,Y,Z ),再用函数)1(373+-=z y x w e 222z y x ---生成数据W , 用三元最近邻内插法方法计算函数w 在插值点x i ,10:5.0:3-= y i ,13:5.0:2-=z i =y i 处拟合数据的值，并作其图形.解输入程序>> X1=-5+5*rand(10,1);Y1=-5+5*rand(10,1);Z1=Y1;[X,Y,Z] = meshgrid(X1,Y1,Z1);W=7-3* X.^3 .* Y.*(Z+1).* exp(-X.^2 - Y.^2- Z.^2);xi=-3:0.5:10;yi=-2:0.5:13;zi=yi;[XI,YI,ZI] = meshgrid(xi,yi,zi);W1=griddata3(X, Y, Z, W, XI, YI, ZI, 'nearest');slice(XI,YI,ZI,W1,[-2 4 9.5],9,[-2 2 9]),%shading flat%lighting flatxlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数W=7-3X^3Y(Z+1)exp(-X^2 - Y^2- Z^2) '); hold oncolorbar('horiz')view([-30 45])运行后屏幕显示三元线性拟合值及其图形（略）.例7.9.2 设节点(X,Y,Z,W )中的X,Y 和Z 分别是在区间]3,3[-和]5.3,5.2[-,Y =Z 上的15个随机数，W 是函数x w +=2e 222z y x ---在(X,Y,Z )的值，拟合点(x i ,y i ,z i )中的x i =-3:0.2:3， y i =-2.5:0.2:3.5，z i =y i ，用'linear '方法计算拟合数据的值，并作其图形.解输入程序>> x=rand(15,1); y=rand(15,1);X1=-3+(3-(-3))*x;Y1=-2.5+(3.5-(-2.5))*y;Z1=Y1;[X,Y,Z] = meshgrid(X1,Y1,Z1);W=2+X.* exp(-X.^2 - Y.^2- Z.^2);xi=-3:0.2:3; yi=-2.5:0.2:3.5; zi=yi;[X2,Y2,Z2]=meshgrid(xi,yi,zi);W1=griddata3(X, Y, Z, W, X2,Y2,Z2,'linear');slice(X2,Y2,Z2,W1,[-1 0 1.5],2,[-2 3]),shading flat,lighting flat,xlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数W=2+X exp(-X^2 - Y^2- Z^2)');hold on,colorbar('horiz'), view([-3 5])运行后屏幕显示三元线性拟合值及其图形（略）.。