数学建模～插值与拟合(课件ppt)

1、插值与拟合一、插值的基本原理二、拟合的基本原理三、插值与拟合的关系四、插值的MATLAB实现五、拟合的Matlab实现 我们经常会遇到大量的数据需要处理，而处理数据的关键就在于这些算法，例如数据拟合、参数估计、插值等数据处理算法。此类问题在MATLAB中有很多现成的函数可以调用，熟悉MATLAB，这些方法都能游刃有余的用好。 一、概述 数据拟合在很多赛题中有应用，与图形处理有关的问题很多与插值和拟合有关系，例如98年美国赛A题，生物组织切片的三维插值处理，94年A题逢山开路，山体海拔高度的插值计算，2003年吵的沸沸扬扬的“非典”问题也要用到数据拟合算法，观察数据的走向进行处理， 2005年的

2、雨量预报的评价的插值计算。2001年的公交车调度拟合问题，2003年的饮酒驾车拟合问题。插值问题雨量预报的评价预测点和实测点的图形插值后的图形拟合问题饮酒驾车喝两瓶酒的拟合曲线喝1-5瓶酒的拟合曲线 在实际中，常常要处理由实验或测量所得到的一些离散数据。插值与拟合方法就是要通过这些数据去确定某一类已知函数的参数或寻求某个近似函数，使所得到的近似函数与已知数据有较高的拟合精度。 如果要求这个近似函数（曲线或曲面）经过所已知的所有数据点，则称此类问题为插值问题。 （不需要函数表达式）二、基本概念 如果不要求近似函数通过所有数据点，而是要求它能较好地反映数据变化规律的近似函数的方法称为数据拟合。（必

3、须有函数表达式） 近似函数不一定（曲线或曲面）通过所有的数据点。 1、联系都是根据实际中一组已知数据来构造一个能够反映数据变化规律的近似函数的方法。2、区别插值问题不一定得到近似函数的表达形式，仅通过插值方法找到未知点对应的值。数据拟合要求得到一个具体的近似函数的表达式。 三、插值与拟合的区别和联系9 区别：运算过程上的区别：拟合：是将数据点用最恰当的曲线描述出来，以反映问题的规律，是特殊到一般的过程。插值：是在知道曲线的形状后得出某些具体点的性质的过程,是从一般到特殊。求解误差上的区别：拟合：考虑观察值的误差（误差不可避免时）。以偏差的某种最小为拟合标准插值：在观测点处无误差如最小二乘法：联

4、系：二者都是函数逼近的主要方法引言2-插值和拟合的联系与区别四、插值的使用及求解 当数据量不够，需要补充，且认定已有数据可信时, 通常利用函数插值方法。 实际问题当中碰到的函数 f (x) 是各种各样的，有的表达式很复杂，有的甚至给不出数学的式子，只提供了一些离散数据，警如，某些点上的函数值和导数值。 4.1 引言 选用不同类型的插值函数，逼近的效果就不同，一般有：（1）拉格朗日插值（lagrange插值）（2）分段线性插值（3）Hermite（4）三次样条插值。 4.2 插值方法数据插值与拟合在工程实践与科学实验中,常常需要从一组试验数据之中找到自变量与因变量之间的关系，一般可用一个近似函数

5、表示。函数产生的办法因观测数据的要求不同而异，数据插值与拟合是两种常用的方法。4.1 MATLAB中的插值函数4.2 拉格朗日插值法4.3 利用均差的牛顿插值法4.4 利用差分的牛顿插值法4.5 Hermite插值4.6 spline三次样条插值4.7 多项式曲线拟合4.8 最小二乘拟合4.1 MATLAB中的插值函数函数插值来源于函数的以下问题：只知道函数在某区间有定义且已得到区间内一些离散点的值，希望用简单的表达式近似给出函数在此区间上的整体描述，并能与已知离散点上的值相等。插值法按插值函数的形式主要分为以下几种形式： （1）代数多项式插值； （2）三角多项式插值； （3）有理分式插值。

