插值与拟合方法.doc

插值与拟合方法 在实际中，常常要处理由实验或测量所得到的一批离散数据插值与拟合方法就是要通过这些数据去确定某一类已知函数的参数或寻找某个近似函数，使所得到的近似函数与已知数据有较高的拟合精度插值问题:要求这个近似函数（曲线或曲面）经过所已知的所有数据点通常插值方法一般用于数据较少的情况数据拟合：不要求近似函数通过所有数据点，而是要求它能较好地反映数据的整体变化趋势。共同点：插值与拟合都是根据实际中一组已知数据来构造一个能够反映数据变化规律的近似函数的方法，由于对近似要求的准则不同，因此二者在数学方法上有很大的差异插值问题的一般提法：已知某函数（未知）的一组观测（或试验）数据，要寻求一个函数，使，则实际中,常常在不知道函数的具体表达式的情况下，对于有实验测量值，寻求另一函数使满足：称此问题为插值问题，并称函数为的插值函数，称为插值节点，称为插值条件，即，则（1） 拉格朗日（Lagrange）插值设函数在个相异点上的函数值为，要求一个次数不超过的代数多项式使在节点上有成立，称之为次代数插值问题，称为插值多项式可以证明次代数插值是唯一的事实上: 可以得到 当时，有二点一次（线性）插值多项式：当=2时，有三点二次（抛物线）插值多项式： （2）牛顿（Newton） 插值牛顿插值的基本思想：由于关于二节点的线性插值为 假设满足插值条件的二次插值多项式一般形式为 由插值条件可得 可以解出 所以 类似的方法，可以得到三次插值多项式等，按这种思想可以得到一般的牛顿插值公式函数的差商及其性质对于给定的函数，用表示关于节点的阶差商，则有一阶差商：，二阶差商：n阶差商：差商有下列性质：（1）差商的分加性：（2）差商的对称性：在中任意调换的次序其值不变牛顿插值公式：一次插值公式为二次插值公式为于是有一般的牛顿插值公式为可以证明：其余项为 实际上，牛顿插值公式是拉格朗日插值公式的一种变形，二者是等价的另外还有著名的埃尔米特（Hermite）插值等 （3）样条函数插值方法 样条，实质上就是由分段多项式光滑连接而成的函数，一般称为多项式样条由于样条函数的特殊性质，决定了样条函数在实际中有着重要的应用 样条函数的一般概念定义 设给定区间的一个分划，如果函数满足条件： (1) 在每个子区间上是次多项式； (2) 及直到-1阶的导数在上连续则称是关于分划的一个次多项式样条函数，称为样条节点，称为内节点，称为边界节点，这类样条函数的全体记作，称为次样条函数空间 若，则是关于分划的次多项式样条函数次多项式样条函数的一般形式为其中和均为任意常数，而在实际中最常用的是和的情况，即为二次样条函数和三次样条函数 二次样条函数：对于上的分划,则其中三次样条函数：对于上的分划,则其中 1 二次样条函数插值 中含有个待定常数，故应需要个插值条件，因此，二次样条插值问题可分为两类： 问题(1)：已知插值节点和相应的函数值，以及端点 (或)处的导数值 (或)，求使得 （5.1） 问题(2)：已知插值节点和相应的导数值，以及端点 (或)处的函数值 (或)，求使得 (5.2)事实上，可以证明这两类插值问题都是唯一可解的对于问题(1)，由条件(5.1)引入记号为未知向量，为已知向量，于是，问题转化为求方程组的解的问题，即可得到二次样条函数的的表达式对于问题(2)的情况类似三次样条函数插值 由于中含有个待定系数，故应需要个插值条件，因此可将三次样条插值问题分为三类： 问题(1)：已知插值节点和相应的函数值，以及两个端点，处的导数值，求使满足条件(5.3) 问题(2)：已知插值节点和相应的函数值，以及两个端点，处的二阶导数值，求使满足条件 (5.4) 问题(3)：类似地，求使满足条件 (5.5)这三类插值问题的条件都是个，可以证明其解都是唯一的8一般的求解方法可以仿照二次样条的情况处理方法，在这里给出一种更简单的方法仅依问题(1)为例，问题(2)和问题(3)的情况类似处理由于在区间上是一个分段光滑，且具有二阶连续导数的三次多项式，则在子区间上是线性函数，记为待定常数由拉格朗日插值公式可得显然在上为常数于是在上有(5.6)则当时，由(5.6)式和问题(1)的条件得故可解得（5.7）将(5.7)式代入(5.6)式得 (5.8) 在上同样的有 (5.9)根据的一阶导数连续性，由(5.9)式得结合(5.7)式整理得引入记号，则（5.10）再由边界条件：得(5.11)联立(5.10),(5.11)式得方程组（5.12）其中，由方程组(6.12)可以唯一解出，代入（5.8）式就可以得三次样条函数的表达式 样条函数插值方法 磨光函数 实际中的许多问题，往往是既要求近似函数（曲线或曲面）有足够的光滑性，又要求与实际函数有相同的凹凸性，一般插值函数和样条函数都不具有这种性质如果对于一个特殊函数进行磨光处理生成磨光函数（多项式），则用磨光函数构造出样条函数作为插值函数，既有足够的光滑性，而且也具有较好的保凹凸性，因此磨光函数在一维插值（曲线）和二维插值（曲面）问题中有着广泛的应用由积分理论可知，对于可积函数通过积分会提高函数的光滑度，因此，我们可以利用积分方法对函数进行磨光处理定义 若为可积函数，对于，则称积分为的一次磨光函数，称为磨光宽度同样的，可以定义的次磨光函数为事实上，磨光函数比的光滑程度要高，且当磨光宽度很小时很接近于 等距样条函数对于任意的函数，定义其步长为的中心差分算子如下：在此取，则是一个单位方波函数（如图5-1），记并取，对进行一次磨光得显然是连续的（如图5-2） o 1 -1/2 0 1/2 x -1 0 1 x 图5-1 图5-2类似地可得到次磨光函数为实际上,可以证明：是分段次多项式，且具有阶连续导数，其阶导数有个间断点，记为从而可知是对应于分划的次多项式样条函数，称之为基本样条函数，简称为次样条由于样条节点为是等距的，故又称为次等距样条函数 对于任意函数的次磨光函数，由归纳法可以得到 4,8 ：特别地，当时，有，从而，且当1时有递推关系 一维等距样条函数插值 等距样条函数与通常的样条如下的关系：定理设有区间的均匀分划，则对任意次样条函数都可以表示为样条函数族的线性组合14根据定理5.1,如果已知曲线上一组点，其中，则可以构造出一条样条磨光曲线（即为样条函数族的线性组合）其中为待定常数用它来逼近曲线，既有较好的精度，又有良好的保凸性实际中，最常用的是的情况，即一般形式为其中个待定系数可以由三类插值条件确定由插值条件(5.3)得(5.13)注意到的局部非零性及其函数值：，当时；且由知，当时则(5.13)中的每一个方程中只有三个非零系数，具体的为(5.14)由方程组(5.14)容易求解出，即可得到三次样条函数表达式类似地，由插值条件(5.4)得待定系数的所满足的方程组为 （5.15）由插值条件(5.5)得待定系数的所满足的方程组为 （5.16） 方程组(5.15),(5.16)也都是容易求解的注：上述等距样条插值公式也适用于近似等距的情形，但在端点和处误差可能较大，实际应用时，为了提高在端点和处的精度，可以适当向左右延拓几个节点二维等距样条函数插值设有空间曲面(未知),如果已知二维等距节点上的值为，则相应的样条磨光曲面的一般形式为其中为待定常数，可以取不同值，常用的也是和的情形这是一种具有良好保凸性的光滑曲面(函数)，在工程设计中是常用的，但只能使用于均匀分划或近似均匀分划的情况（4） 最小二乘拟合方法 最小二乘拟合方法的思想：由于一般插值问题并不总是可解的（即当插值条件多于待定系数的个数时，其问题无解），同时，问题的插值条件本身一般是近似的，为此，只要求在节点上近似地满足插值条件，并使它们的整体误差最小，这就是最小二乘拟合法最小二乘拟合方法可以分为线性最小二乘拟合方法和非线性最小二乘拟合方法 线性最小二乘拟合方法 设是一个线性无关的函数系，则称线性组合为广义多项式如三角多项式：设由给定的一组测量数据和一组正数，求一个广义多项式使得目标函数 (5.17)达到最小，则称函数为数据关于权系数的最小二乘拟合函数，由于关于待定系数是线性的，故此问题又称为线性最小二乘问题注意：这里可根据实际来选择，权系数的选取更是灵活多变的，有时可选取，或，对于，则相应问题称为均方差的极小化问题 最小二乘拟合函数的求解 要使最小二乘问题的目标函数(5.17)达到最小，则由多元函数取得极值的必要条件得即亦即 (5.18)是未知量为的线性方程组，称(5.18)式为正规方程组实际中可适当选择函数系，由正规方程组解出，于是可得最小二乘拟合函数 一般线性最小二乘拟合方法 将上面一元函数的最小二乘拟合问题推广到多元函数，即为多维线性最小二乘拟合问题 假设已知多元函数的一组测量数据和一组线性无关的函数系，求函数 对于一组正数，使得目标函数达到最小其中待定系数由正规方程组确定，此处 注：上面的函数关于都是线性的，这就是线性最小二乘拟合问题，对于这类问题的正规组总是容易求解的如果关于是非线性的，则相应的问题称为非线性最小二乘拟合问题 非线性最小二乘拟合方法 假设已知多元函数的一组测量数据，要求一个关于参数是非线性的函数对一组正数使得目标函数达到最小，则称之为非线性最小二乘问题这类问题属于无约束的最优化问题，一般问题的求解是很复杂的，通常情况下，可以采用共轭梯度法、最速下降法、拟牛顿法和变尺度法等方法求解 实例： 黄河小浪底调水调沙问题 问题的提出2004年月至月黄河进行了第三次调水调沙试验，特别是首次由小浪底、三门峡和万家寨三大水库联合调度，采用接力式防洪预泄放水，形成人造洪峰进行调沙试验获得成功整个试验期为20多天，小浪底从月19日开始预泄放水，直到月13日恢复正常供水结束小浪底水利工程按设计拦沙量为75.5亿立方米，在这之前，小浪底共积泥沙达14.15亿吨这次调水调试验一个重要目的就是由小浪底上游的三门峡和万家寨水库泄洪，在小浪底形成人造洪峰，冲刷小浪底库区沉积的泥沙在小浪底水库开闸泄洪以后，从6月27日开始三门峡水库和万家寨水库陆续开闸放水，人造洪峰于29日先后到达小浪底，7月3日达到最大流量2700立方米/每秒，使小浪底水库的排沙量也不断地增加下面是由小浪底观测站从6月29日到7月10日检测到的试验数据：表5-1: 试验观测数据单位：水流为立方米/每秒，含沙量为公斤/立方米日期6.296.307.17.27.37.4时间8:0020:008:0020:008:0020:008:0020:008:0020:008:0020:00水流量180019002100220023002400250026002650270027202650含沙量326075859098100102108112115116日期7.57.67.77.87.97.10时间8:0020:008:0020:008:0020:008:0020:008:0020:008:0020:00水流量26002500230022002000185018201800175015001000900含沙量11812011810580605030262085注：以上数据主要是根据媒体公开报道的结果整理而成的，不一定与真实数据完全相符现在，根据试验数据建立数学模型研究下面的问题： (1) 给出估算任意时刻的排沙量及总排沙量的方法；(2) 确定排沙量与水流量的变化关系 模型的建立与求解 对于问题（），根据所给问题的试验数据，要计算任意时刻的排沙量，就要确定出排沙量随时间变化的规律，可以通过插值来实现考虑到实际中排沙量应该是随时间连续变化的，为了提高精度，我们采用三次样条函数进行插值 下面构造三次样条函数由试验数据，时间是每天的早点和晚点，间隔都是12个小时，共24个点为了计算方便，令 (5.19)则对应于于是以为插值节点（等距），步长其相应的排沙量为对应关系如下表：表5-2: 插值数据对应关系 单位：排沙量为公斤日期6.296.307.17.27.37.4时间8:0020:008:0020:008:0020:008:0020:008:0020:008:0020:00节点123456789101112排沙量57600 114000157500187000207000235200250000 265200286200302400312800307400日期7.57.67.77.87.97.10时间8:0020:008:0020:008:0020:008:0020:008:0020:008:0020:00节点131415161718192021222324排沙量306800 30000027140023100016000011100091000 54000455003000080004500函数所满足的条件为(1)；(2) .取的三次样条函数一般形式为其中为待定常数，在这里且易知和根据样条函数的性质，在上连续，则有由插值条件(1),(2)可得到下列方程组即 将代入前24个方程中的第一个和最后一个，便可得到方程组，其中， 显然为满秩阵，方程组一定有解，用消元法求解可得问题的解为, , , , , , , , , , ,.将代入 (5.20)即得排沙量的变化规律由（5.19）和(5.20)式可得到第时间段（12小时为一段）内，任意时刻的排沙量则总的排沙量为经计算可得吨，即从6月29日至7月10日小浪底水库排沙总量大约为1.844亿吨,此与媒体报道的排沙量基本相符.对于问题(2),研究排沙量与水量的关系，从试验数据可以看出，开始排沙量是随着水流量的增加而增长，而后是随着水流量的减少而减少显然，变化规律并非是线性的关系，为此，我们问题分为两部分，从开始水流量增加到最大值2720立方米/每秒（即增长的过程）为一段，从水流量的最大值到结束为第二段，分别来研究水流量与排沙量的关系具体数据如表5-3和5-4.表5-3: 第一阶段试验观测数据单位：水流为立方米/每秒，含沙量为公