6、代数多项式插值是最常用的插值方式，其内容也是最丰富的，它又可分为以下几种插值方式： （1）非等距节点插值，包括拉格朗日插值、利用均差的牛顿插值和埃特金插值； （2）非等距节点插值，包括利用差分的牛顿插值和高斯插值等； （3）在插值中增加了导数的Hermite（埃尔米特）插值； （4）分段插值，包括分段线性插值、分段Hermite（埃尔米特）插值和样条函数插值；（5）反插值。按被插值函数的变量个数还可把插值法分为一元插值和多元插值。154.1 函数插值 在化工领域中，通过实验可获得有限个离散点的数据表（如表4-1所示）。 这种用表格形式给出的函数 通常称为列表函数 表4-1 列表函数 0 导入1

7、6 如果可以寻找出与已测得的实验数据相适应的近似解析函数式，则可根据近似解析表达式求出未列点的函数值。 列表函数虽然可以对实验数据点的函数变化规律有一定的反映，但它不能给出数据点外的函数值，因此使用起来往往很不方便。170.1 函数插值的定义图41 n次插值多项式的几何表示插值就是定义一个在特定点取给定值的函数的过程。当选用不同的插值函数时，相应的得到不同的插值方法。180.2 代数多项式插值法 为给实验数据构造相适应的近似函数表达式，最常用的是代数多项式插值法。 代数多项式插值法的基本思想: 设法构造某个简单函数 （ 是代数多项式形式，称插值函数）作为 的近似表达式，然后计算 未列点的自变量

8、 的函数值 值作为 的近似值。19对于多项式插值，我们主要讨论:插值多项式的常用构造方法有哪些？根据化工计算的要求，本节着重介绍拉格朗日（Lagrange）多项式插值法、分段线性插值法和三次样条插值法。20 4.1.1 拉格朗日多项式插值法1 插值多项式模型 已知函数y=f(x)在n个点 上的值f( ),求一个项数低于n的插值多项式（k=1,2,n）表示插值多项式的最紧凑的方式是拉格朗日（Lagrange）形式（40）21（41）（42）22拉格朗日插值法求多项式模型为：（44）23例1 求二次插值多项式。 解 : 按拉格朗日方法，有：242 拉格朗日多项式插值法的Matlab实现对于拉格朗日

9、多项式插值法，Matlab中没有专用指令，可根据拉格朗日多项式插值法的数学公式来编程实现。程序如下:25例4.1： 从手册中查到水在不同温度（T）时的导热系数（ ）的数据如表4-2所示。 表4-2 水在不同温度时的导热系数 试用拉格朗日多项式插值法分别计算水在80C、250C、560C、730C、960C时热导率数值。T/0C /W/(mK)T/0C /W/(mK)00.553600.659200.599800.675400.6341000.68326解：Matlab计算程序如下： T0=0,20,40,60,80,100;K0=0.553,0.599,0.634,0.659,0.675,0.

10、683;T=8,25,56,73,96;K=lagrange(T0,K0,T);T;K; 执行结果：ans = 8.0000 25.0000 56.0000 73.0000 96.0000 0.5728 0.6087 0.6548 0.6704 0.6820这里:T0和K0分别表示样本点温度和导热系数，T和K分别表示插值点温度和导热系数。27 Lagrange插值多项式的优点在于不要求数据点是等间隔的，其缺点是数据点数不宜过大，通常不超过7个，否则计算工作量大且误差大，计算不稳定。6 分段低次插值 /* piecewise polynomial approximation */例：在5, 5上

11、考察 的Ln(x)。取 -5 -4 -3 -2 -1 0 1 2 3 4 5 -0.5 0 0.5 1 1.5 2 2.5 n 越大，端点附近抖动越大，称为Runge 现象Ln(x) f (x)分段低次插值29 原因： 当 n 很大时，数据的误差可能对插值多项式计算结果带来很大误差，引起“龙格现象” ，也可以说是拉格朗日插值法的不稳定现象。 为避免“龙格现象” ，通常限定 n x=0:2*pi; y=sin(x); xx=0:0.5:2*pi yy=interp1(x,y,xx); plot(x,y,s,xx,yy)注：例4-1中用默认的(分段线性插值的linear)对已知的7个sin函数的数

12、据点进行插值，用plot画出插值结果。从图中可以看出分段线性就是联结两个邻近的已知点的线性函数插值计算该区间内插值点上的函数值。例4-2 用其他一维插值方法对以下7个离散数据点(1,3.5)、(2,2.1)、(3,1.3)、(4.0.8)、(5,2.9)、(6,4.2)、(7,5.7)进行一维插值方法。解：在MATLAB命令窗口中输入以下命令： x=1 2 3 4 5 6 7; y=3.5 2.1 1.3 0.8 2.9 4.2 5.7; xx=1:0.5:7; y1=interp1(x,y,xx,nearest); y2=interp1(x,y,xx,spline); y3=interp1(

13、x,y,xx,cubic); plot(x,y,o,xx,y1,-,xx,y2,-.,xx,y3,:)41例4.2： 已知某转子流量计在1001000mL/min流量范围内，刻度值与校正值有如下关系（见表4-3）。试用线性插值法计算流量计的刻度值为785时，实际流量为多少？刻度值校正值刻度值校正值100105.3600605.8200207.2700707.4300308.1800806.7400406.9900908.0500507.51000107.9表4-3 例4.2数据42解：Matlab计算程序如下：X=100,200,300,400,500,600,700,800,900,1000

14、;Y=105.3,207.2,308.1,406.9,507.5,605.8,707.4,806.7,908.0,107.9;Xk=780;Yk=interp1(X,Y,Xk)执行结果：Yk = 786.8400这里：X和Y分别表示样本点的刻度值和校正值； Xk和Yk分别表示插值点的刻度值和校正值。43什么是样条:是指在飞机或轮船等的设计制造过程中为描绘出光滑的外形曲线(放样)所用的工具三次样条本质上是一段一段的三次多项式拼合而成的曲线，在拼接处,不仅函数是连续的,且一阶和二阶导数也是连续的1946年,Schoenberg将样条引入数学,即所谓的样条函数1）三次样条插值函数4.1.3 三次样条

15、插值法1 方法概述4445这称为三次样条函数的M表达式。 由三次样条函数的 表达式可见，只需确定利用 导数的连续性来推出 满足的方程。推导过程（略）就可确定。46三弯矩方程组可由：求得。472 三次样条插值法的Matlab实现 Matlab中设有专用函数实现三次样条插值法。可调用interp1和spline两个函数实现。interp1函数调用格式为：yi = interp1(x,y,xi,spline)spline函数调用格式为：yi =spline(x,y,xi)两种函数的插值效果等价，其中参数x、y、xi，yi的意义及要求与线性插值interp1中的完全一致。48例4.3：为了绘制乙醇水二

16、元溶液的气液平衡曲线，查到24点实验数据如表4-4所示。表中x和y分别表示液相和气相中乙醇摩尔分数。但表中缺少低浓度点的数据，且数据分布不均匀。试用三次样条插值法对以下插值点插值并画出乙醇水二元溶液的气液平衡曲线。 表4-4 乙醇水二元溶液的气液平衡曲线实验值x00.0180.0540.0780.1110.1470.1720.2030.2100.2410.2970.320y00.1630.3280.3860.4520.4910.5120.5240.5200.5460.5730.577x0.3810.3910.4230.4250.4970.5070.5670.6850.7670.8380.882

17、0.894y0.6000.6080.6210.6230.6530.6560.6840.7430.7960.8460.8840.89449解：用interp1或spline函数，分别求出插值点处的对应函数值并画出乙醇水二元溶液的气液平衡曲线。Matlab计算程序如下：x=0,0.018,0.054,0.078, 0.111,0.147,0.172,0.203,0.210,0.241,0.297,0.320,0.381,0.391,0.423,0.425,0.497,0.507,0.567,0.685,0.767,0.838,0.882,0.894;y=0,0.163,0.328,0.386,0.

18、452,0.491,0.512,0.524,0.520,0.546,0.573,0.577,0.600,0.608,0.621,0.623,0.653,0.656,0.684,0.743,0.796,0.846,0.884,0.894;xk=0.010,0.018,0.020,0.050,0.078,0.100,0.118,0.200,0.300,0.400,0.425,0.500,0.600,0.700,0.800,0.838,0.866,0.899;yk=spline(x,y,xk);plot(x,y,*,xk,yk,r);grid;legend(样本点,三次样条插值)这里：x和y分别表示

19、样本点的液相和气相中乙醇摩尔分数， xk和yk分别表示插值点的液相和气相中乙醇摩尔分数。注：用yk=interp1(x,y,xk,spline); 求插值可得到同样的结果。50图4-1 三次样条插值乙醇水二元溶液的气液平衡曲线514.1.4 分段三次Hermite插值法 已知插值点的函数值和一阶导数值，利用它们可以构造一个三次函数，使其在插值点的函数值和一阶导数值与已知值相等，这样的分段三次函数成为分段三次Hermite插值Matlab函数 pchip P99自学yi=pchip(x,y,xi)或者yi=interp1(x,y,xi,pchip)插值法小结 Lagrange : 给出 y0 y

20、n，选基函数 li(x)，其次数为 节点数 1。 Newton Ln(x)，只是形式不同；节点等距或渐增节点时方便处理。 Hermite： 给出 yi 及 yi ，选 i(x) 及 i(x) 。 Spline：分段低次, 自身光滑, f 的导数只在边界给出。 Matlab 实现：实现分段线性插值不需要编制函数程序，它自身提供了内部的功能函数interp1(一维插值) intep2(二维) interp3(三维) intern(n维) 4.3 MATLAB实现插值用MATLAB作插值计算一维插值函数：yi=interp1(x，y，xi，method)插值方法被插值点插值节点xi处的插值结果nea

21、rest 最邻近插值；linear 线性插值；spline 三次样条插值；cubic 立方插值； 缺省时 分段线性插值 注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围例：从1点12点的11小时内，每隔1小时测量一次温度，测得的温度的数值依次为：5，8，9，15，25，29，31，30，22，25，27，24试估计每隔1/10小时的温度值To MATLAB (temp)hours=1:12;temps=5 8 9 15 25 29 31 30 22 25 27 24;h=1:0.1:12;t=interp1(hours,temps,h,spline); plot(hours,t

22、emps,+,h,t,hours,temps,r:) %作图xlabel(Hour),ylabel(Degrees Celsius)xy机翼下轮廓线例 已知飞机下轮廓线上数据如下，求x每改变0.1时的y值To MATLAB(plane)返回 要求x0,y0单调；x，y可取为矩阵，或x取行向量，y取为列向量，x,y的值分别不能超出x0,y0的范围z=interp2(x0,y0,z0,x,y,method)被插值点插值方法用MATLAB作网格节点数据的插值插值节点被插值点的函数值nearest 最邻近插值；linear 双线性插值；cubic 双三次插值；缺省时 双线性插值.例：测得平板表面35网

23、格点处的温度分别为： 82 81 80 82 84 79 63 61 65 81 84 84 82 85 86 试作出平板表面的温度分布曲面z=f(x,y)的图形输入以下命令：x=1:5;y=1:3;temps=82 81 80 82 84;79 63 61 65 81;84 84 82 85 86;mesh(x,y,temps)1.先在三维坐标画出原始数据，画出粗糙的温度分布曲线图.2以平滑数据,在 x、y方向上每隔0.2个单位的地方进行插值.再输入以下命令:xi=1:0.2:5;yi=1:0.2:3;zi=interp2(x,y,temps,xi,yi,cubic);mesh(xi,yi,

24、zi)画出插值后的温度分布曲面图. To MATLAB (wendu) 通过此例对最近邻点插值、双线性插值方法和双三次插值方法的插值效果进行比较To MATLAB (moutain)返回 插值函数griddata格式为: cz =griddata（x，y，z，cx，cy，method）用MATLAB作散点数据的插值计算 要求cx取行向量，cy取为列向量被插值点插值方法插值节点被插值点的函数值nearest最邻近插值linear 双线性插值cubic 双三次插值v4- MATLAB提供的插值方法缺省时, 双线性插值 例 在某海域测得一些点(x,y)处的水深z由下表给出，船的吃水深度为5英尺，在矩

25、形区域（75，200）（-50，150）里的哪些地方船要避免进入To MATLAB hd1返回4.作出水深小于5的海域范围,即z=5的等高线.2.在矩形区域(75,200)(-50,150)进行插值。 1.输入插值基点数据 3. 作海底曲面图 %程序一：插值并作海底曲面图 x =129.0 140.0 103.5 88.0 185.5 195.0 105.5 157.5 107.5 77.0 81.0 162.0 162.0 117.5 ;y = 7.5 141.5 23.0 147.0 22.5 137.5 85.5 -6.5 -81 3.0 56.5 -66.5 84.0 -33.5 ;z

26、 = 4 8 6 8 6 8 8 9 9 8 8 9 4 9 ;x1=75:1:200;y1=-50:1:150;x1,y1=meshgrid(x1,y1);z1=griddata(x,y,z,x1,y1,v4);meshc(x1,y1,z1) 海底曲面图%程序二：插值并作出水深小于5的海域范围。x1=75:1:200;y1=-50:1:150;x1,y1=meshgrid(x1,y1);z1=griddata(x,y,z,x1,y1,v4); %插值z1(z1=5)=nan; %将水深大于5的置为nan，这样绘图就不会显示出来meshc(x1,y1,z1) 水深小于5的海域范围实验作业1 山

27、区地貌：在某山区测得一些地点的高程如下表：(平面区域1200 x 4000,1200y 3600)，试作出该山区的地貌图和等高线图，并对几种插值方法进行比较返回69 2 最小二乘法曲线拟合 化工实验和工程实践中，可测得许多离散的实验数据和工业数据，通常需要寻找一条连续光滑曲线 来近似反映已知数据组间存在的某种关系的一般趋势，所得近似函数 可以很好地逼近离散数据（ ），这个函数逼近的过程称为曲线拟合或经验建模。0 导入常用最小二乘法曲线拟合。70710.1 最小二乘法曲线拟合的原理 如果观测数据存在较大误差，通常采用“近似函数在各实验点的计算结果与实验结果的偏差平方和最小”的原则建立近似函数。最

28、小若称此曲线拟合法为最小二乘法曲线拟合定义:72经验建模经验建模又分为两种情况：一是无任何理论依据，但有经验公式可供选择，例如很多物性数据（热容、密度、饱和蒸气压）与温度的关系常表示为:二是没有任何经验可循的情况，只能将实验数据画出图形与已知函数图形进行比较，选择图形接近的函数形式作拟合模型。最小二乘法的优点是函数形式多种多样，根据其来源不同，可分为半经验建模和经验建模两种。半经验建模如果建模过程中先由一定的理论依据写出模型结构，再由实验数据估计模型参数，这时建立的模型为半经验模型。例如，描述反应速率常数与温度的关系可用阿仑纽斯方程，即 这种情况下，工作要点在于如何确定函数中的各未知系数 ，

29、73 不论何种建模情况，在选定关联函数的形式之后，就是如何根据实验数据去确定所选关联函数中的待定系数。 最小二乘法按计算方法特点又分为线性最小二乘法和非线性最小二乘法。74对于一元线性函数：测定了m个自变量值：和m个应变量值：计算出m个应变量值：定义误差：2.1 线性最小二乘法 线性最小二乘法是常用的曲线拟合方法。线性最小二乘法又分为一元和多元等不同情况。1 一元线性最小二乘法的方法概述75欲使Q最小，按极值的必要条件，要满足:由最小二乘法：设76可推导出解此方程组，就可求出参数a,b式（4-7）称为一元线性最小二乘法的法方程-(4-7)p = polyfit(x,y,n)772 多元线性最小

30、二乘法的方法概述设系统共有n个影响因子，得到m次实验数据。若可用多元线性函数拟合时，形式如下： （4-8）若k代表第k次实验的数据，则相应的预测值表示为： （4-9）由最小二乘法设： 欲使Q最小，按极值的必要条件，要满足: （4-10）78共n个影响因子，有m次实验数据，若k代表第k次实验的数据，则：设多元线性函数：2 多元线性最小二乘法的方法概述79 则（4-6）转化为以 为未知数的方程组： 式（4-7）称为多元线性最小二乘法的法方程。解此方程组，可求出参数 ，因此拟合方程 便可确定。因此，要求：根据最小二乘原则：要使达到最小令： -(4-6)-(4-7)80 在应用最小二乘法曲线拟合时，通

31、常遇到更多的是非线性函数。对比线性模型拟合，非线性模型拟合要困难的多。4.2.2 非线性最小二乘法最好设法使模型转化为线性形式。有些非线性模型是不能变换成线性模型的，这时应该用直接非线性最小二乘法进行处理 非线性模型拟合的二个途径直接采用非线性拟合通过代换转化为线性关系811 非线性函数的的线性化处理 化工中常见的函数双曲线幂函数指数函数负指数函数对数函数S型曲线n次多项式821） 双曲线令：代换方程为：832） 幂函数两边取对数令：代换方程为：yxb1b=1b11a 843） 指数函数两边取对数令：代换方程为：yxb0b0a85xyb0b0a4）负指数函数 令则代换方程为：86yxb0b00

32、5）对数函数 令则代换方程为：876） S型曲线变形后，令：代换方程为：yx1/a887） n次多项式令：代换方程为：892 非线性直接拟合有些非线性方程无法通过代换法转换成线性方程，则需要采用直接非线性最小二乘法 如：下式是一种常用的饱和蒸汽压计算公式这里介绍非线性直接拟合的常用方法之一-高斯牛顿法90 理论基础：泰勒展开 对于非线性函数若 的近似值为 ，误差为 ，则当初值给定时对非线性函数在初值 附近作泰勒展开，并略去 的二次以上的高次项，可以得到：-(4-12)91其中：92由最小二乘法设 欲使Q最小，按极值的必要条件，要满足: （4-13）则（4-13）转化为以为未知数的方程组： 令9

33、3 将解此法方程所得到的第一套修正值 带入式4-12可求得 ，再用上述方法求得第二套修正值 ，并求得 ，这样经过n此迭代后，若 小到一定程度，就逼近了真值。 即可确定。944.2.3 最小二乘法曲线拟合的Matlab实现 在Matlab中有专用的拟合指令，可以方便地给出求解运算。1 一元多项式拟合函数 polyfit是多项式函数拟合命令，它只能拟合一元多项式。可化为线性的非线性函数拟合，必须变换变量后,然后才可以用此函数。调用格式：p = polyfit(x,y,n)其中: x，y为输入的拟合数据，通常用数组方式输入； n表示多项式的最高阶数，当n=1时不可省略； 输出参数p为拟合多项式的系数

34、注：在Matlab中，多项式是用它的系数构成的行向量表示的。该向量的分量是自左至右，依次表示多项式高次幂项到低次幂项的系数，缺少的幂次项，其系数必须用零填补。如 ，可用表示成 。952 多元线性拟合函数Matlab统计工具箱中使用函数regress实现多元线性拟合 Y = X*B常用调用格式为：b=regress(y，x)其中：x,y代表自变量数据矩阵x和因变量数据向量y； b为多元线性拟合的参数结果向量。注：因变量数据向量y和自变量数据矩阵x按以下排列方式输入 这里：n为变量个数， m为实验次数。963 直接非线性拟合函数Matlab求解非线性拟合的函数较多，多采用最优化方法解决其中nlin

35、fit函数的简单调用格式为：beta = nlinfit(x,y,fun,beta0) 这里：输入数据x，y分别为mn矩阵和m维列向量（n为变量个数， m为实验次数）因变量数据向量y和自变量数据矩阵x排列方式同regress函数fun是事先用 m-文件定义的非线性函数；beta0是回归系数的初值；beta是估计出的回归系数。 974 用最小二乘法拟合生成样条曲线 与样条插值不同，样条拟合并不要求曲线通过全部的数据点。在工程实践与科学研究中，由实验观测得到的一组离散数据一般包含实验误差，因此，样条拟合比样条插值应用更广泛。 Matlab样条工具箱提供的函数较多，常用的三个拟合函数如表4-5所示。

36、 表4.5 常用的三个拟合函数函数名曲线类型拟合准则是否平滑处理csaps()三次样条曲线最小二乘法是spap2()B样条曲线最小二乘法否spaps()B样条曲线最小二乘法是981) 函数csaps( )的用法功能：平滑生成三次样条函数，即对于数据(xi,yi)，所求的三次样条函数y=f(x)满足 调用格式：sp=csaps(x, y, p) ys=csaps(x, y, p, xx, w)输入参数：x,y - 要处理的离散数据(xi, yi)； p - 平滑参数，取值区间为0,1。当p=0时，相当于最小二乘直线拟合；当p=1时，相当于“自然的”三次样条插值，即相当于csapi或spline；

37、 xx - 用于指定在给定点xx上计算其三次样条函数值ys； w - 权值（权重），默认为1 ；输出参数：sp - 拟合得到的样条函数； ys - 在给定点xx上的三次样条函数值。992) 函数spap2( )的用法功能：用最小二乘法拟合生成B样条曲线，即：对于离散数据(xi,yi)，所求的k次样条函数y=f(x)满足 调用格式：sp=spap2(knots,k,x,y) sp=spap2(knots,k,x,y,w)输入参数：knots - 节点序数(knot sequence) k - 样条函数的阶次，一般取k=3，有时取k=4 x,y - 要处理的离散数据(xi,yi) w - 权值（权

38、重），默认为1输出参数：sp - 拟合得到的样条函数1003)函数spaps( )用法功能： 平滑生成B样条函数，所用最小二乘拟合准则同函数spap2。调用格式：sp=spaps(x,y,tol) sp,ys=spaps(x,y,tol,m,w)输入参数：x,y- 要处理的离散数据(xi,yi) tol - 光滑时的允许精度 m - 默认值是2，即平滑生成三次B样条曲线 w - 权值（权重），默认为1输出参数：sp - 拟合得到的样条曲线 ys - 在x上经平滑处理的B样条函数值4) 函数fnval( )的用法功能：计算函数f在给定点x处的函数值values。调用格式：values=fnval

39、(f,x)101表4-7 例4.5数据时间t510152025303540455055浓度y1.272.162.863.443.874.154.374.514.584.624.64解：（1）：建立数学模型，首先对实验点绘图，得：例4.5：某化学反应其反应产物的浓度随时间变化的数据如下：试用最小二乘法关联y与t的关系。102由曲线形状特点，可用 拟合实验数据。103（2）线性化变换令： 则变换为：（3）计算参数，程序如下： t=5,10,15,20,25,30,35,40,45,50,55;y=1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.62,4

40、.64;Y=log(y);X=1./t;p=polyfit(X,Y,1)执行结果： p = -7.4962 1.6516即：b=-7.4962, =1.6516 则：a=exp(1.6516)= 5.2153所以y与t的关系可表示为104例4.6：设某气体反应可表示为：其反应动力学方程可用下列非线性方程表示：式中，V为反应速度，K为反应速度常数，依次为A、B、C的分压，表4-8为实验测定的不同分压下的V值 。试确定此气体反应动力学方程。 VV8.998 8.2982.6998.587.0013.9009.8952.188.1997.0014.4026.053.3103.4019.7962.11

41、7.9016.2035.9004.736.5012.60110.9031.887.0014.3028.1993.357.9972.19917.7971.04表4-8 不同分压下的反应速度105解：（1）首先将动力学方程的两边取对数： 令： 则上式变换为：（2）编程计算参数，程序如下：PA=8.998, 8.199,7.901, 7.001, 7.001, 3.310, 6.501, 7.997;PB=8.298,7.001, 6.203, 4.302, 3.900, 3.401, 2.601, 2.199;PC=2.699, 4.402,5.900 ,8.199, 9.895 ,9.796,1

42、0.903,17.797;V=8.58, 6.05, 4.73, 3.35, 2.18, 2.11, 1.88 ,1.04; x1=log(PA);x2=log(PB);x3=log(PC);y=logV; x=ones(8,1), x1, x2, x3; 将x表示成要求的x形式 b=regress(y,x)106执行结果: b = 1.4914 0.0097 0.6454 -0.6611则：1.4914（K =exp(1.4914)4.4433）；所以此气体反应动力学方程可表示为：0.0097；0.6454；-0.6611 107例4.7：某催化剂活性Y与工作持续时间t的关系为： 将表4-9

43、所列的实验数据通过曲线拟合求系数A、B、C。t(h)02740527089106Y(%)10082.276.371.866.463.361.3表4-9 某催化剂活性Y与工作持续时间t的实验数据108解：本题可以用两种方法解决：（1）先用非线性函数的线性变换法将模型转化为线性形式，然后编程用polyfit 函数得计算结果；（2）直接用非线性拟合函数nlinfit求解。这里用第二种方法，其中系数A、B、C的初始值分别估计为90,0.001,0.001。109t=0 27 24 52 70 89 106;y=100 82.2 76.3 71.8 66.4 63.3 61.3;beta0=90,0.0

44、01,0.001;beta=nlinfit(t,y,model,beta0)其中：model函数的m-文件为：function y=model(beta,t) y=beta(1)*exp(beta(2)*t+beta(3)*t.2)执行结果：beta = 98.5783 -0.0088 0.0000即：A、B、C的结果分别为：98.5783，-0.0088，0.0000。其程序如下：110习题4 3,5(P110)5 3(P126)5.1 引言 对于情况较复杂的实际问题（因素不易化简，作用机理不详）可直接使用数据组建模，寻找简单的因果变量之间的数量关系， 从而对未知的情形作预报。这样组建的模型

45、为拟合模型。 拟合模型的组建主要是处理好观测数据的误差，使用数学表达式从数量上近似因果变量之间的关系。拟合模型的组建是通过对有关变量的观测数据的观察、分析和选择恰当的数学表达方式得到的。 五、拟合的使用及求解5.2 拟合模型的分类 5.2.1 直线拟合5.2.2 曲线拟合5.2.3 观察数据修匀 对于已给一批实测数据，由于实测方法、实验环境等一些外界因素的影响，不可避免地会产生随机干扰和误差。我们自然希望根据数据分布的总趋势去剔除观察数据中的偶然误差，这就是所谓的数据修匀（或称数据平滑）问题。 直 线 拟 合 问 题 引 例 1温度t(C) 20.5 32.7 51.0 73.0 95.7电阻

46、R() 765 826 873 942 1032已知热敏电阻数据：求60C时的电阻R 设 R=at+ba,b为待定系数曲 线 拟 合 问 题 引 例 2 t (h) 0.25 0.5 1 1.5 2 3 4 6 8c (g/ml) 19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01已知一室模型快速静脉注射下的血药浓度数据(t=0注射300mg)求血药浓度随时间的变化规律c(t).在直角坐标系下作图如下(plot)MATLAB(aa1)曲 线 拟 合 问 题 的 提 法已知一组（二维）数据，即平面上 n个点（xi,yi) i=1,n, 寻求一个函数

47、（曲线）y=f(x), 使 f(x) 在某种准则下与所有数据点最为接近，即曲线拟合得最好 +xyy=f(x)(xi,yi)ii 为点（xi,yi) 与曲线 y=f(x) 的距离曲线拟合问题最常用的解法线性最小二乘法的基本思路第一步:先选定一组函数 r1(x), r2(x), ,rm(x), mn, 令 f(x)=a1r1(x)+a2r2(x)+ +amrm(x) （1）其中 a1,a2, ,am 为待定系数 第二步: 确定a1,a2, ,am 的准则（最小二乘准则）：使n个点（xi,yi) 与曲线 y=f(x) 的距离i 的平方和最小 记 问题归结为，求 a1,a2, ,am 使 J (a1,a2, ,am) 最小用MATLAB作线性最小二乘拟合1. 作多项式f(x)=a1xm+ +amx+am+1拟合,可利用已有程序:a=polyfit(x,y,m)输入同长度的数组x，y拟合多项式次数2.多项式在x处的值y可用以下命令计算： y=polyval（