文档库 最新最全的文档下载
当前位置:文档库 › 曲线拟合与插值理论与实例

曲线拟合与插值理论与实例

曲线拟合与插值理论与实例
曲线拟合与插值理论与实例

第11章曲线拟合与插值

在大量的应用领域中,人们经常面临用一个解析函数描述数据(通常是测量值)的任务。对这个问题有两种方法。在插值法里,数据假定是正确的,要求以某种方法描述数据点之间所发生的情况。这种方法在下一节讨论。这里讨论的方法是曲线拟合或回归。人们设法找出某条光滑曲线,它最佳地拟合数据,但不必要经过任何数据点。图11.1说明了这两种方法。标有'o'的是数据点;连接数据点的实线描绘了线性内插,虚线是数据的最佳拟合。

11.1 曲线拟合

曲线拟合涉及回答两个基本问题:最佳拟合意味着什么?应该用什么样的曲线?可用许多不同的方法定义最佳拟合,并存在无穷数目的曲线。所以,从这里开始,我们走向何方?正如它证实的那样,当最佳拟合被解释为在数据点的最小误差平方和,且所用的曲线限定为多项式时,那么曲线拟合是相当简捷的。数学上,称为多项式的最小二乘曲线拟合。如果这种描述使你混淆,再研究图11.1。虚线和标志的数据点之间的垂直距离是在该点的误差。对各数据点距离求平方,并把平方距离全加起来,就是误差平方和。这条虚线是使误差平方和尽可能小的曲线,即是最佳拟合。最小二乘这个术语仅仅是使误差平方和最小的省略说法。

图11.1 2阶曲线拟合

在MATLAB中,函数polyfit求解最小二乘曲线拟合问题。为了阐述这个函数的用法,让我们以上面图11.1中的数据开始。

? x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1];

? y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];

为了用polyfit,我们必须给函数赋予上面的数据和我们希望最佳拟合数据的多项式的阶次或度。如果我们选择n=1作为阶次,得到最简单的线性近似。通常称为线性回归。相反,如果我们选择n=2作为阶次,得到一个2阶多项式。现在,我们选择一个2阶多项式。

? n=2; % polynomial order

? p=polyfit(x, y, n)

p =

-9.8108 20.1293 -0.0317

polyfit的输出是一个多项式系数的行向量。其解是y = -9.8108x2+20.1293x-0.0317。为了将曲线拟合解与数据点比较,让我们把二者都绘成图。

? xi=linspace(0, 1, 100); % x-axis data for plotting

? z=polyval(p, xi);

为了计算在xi数据点的多项式值,调用MATLAB的函数polyval。

? plot(x, y, ' o ' , x, y, xi, z, ' : ' )

画出了原始数据x和y,用'o'标出该数据点,在数据点之间,再用直线重画原始数据,并用点' : '线,画出多项式数据xi和z。

? xlabel(' x '), ylabel(' y=f(x) '), title(' Second Order Curve Fitting ')

将图作标志。这些步骤的结果表示于前面的图11.1中。

多项式阶次的选择是有点任意的。两点决定一直线或一阶多项式。三点决定一个平方或2阶多项式。按此进行,n+1数据点唯一地确定n阶多项式。于是,在上面的情况下,有11个数据点,我们可选一个高达10阶的多项式。然而,高阶多项式给出很差的数值特性,人们不应选择比所需的阶次高的多项式。此外,随着多项式阶次的提高,近似变得不够光滑,因为较高阶次多项式在变零前,可多次求导。例如,选一个10阶多项式

? pp=polyfit(x, y, 10) ;

? format short e % change display format

? pp.' % display polynomial coefficients as a column

ans =

-4.6436e+005

2.2965e+006

-4.8773e+006

5.8233e+006

-4.2948e+006

2.0211e+006

-6.0322e+005

1.0896e+005

-1.0626e+004

4.3599e+002

-4.4700e-001

要注意在现在情况下,多项式系数的规模与前面的2阶拟合的比较。还要注意在最小(-4.4700e-001)和最大(5.8233e+006)系数之间有7个数量级的幅度差。将这个解作图,并把此图与原始数据及2阶曲线拟合相比较,结果如何呢?

? zz=polyval(pp, xi); % evaluate 10th order polynomial

? plot(x, y, ' o ' , xi, z, ' : ' , xi, zz) % plot data

? xlabel(' x '), ylabel(' y=f(x) '), title(' 2nd and 10th Order curve Fitting ')

在下面的图11.2中,原始数据标以'o',2阶曲线拟合是虚线,10阶拟合是实线。注意,在10阶拟合中,在左边和右边的极值处,数据点之间出现大的纹波。当企图进行高阶曲线拟合时,这种纹波现象经常发生。根据图11.2,显然,‘越多就越好’的观念在这里不适用。

图11.2 2阶和10阶曲线拟合

11.2 一维插值

正如在前一节对曲线拟合所描述的那样,插值定义为对数据点之间函数的估值方法,这些数据点是由某些集合给定。当人们不能很快地求出所需中间点的函数值时,插值是一个有价值的工具。例如,当数据点是某些实验测量的结果或是过长的计算过程时,就有这种情况。

或许最简单插值的例子是MATLAB的作图。按缺省,MATLAB用直线连接所用的数据点以作图。这个线性插值猜测中间值落在数据点之间的直线上。当然,当数据点个数的增加和它们之间距离的减小时,线性插值就更精确。例如,

? x1=linspace(0, 2*pi, 60);

? x2=linspace(0, 2*pi, 6);

? plot(x1, sin(x1), x2, sin(x2), ' - ')

? xlabel(' x '), ylabel(' sin(x) '), title(' Linear Interpolation ')

图11.3 线性插值

图11.3是sine函数的两个图,一个在数据点之间用60个点,它比另一个只用6个点更

光滑和更精确。

如曲线拟合一样,插值要作决策。根据所作的假设,有多种插值。而且,可以在一维以上空间中进行插值。即如果有反映两个变量函数的插值,z=f(x, y),那么就可在x之间和在y之间,找出z的中间值进行插值。MATLAB在一维函数interp1和在二维函数interp2中,提供了许多的插值选择。其中的每个函数将在下面阐述。

为了说明一维插值,考虑下列问题,12小时内,一小时测量一次室外温度。数据存储在两个MATLAB变量中。

? hours=1:12; % index for hour data was recorded

?temps=[5 8 9 15 25 29 31 30 22 25 27 24]; % recorded temperatures

? plot(hours, temps, hours, temps,' + ') % view temperatures

? title(' Temperature ')

? xlabel(' Hour '), ylabel(' Degrees Celsius ')

图11.4 在线性插值下室外温度曲线

正如图11.4看到的,MATLAB画出了数据点线性插值的直线。为了计算在任意给定时间的温度,人们可试着对可视的图作解释。另外一种方法,可用函数interp1。

? t=interp1(hours, temps, 9.3) % estimate temperature at hour=9.3

t =

22.9000

? t=interp1(hours, temps, 4.7) % estimate temperature at hour=4.7

t =

22

? t=interp1(hours, temps, [3.2 6.5 7.1 11.7]) % find temp at many points!

t =

10.2000

30.0000

30.9000

24.9000

interp1的缺省用法是由interp1(x, y, xo)来描述,这里x是独立变量(横坐标),y是应变量(纵坐标),xo是进行插值的一个数值数组。另外,该缺省的使用假定为线性插值。

若不采用直线连接数据点,我们可采用某些更光滑的曲线来拟合数据点。最常用的方法是用一个3阶多项式,即3次多项式,来对相继数据点之间的各段建模,每个3次多项式的头两个导数与该数据点相一致。这种类型的插值被称为3次样条或简称为样条。函数interp1也能执行3次样条插值。

? t=interp1(hours, temps, 9.3, ' spline ') % estimate temperature at hour=9.3

t =

21.8577

? t=interp1(hours, temps, 4.7, ' spline ') % estimate temperature at hour=4.7

t =

22.3143

? t=interp1(hours, temps, [3.2 6.5 7.1 11.7], ' spline ')

t =

9.6734

30.0427

31.1755

25.3820

注意,样条插值得到的结果,与上面所示的线性插值的结果不同。因为插值是一个估计或猜测的过程,其意义在于,应用不同的估计规则导致不同的结果。

一个最常用的样条插值是对数据平滑。也就是,给定一组数据,使用样条插值在更细的间隔求值。例如,

? h=1:0.1:12; % estimate temperature every 1/10 hour

? t=interp1(hours, temps, h, ' spline ') ;

? plot(hours, temps, ' - ' , hours, temps, ' + ' , h, t) % plot comparative results

? title(' Springfield Temperature ')

? xlabel(' Hour '), ylabel(' Degrees Celsius ')

在图11.5中,虚线是线性插值,实线是平滑的样条插值,标有' + '的是原始数据。如要求在时间轴上有更细的分辨率,并使用样条插值,我们有一个更平滑、但不一定更精确地对温度的估计。尤其应注意,在数据点,样条解的斜率不突然改变。作为这个平滑插值的回报,3次样条插值要求更大量的计算,因为必须找到3次多项式以描述给定数据之间的特征。关于样条的更详细信息可见下一章。

图11.5 在不同插值下室外温度曲线

在讨论二维插值之前,了解interp1所强制的二个强约束是很重要的。首先,人们不能要求有独立变量范围以外的结果,例如,interp1(hours, temps, 13.5)导致一个错误,因为hours在1到12之间变化。其次,独立变量必须是单调的。即独立变量在值上必须总是增加的或总是减小的。在我们的例子里,hours是单调的。然而,如果我们已经定义独立变量为一天的实际时间,

? time_of_day=[7:12 1:6] % start at 7AM,end at 6PM

time_of_day =

7 8 9 10 11 12 1 2 3 4 5 6

则独立变量将不是单调的,因为time_of_day增加到12,然后跌到1,再然后增加。如果用time_of_day代替interp1中的hours,将会返回一个错误。同样的理由,人们不能对temps 插值来找出产生某温度的时间(小时),因为temps不是单调的。

11.3 二维插值

二维插值是基于与一维插值同样的基本思想。然而,正如名字所隐含的,二维插值是对两变量的函数z=f(x, y)进行插值。为了说明这个附加的维数,考虑一个问题。设人们对平板上的温度分布估计感兴趣,给定的温度值取自平板表面均匀分布的格栅。

采集了下列的数据:

? width=1:5; % index for width of plate (i.e.,the x-dimension)

? depth=1:3; % index for depth of plate (i,e,,the y-dimension)

?temps=[82 81 80 82 84; 79 63 61 65 81; 84 84 82 85 86] % temperature data

temps =

82 81 80 82 84

79 63 61 65 81

84 84 82 85 86

如同在标引点上测量一样,矩阵temps表示整个平板的温度分布。temps的列与下标depth或y-维相联系,行与下标width或x-维相联系(见图11.6)。为了估计在中间点的温度,我们必须对它们进行辨识。

? wi=1:0.2:5; % estimate across width of plate

? d=2; % at a depth of 2

? zlinear=interp2(width, depth, temps, wi, d) ; % linear interpolation

? zcubic=interp2(width, depth, temps, wi,d, ' cubic ') ; % cubic interpolation

? plot(wi, zlinear, ' - ' , wi, zcubic) % plot results

? xlabel(' Width of Plate '), ylabel(' Degrees Celsius ')

? title( [' Temperature at Depth = ' num2str(d) ] )

另一种方法,我们可以在两个方向插值。先在三维坐标画出原始数据,看一下该数据的粗糙程度(见图11.7)。

? mesh(width, depth, temps) % use mesh plot

? xlabel(' Width of Plate '), ylabel(' Depth of Plate ')

? zlabel(' Degrees Celsius '), axis(' ij '), grid

图11.6 在深度d=2处的平板温度

图11.7 平板温度

然后在两个方向上插值,以平滑数据。

? di=1:0.2:3; % choose higher resolution for depth

? wi=1:0.2:5; % choose higher resolution for width

? zcubic=interp2(width, depth, temps, wi, di, ' cubic ') ; % cubic

? mesh(wi, di, zcubic)

? xlabel(' Width of Plate '), ylabel(' Depth of Plate ')

? zlabel(' Degrees Celsius '), axis(' ij '), grid

上面的例子清楚地证明了,二维插值更为复杂,只是因为有更多的量要保持跟踪。interp2的基本形式是interp2(x, y, z, xi, yi, method)。这里x和y是两个独立变量,z是一个应变量矩阵。x和y对z的关系是

z(i, :) = f(x, y(i)) 和z(:, j) = f(x(j), y).

也就是,当x变化时,z的第i行与y的第i个元素y(i)相关,当y变化时,z的第j列与x的第j个元素x(j)相关,。xi是沿x-轴插值的一个数值数组;yi是沿y-轴插值的一个数值数组。

图11.8 二维插值后的平板温度

可选的参数method可以是'linear','cubic'或'nearest'。在这种情况下,cubic不意味着3次样条,而是使用3次多项式的另一种算法。linear方法是线性插值,仅用作连接图上数据点。nearest方法只选择最接近各估计点的粗略数据点。在所有的情况下,假定独立变量x和y是线性间隔和单调的。关于这些方法的更多的信息,可请求在线帮助,例如,? help interp2,或参阅MATLAB参考手册。

11.4M文件举例

虽然对于许多应用,函数interp1和interp2是很有用的,但它们限制为对单调向量进行插值。在某些情况,这个限制太严格。例如,考虑下面的插值:

? x=linspace(0, 5);

? y=1-exp(-x).*sin(2*pi*x);

? plot(x, y)

图11.9 函数1-exp(-x).*sin(2*pi*x)的曲线

函数interp1可用来在任何值或x的值上估计y值。

? yi=interp1(x, y, 1.8)

yi =

1.1556

然而,interp1不能找出对应于某些y值的x值。例如,如在图11.9上所示,考虑寻找y=1.1处的x值:

图11.10 给y值在函数曲线上求x的值

? plot(x, y, [0, 5], [1.1 1.1] )

从图11.10上,我们看到有四个交点。使用interp1,我们得到:

? xi=interp1(y, x, 1.1)

??? Error using ==> table1

First column of the table must be monotonic.

这个函数interp1失败,由于y不是单调的。

在本章精通MATLAB工具箱所说明的M文件例子,消除了单调性的要求。? table=[x; y].' ; % create column oriented table from data

? xi=mminterp(table, 2, 1.1)

xi =

0.5281 1.1000

0.9580 1.1000

1.5825 1.1000

1.8847 1.1000

这里使用了线性插值,函数mminterp估计了y=1.1处的四个点。由于函数mminterp 的一般性质,要插值的数据是由面向列矩阵给出,在上面的例子中称作为表(table)。第二个输入参量是被搜索矩阵table的列,第三个参量是要找的值。

这个精通MATLAB工具箱函数的主体由下面给出:

function y=mminterp(tab, col, val)

% MMINTERP 1-D Table Search by Linear Interpolation.

% Y=MMINTERP(TAB,COL,VAL) linearly interpolates the table

% TAB searching for the scalar value VAL in the column COL.

% All crossings are found and TAB(:,COL) need not be monotonic.

% Each crossing is returned as a separate row in Y and Y has as

% many columns as TAB.Naturally,the column COL of Y contains

% the value VAL. If VAL is not found in the table,Y=[].

% Copyright (c) 1996 by Prentice-Hall,Inc.

[rt, ct]=size(tab);

if length(val) > 1, error(' VAL must be a scalar. '), end

if col>ct|col < 1, error(' Chosen column outside table width. '), end

if rt < 2, error(' Table too small or not oriented in columns. '), end

above=tab(: , col) > val; % True where > VAL

below=tab(: , col) < val; % True where < VAL

equal=tab(: , col) = = val; % True where = VAL

if all(above = = 0) | all(below = = 0), % handle simplest case

y=tab(find(equal), : ); return

end

pslope=find(below(1:rt-1)&above(2:rt)); % indices where slope is +

nslope=find(below(2:rt)&above(1:rt-1)); % indices where slope is -

ib=sort([pslope; nslope+1]); % put indices below in order

ia=sort([nslope; pslope+1]); % put indices above in order

ie=find(equal); % indices where equal to val

[tmp,ix]=sort( [ib, ie] ); % find where equals fit in result

ieq=ix > length(ib); % True where equals values fit

ry=length(tmp); % # of rows in result y

y=zeros(ry, ct); % poke data into a zero matrix

alpha=(val-tab(ib,col))./(tab(ia,col)-tab(ib,col));

alpha=alpha(: , ones(1, ct)); % duplicate for all columns

y(~ieq, : )=alpha.*tab(ia, : )+(1-alpha).*tab(ib, : ); % interpolated values

y(ieq, : )=tab(ie, : ); % equal values

y( : , col)=val*ones(ry, 1); % remove roundoff error

正如所见的,mminterp利用了find和sort函数、逻辑数组和数组操作技术。没有For 循环和While循环。不论用其中哪一种技术来实现将使运行变慢,尤其对大的表。注意mminterp与含有大于或等于2的任意数列的表一起工作,如同函数interp1一样。而且,在这种情况下,插值变量可以是任意的列。例如,

? z=sin(pi*x); % add more data to table

? table=[x; y; z].' ;

? t=mminterp(table, 2, 1.1) % same interpolation as earlier

t =

0.5281 1.1000 0.9930

0.9580 1.1000 0.1314

1.5825 1.1000 -0.9639

1.8847 1.1000 -0.3533

? t=mminterp(table, 3, -.5) % second third column now

t =

1.1669 0.7316 -0.5000

1.8329 1.1377 -0.5000

3.1671 0.9639 -0.5000

3.8331 1.0187 -0.5000

这些最后的结果估计了x和y在z= -.5处的值。

尽管逐条地对函数mminterp解释如何工作是很有帮助的,但这样做要求有更多的篇幅和时间。解释mminterp如何工作最容易的方法是创建一个小表格,然后,在重要的语句末尾删除分号以后,调用函数。这样,中间值将帮助用户理解函数是如何找到与所需值相符的数据值以及如何执行插值。

前面已阐述了interp1的用法。当用于线性插值时,只要所要求的插值点的个数少,interp1工作很好。在要求许多插值点情况下,由于所用的算法,interp1工作较慢。为了克服这个问题,精通MATLAB工具箱包括了函数mmtable,它的帮助文本是:

?help table

MMTABLE 1-D Monotonic Table Search by Linear Interpolation.

YI=MMTABLE(TAB,COL,VALS) linearly interpolates the table TAB

searching for values VALS in the column COL.

TAB(:,COL) must be monotonic, but need NOT be equally spaced.

YI has as many rows as VALS and as many columns TAB

NaNs are returned where VALS are outside the range of TAB(:,COL).

YI=MMTABLE(TAB,VALS) interpolates using COL=1 and does not return

TAB(:,1) in Y. This matches the usage of TABLE1(TAB,X0).

YI=MMTABLE(X,Y,XI) interpolates the vector X to find YI associated

with XI. This match the usage of INTERP1(X,Y,XI)

This routine is 10X faster than TABLE1which is called by INTERP1.

MMTABLE由线性插值实现一维单调表搜索

YI=MMTABLE(TAB,COL,VALS) 线性地对表TAB进行插值,在列COL中搜索值为VALS

TAB(:,COL)必须是单调的,但不必等价地生成空间。

YI与VALS有同样的行和与TAB有同样的列。

当VALS超出TAB(:,COL)的范围,返回NaNs.

YI=MMTABLE(TAB,VALS) 使用COL=1进行插值,不返回在Y中的TAB(:,1)

这和TABLE1(TAB,X0)的用法匹配。

YI=MMTABLE(X,Y,XI) 为了找出YI和XI的关系,对向量X进行插值。

这和INTERP1(X,Y,XI)的用法匹配。

这个例程比由INTERP1调用TABLE1快10倍。

正如前面描述的,可以用几种方式调用mmtable。此外,要插值的列或向量不需要线性间隔。由于这个原因,mmtable比ilinear函数更普遍。在MATLAB版本5中,interp1将用ilinear来实现线性插值。

11.5 小结

下面的表11.1总结了在MATLAB中所具有的曲线拟合和插值函数。

表11.1

曲线拟合和插值函数

polyfit(x, y, n) 对描述n阶多项式y=f(x)的数据

进行最小二乘曲线拟合

interp1(x, y, xo) 1维线性插值

interp1(x, y, xo, ' spline ') 1维3次样条插值

interp1(x, y, xo, ' cubic ') 1维3次插值

interp2(x, y, Z, xi, yi) 2维线性插值

interp2(x, y, Z, xi, yi, ' cubic ') 2维3次插值

interp2(x, y, Z, xi, yi, ' nearest ') 2维最近邻插值

插值与数据拟合模型

第二讲 插值与数据拟合模型 函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似,由于近似的要求不同,二者的数学方法上是完全不同的。而面对一个实际问题,究竟用插值还是拟合,有时容易确定,有时则并不明显。 在数学建模过程中,常常需要确定一个变量依存于另一个或更多的变量的关系,即函数。但实际上确定函数的形式(线性形式、乘法形式、幂指形式或其它形式)时往往没有先验的依据。只能在收集的实际数据的基础上对若干合乎理论的形式进行试验,从中选择一个最能拟合有关数据,即最有可能反映实际问题的函数形式,这就是数据拟合问题。 一、插值方法简介 插值问题的提法是,已知1+n 个节点n j y x j j ,,2,1,0),,( =,其中j x 互不相同,不妨设b x x x a n =<<<= 10,求任一插值点)(*j x x ≠处的插值*y 。),(j j y x 可以看成是由某个函数)(x g y =产生的,g 的解析表达式可能十分复杂,或不存在封闭形式。也可以未知。 求解的基本思路是,构造一个相对简单的函数)(x f y =,使f 通过全部节点,即),,2,1,0()(n j y x f j j ==,再由)(x f 计算插值,即*)(*x f y =。 1.拉格朗日多项式插值 插值多项式 从理论和计算的角度看,多项式是最简单的函数,设)(x f 是n 次多项式,记作 0111)(a x a x a x a x L n n n n n ++++=-- (1) 对于节点),(j j y x 应有 n j y x L j j n ,,2,1,0,)( == (2) 为了确定插值多项式)(x L n 中的系数011,,,,a a a a n n -,将(1)代入(2),有 ???????=++++=++++=++++---n n n n n n n n n n n n n n n n y a x a x a x a y a x a x a x a y a x a x a x a 01110111110001010 (3) 记 T n T n n n n n n n n n n y y y Y a a a A x x x x x x X ),,,(,),,,(,11110011111 100 ==?????? ? ??=---- 方程组(3)简写成 Y XA = (4) 注意X det 是Vandermonde 行列式,利用行列式性质可得 ∏≤<≤-= n k j j k x x X 0)(det 因j x 互不相同,故0det ≠X ,于是方程(4)中A 有唯一解,即根据1+n 个节点可以确定唯一的n 次插值多项式。 拉格朗日插值多项式 实际上比较方便的做法不是解方程(4)求A ,而是先构造一组基函数: n i x x x x x x x x x x x x x x x x x l n i i i i i i n i i i ,,2,1,0,) ())(()()())(()()(110110 =--------=+-+- (5) )(x l i 是n 次多项式,满足

曲线拟合与插值1一维插值

一维插值 正如在前一节对曲线拟合所描述的那样,插值定义为对数据点之间函数的估值方法,这些数据点是由某些集合给定。当人们不能很快地求出所需中间点的函数值时,插值是一个有价值的工具。例如,当数据点是某些实验测量的结果或是过长的计算过程时,就有这种情况。 或许最简单插值的例子是MATLAB的作图。按缺省,MATLAB用直线连接所用的数据点以作图。这个线性插值猜测中间值落在数据点之间的直线上。当然,当数据点个数的增加和它们之间距离的减小时,线性插值就更精确。例如, ? x1=linspace(0, 2*pi, 60); ? x2=linspace(0, 2*pi, 6); ? plot(x1, sin(x1), x2, sin(x2), ' - ') ? xlabel(' x '),ylabel(' sin(x) '),title(' Linear Interpolation ') 图11.3线性插值 图11.3是sine函数的两个图,一个在数据点之间用60个点,它比另一个只用6个点更光滑和更精确。 如曲线拟合一样,插值要作决策。根据所作的假设,有多种插值。而且,可以在一维以上空间中进行插值。即如果有反映两个变量函数的插值,z=f(x, y),那么就可在x之间和在y之间,找出z的中间值进行插值。MATLAB在一维函数interp1和在二维函数interp2中,提供了许多的插值选择。其中的每个函数将在下面阐述。 为了说明一维插值,考虑下列问题,12小时内,一小时测量一次室外温度。数据存储在两个MATLAB变量中。 ? hou rs=1:12;%index for hour data was recorded ? temps=[5 89152529313022252724]; %recorded temperatures ? plot(hours, temps, hours, temps,' + ')%view temperatures

实验四 插值法与曲线拟合

计算方法实验报告 专业班级:医学信息工程一班姓名:陈小芳学号:201612203501002 实验成绩: 1.【实验题目】 插值法与曲线拟合 2.【实验目的】 3.【实验内容】 4. 【实验要求】

5. 【源程序(带注释)】 (1)拉格朗日插值 #include #include #include #include #include #define n 4 //插值节点的最大下标 main() { double x1[n+1]={0.4,0.55,0.65,0.8,0.9}; double y1[n+1]={0.4175,0.57815,0.69657,0.88811,1.02652}; double Lagrange(double x1[n+1],double y1[n+1],float t); int m,k;float x,y;float X;double z; printf("\n The number of the interpolation points is m ="); //输入插值点的个数 while(!scanf("%d",&m)) { fflush(stdin); printf("\n输入错误,请重新输入:\n"); printf("\n The number of the interpolation points is m ="); } for(k=1;k<=m;k++) { printf("\ninput X%d=",k); while(!scanf("%f",&X)) { fflush(stdin); printf("\n输入错误,请重新输入:\n"); printf("\ninput X%d=",k); } z=Lagrange(x1,y1,X); printf("P(%f)=%f\n",X,z); } getch(); return (0); } double Lagrange(double x[n+1],double y[n+1],float X) { int i,j;

第五章插值与拟合答案—牟善军

习题5.1: Matlab程序如下: clc,clear x=1:0.5:10; y=x.^3-6*x.^2+5*x-3; y0=y+rand; f1=polyfit(x,y0,1) y1=polyval(f1,x); plot(x,y,'+',x,y1); grid on title('一次拟合曲线'); figure(2); f2=polyfit(x,y0,2) y2=polyval(f2,x); plot(x,y,'+',x,y2); grid on title('二次拟合曲线'); figure(3); f4=polyfit(x,y0,4) y3=polyval(f4,x); plot(x,y,'+',x,y3); grid on title('四次拟合曲线'); figure(4); f6=polyfit(x,y0,6) y4=polyval(f6,x); plot(x,y,'+',x,y4); grid on title('六次拟合曲线'); 计算结果及图如下 f1 = 43.2000 -148.8307 f2 = 10.5000 -72.3000 90.0443

f4 = 0.0000 1.0000 -6.0000 5.0000 -2.3557 f6 = -0.0000 0.0000 -0.0000 1.0000 -6.0000 5.0000 -2.3557 5.2高程数据问题解答如下:matlab程序: clc,clear x0=0:400:5600 y0=0:400:4800 z0=[1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 620 460 370 350 1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500

曲线拟合与插值理论与实例

第11章曲线拟合与插值 在大量的应用领域中,人们经常面临用一个解析函数描述数据(通常是测量值)的任务。对这个问题有两种方法。在插值法里,数据假定是正确的,要求以某种方法描述数据点之间所发生的情况。这种方法在下一节讨论。这里讨论的方法是曲线拟合或回归。人们设法找出某条光滑曲线,它最佳地拟合数据,但不必要经过任何数据点。图11.1说明了这两种方法。标有'o'的是数据点;连接数据点的实线描绘了线性内插,虚线是数据的最佳拟合。 11.1 曲线拟合 曲线拟合涉及回答两个基本问题:最佳拟合意味着什么?应该用什么样的曲线?可用许多不同的方法定义最佳拟合,并存在无穷数目的曲线。所以,从这里开始,我们走向何方?正如它证实的那样,当最佳拟合被解释为在数据点的最小误差平方和,且所用的曲线限定为多项式时,那么曲线拟合是相当简捷的。数学上,称为多项式的最小二乘曲线拟合。如果这种描述使你混淆,再研究图11.1。虚线和标志的数据点之间的垂直距离是在该点的误差。对各数据点距离求平方,并把平方距离全加起来,就是误差平方和。这条虚线是使误差平方和尽可能小的曲线,即是最佳拟合。最小二乘这个术语仅仅是使误差平方和最小的省略说法。

图11.1 2阶曲线拟合 在MATLAB中,函数polyfit求解最小二乘曲线拟合问题。为了阐述这个函数的用法,让我们以上面图11.1中的数据开始。 ? x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1]; ? y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2]; 为了用polyfit,我们必须给函数赋予上面的数据和我们希望最佳拟合数据的多项式的阶次或度。如果我们选择n=1作为阶次,得到最简单的线性近似。通常称为线性回归。相反,如果我们选择n=2作为阶次,得到一个2阶多项式。现在,我们选择一个2阶多项式。 ? n=2; % polynomial order ? p=polyfit(x, y, n) p = -9.8108 20.1293 -0.0317 polyfit的输出是一个多项式系数的行向量。其解是y = -9.8108x2+20.1293x-0.0317。为了将曲线拟合解与数据点比较,让我们把二者都绘成图。

曲线拟合与插值3二维插值

二维插值 二维插值是基于与一维插值同样的基本思想。然而,正如名字所隐含的,二维插值是对两变量的函数 z=f(x, y)进行插值。为了说明这个附加的维数,考虑一个问题。设人们对平板上的温度分布估计感兴趣,给定的温度值取自平板表面均匀分布的格栅。 采集了下列的数据: ? width=1:5;%index for width of plate (i.e.,the x-dimension) ? depth=1:3;%index for depth of plate (i,e,,the y-dimension) ? t emps=[8281808284; 7963616581; 8484828586]%temperature data temps = 8281808284 7963616581 8484828586 如同在标引点上测量一样,矩阵temps表示整个平板的温度分布。temps的列与下标depth或y-维相联系,行与下标width或x-维相联系(见图11.6)。为了估计在中间点的温度,我们必须对它们进行辨识。 ? wi=1:0.2:5;%estimate across width of plat e ? d=2;%at a depth of 2 ? zlinear=interp2(width, depth, temps, wi, d) ;%linear interpolation ? zcubic=interp2(width, depth, temps, wi,d, ' cubic ') ;%cubic interpolation ? plot(wi, zlinear, ' - ' , wi, zcubic)%plot results ? xlabel(' Width of Plate '),y label(' Degrees Celsius ') ? title( [' Temperature at Depth ='num2str(d) ] ) 另一种方法,我们可以在两个方向插值。先在三维坐标画出原始数据,看一下该数据的粗糙程度(见图11.7)。? mesh(width, depth, temps)%use mesh plot ? xlabel(' Width of Plate '),ylabel(' Depth of Plate ') ? zlabel(' Degrees Celsius '),axis(' ij '),grid

(完整版)Matlab学习系列13.数据插值与拟合

13. 数据插值与拟合 实际中,通常需要处理实验或测量得到的离散数据(点)。插值与拟合方法就是要通过离散数据去确定一个近似函数(曲线或曲面),使其与已知数据有较高的拟合精度。 1.如果要求近似函数经过所已知的所有数据点,此时称为插值问 题(不需要函数表达式)。 2.如果不要求近似函数经过所有数据点,而是要求它能较好地反 映数据变化规律,称为数据拟合(必须有函数表达式)。 插值与拟合都是根据实际中一组已知数据来构造一个能够反映数据变化规律的近似函数。区别是:【插值】不一定得到近似函数的表达形式,仅通过插值方法找到未知点对应的值。【拟合】要求得到一个具体的近似函数的表达式。 因此,当数据量不够,但已知已有数据可信,需要补充数据,此时用【插值】。当数据基本够用,需要寻找因果变量之间的数量关系(推断出表达式),进而对未知的情形作预测,此时用【拟合】。

一、数据插值 根据选用不同类型的插值函数,逼近的效果就不同,一般有:(1)拉格朗日插值(lagrange插值) (2)分段线性插值 (3)Hermite (4)三次样条插值 Matlab 插值函数实现: (1)interp1( ) 一维插值 (2)intep2( ) 二维插值 (3)interp3( ) 三维插值 (4)intern( ) n维插值 1.一维插值(自变量是1维数据) 语法:yi = interp1(x0, y0, xi, ‘method’) 其中,x0, y0为原离散数据(x0为自变量,y0为因变量);xi为需要插值的节点,method为插值方法。 注:(1)要求x0是单调的,xi不超过x0的范围; (2)插值方法有‘nearest’——最邻近插值;‘linear’——线性插值;‘spline’——三次样条插值;‘cubic’——三次插值;

第四章曲线拟合和多项式插值 - Hujiawei Bujidao

数值计算之美SHU ZHI JI SUAN ZHI MEI 胡家威 http://hujiaweibujidao.github.io/ 清华大学逸夫图书馆·北京

内容简介 本书是我对数值计算中的若干常见的重要算法及其应用的总结,内容还算比较完整。 本人才疏学浅,再加上时间和精力有限,所以本书不会详细介绍很多的概念,需要读者有一定的基础或者有其他的参考书籍,这里推荐参考文献中的几本关于数值计算的教材。 本书只会简单介绍下算法的原理,对于每个算法都会附上我阅读过的较好的参考资料以及算法的实现(Matlab或者其他语言),大部分代码是来源于参考文献[1]或者是经过我改编而成的,肯定都是可以直接使用的,需要注意的是由于Latex对代码的排版问题,导致中文注释中的英文字符经常出现错位,对于这种情况请读者自行分析,不便之处还望谅解。写下这些内容的目的是让自己理解地更加深刻些,顺便能够作为自己的HandBook,如有错误之处,还请您指正,本人邮箱地址是:hujiawei090807@https://www.wendangku.net/doc/0a12841811.html,。

目录 第四章曲线拟合和多项式插值1 4.1曲线拟合 (1) 4.1.1使用线性方程进行曲线拟合 (1) 4.1.2非线性方程进行曲线拟合 (2) 4.1.3使用二次或者高次多项式进行曲线拟合[最小二乘问题].3 4.2多项式插值 (4) 4.2.1拉格朗日插值多项式 (4) 4.2.2牛顿插值多项式 (5) 4.2.3分段线性插值 (7) 4.2.4保形分段三次插值 (8) 4.2.5三次样条插值 (10) 4.3Matlab函数解析 (13) 参考文献14

插值与拟合(使用插值还是拟合)

利用matlab实现插值与拟合实验 张体强1026222 张影 晁亚敏 [摘要]:在测绘学中,无论是图形处理,还是地形图处理等,大多离不开插值与拟合的应用,根据插值与拟合原理,构造出插值和拟合函数,理解其原理,并在matlab平台下,实现一维插值,二维插值运算,实现多项式拟合,非线性拟合等,并在此基础上,联系自己所学专业,分析其生活中特殊例子,提出问题,建立模型,编写程序,以至于深刻理解插值与拟合的作用。 [关键字]: 测绘学插值多项式拟合非线性拟合 [ Abstract]: in surveying and mapping, whether the graphics processing, or topographic map processing and so on, are inseparable from the interpolation and fitting application, according to the interpolation and fitting theory, construct the fitting and interpolation function, understanding its principle, and MATLAB platform, achieve one-dimensional interpolation, two-dimensional interpolation, polynomial fitting, non-linear fitting, and on this basis, to contact their studies, analysis of their living in a special example, put forward the question, modeling, programming, so that a deep understanding of interpolation and fitting function. [ Key words]: Surveying and mapping interpolation polynomial fitting nonlinear

数值分析法 曲线拟合法插值建模法

数值分析法 相关知识 在生产和科学实验中,自变量x 与因变量y 间的函数关系()y f x =有时不能写出解析表达式,而只能得到函数在若干点的函数值或导数值,或者表达式过于复杂而需要较大的计算量。当要求知道其它点的函数值时,需要估计函数值在该点的值。 为了完成这样的任务,需要构造一个比较简单的函数()y x ?=,使函数在观测点的值等于已知的值,或使函数在该点的导数值等于已知的值,寻找这样的函数()y x ?=有很多方法。根据测量数据的类型有以下两类处理观测数据的方法。 (1)测量值是准确的,没有误差,一般用插值。 (2)测量值与真实值有误差,一般用曲线拟合。 曲线拟合法 已知离散点上的数据集1122{(,),(,),,(,)}n n x y x y x y ,即已知在点集12{,,,}n x x x 上的函数值12{,,,}n y y y ,构造一个解析函数(其图形为一曲线)使()f x 在原离散点 i x 上尽可能接近给定的i y 值,这一过程称为曲线拟合。 曲线拟合的一般步骤是先根据实验数据,结合相关定律,将要寻求的最恰当的拟合曲线方程形式预测出来,再用其他的数学方法确定经验公式中的参数。 对于事先给定的一组数据,确定经验公式一般可分为三步进行: (1)、确定经验公式的形式:根据系统和测定的数据的特点,并参照已知图形的特点确定经验公式的形式。 (2)、确定经验公式中的待定系数:计算待定系数的方法有许多常用的法有图示法、均值法、差分法、最小二乘法、插值法等。 (3)、检验:求出经验公式后,还要将测定的数据与用经验公式求出的理论

数据作比较,验证经验公式的正确性,必要时还要修正经验公式。 关于确定经验公式的形式,可从以下几个方面入手: (1)、利用已知的结论确定经验公式形式,如由已知的胡克定律可以确定在一定条件下,弹性体的应变与应力呈线性关系等。 (2)、从分析实验数据的特点入手,将之与已知形式的函数图形相对照,确定经验公式的形式。 (3)、描点作图法:将已知的点用光滑的曲线连接起来,寻找曲线的形式。 (4)、多项式近似、线性插值或样条插值等。多项式近似是工程中十分常见的方法,它首先需要确定多项式的次数,一般可以用差分法、差商法来估计。 <一>、差分方程法 <1>、差分方程:差分方程反映的是关于离散变量的取值与变化规律。通过建立一个或几个离散变量取值所满足的平衡关系,从而建立差分方程。 (1)、说明:差分方程就是针对要解决的目标,引入系统或过程中的离散变量,根据实际背景的规律、性质、平衡关系,建立离散变量所满足的平衡关系等式,从而建立差分方程。通过求出和分析方程的解,或者分析得到方程解的特别性质(平衡性、稳定性、渐近性、振动性、周期性等),从而把握这个离散变量的变化过程的规律,进一步再结合其他分析,得到原问题的解。 (2)、应用:差分方程模型有着广泛的应用。实际上,连续变量可以用离散变量来近似和逼近,从而微分方程模型就可以近似于某个差分方程模型。差分方程模型有着非常广泛的实际背景。在经济金融保险领域、生物种群的数量结构规律分析、疾病和病虫害的控制与防治、遗传规律的研究等许许多多的方面都有着非常重要的作用。可以这样讲,只要牵涉到关于变量的规律、性质,就可以适当地用差分方程模型来表现与分析求解。

曲线拟合和插值运算原理和方法

实验10 曲线拟合和插值运算 一. 实验目的 学会MATLAB 软件中软件拟合与插值运算的方法。 二. 实验内容与要求 在生产和科学实验中,自变量x 与因变量y=f(x)的关系式有时不能直接写出表达式,而只能得到函数在若干个点的函数值或导数值。当要求知道观测点之外的函数值时,需要估计函数值在该点的值。 要根据观测点的值,构造一个比较简单的函数y=t (x),使函数在观测点的值等于已知的数值或导数值,寻找这样的函数t(x),办法是很多的。 根据测量数据的类型有如下两种处理观测数据的方法。 (1) 测量值是准确的,没有误差,一般用插值。 (2) 测量值与真实值有误差,一般用曲线拟合。 MATLAB 中提供了众多的数据处理命令,有插值命令,拟合命令。 1.曲线拟合 已知离散点上的数据集[(1x ,1y ),………(n x ,n y )],求得一解析函数y=f (x),使f(x)在原离散点i x 上尽可能接近给定i y 的值,之一过程叫曲线拟合。最常用的的曲线拟合是最小二乘法曲线拟合,拟合结果可使误差的平方和最小,即使求使21|()|n i i i f x y =-∑ 最小的f(x). 格式:p=polyfit(x,Y ,n). 说明:求出已知数据x,Y 的n 阶拟合多项式f(x)的系数p ,x 必须是单调的。 [例 1.9] >>x=[0.5,1.0,1.5,2.0,2.5,3.0]; %给出数据点的x 值 >>y=[1.75,2.45,3.81,4.80,7.00,8.60]; %给出数据点的y 值 >>p=polyfit (x,y,2); %求出二阶拟合多项式f(x)的系数 >>x1=0.5:0.05:3.0; %给出x 在0.5~3.0之间的离散值 >>y1=polyval(p,1x ); %求出f(x)在1x 的值 >>plot(x,y,?*r ?, 11,x y ?-b ?) %比较拟合曲线效果 计算结果为: p= 0.5614 0.8287 1.1560 即用f(x)=0.56142 x +0.8287x+1.1560拟合已知数据,拟合曲线效果如图所示。

数学建模案例分析-- 插值与拟合方法建模1数据插值方法及应用

第十章 插值与拟合方法建模 在生产实际中,常常要处理由实验或测量所得到的一批离散数据,插值与拟合方法就是要通过这些数据去确定某一类已经函数的参数,或寻求某个近似函数使之与已知数据有较高的拟合精度。插值与拟合的方法很多,这里主要介绍线性插值方法、多项式插值方法和样条插值方法,以及最小二乘拟合方法在实际问题中的应用。相应的理论和算法是数值分析的内容,这里不作详细介绍,请参阅有关的书籍。 §1 数据插值方法及应用 在生产实践和科学研究中,常常有这样的问题:由实验或测量得到变量间的一批离散样点,要求由此建立变量之间的函数关系或得到样点之外的数据。与此有关的一类问题是当原始数据 ),(,),,(),,(1100n n y x y x y x 精度较高,要求确定一个初等函数)(x P y =(一般用多项式或分段 多项式函数)通过已知各数据点(节点),即n i x P y i i ,,1,0,)( ==,或要求得函数在另外一些点(插值点)处的数值,这便是插值问题。 1、分段线性插值 这是最通俗的一种方法,直观上就是将各数据点用折线连接起来。如果 b x x x a n =<<<= 10 那么分段线性插值公式为 n i x x x y x x x x y x x x x x P i i i i i i i i i i ,,2,1,,)(11 1 11 =≤<--+--= ----- 可以证明,当分点足够细时,分段线性插值是收敛的。其缺点是不能形成一条光滑曲线。 例1、已知欧洲一个国家的地图,为了算出它的国土面积,对地图作了如下测量:以由西向东方向为x 轴,由南向北方向为y 轴,选择方便的原点,并将从最西边界点到最东边界点在x 轴上的区间适当的分为若干段,在每个分点的y 方向测出南边界点和北边界点的y 坐标y1和y2,这样就得到下表的数据(单位:mm )。 根据地图的比例,18 mm 相当于40 km 。

第4、5讲 插值与拟合 作业参考答案

第四、五讲作业题参考答案 一、填空题 1、拉格朗日插值基函数在节点上的取值是( 0或1 )。 2、当1,1,2x =-,时()034f x =-, ,,则()f x 的二次插值多项式为 ( 2527 633 x x +- )。 3、由下列数据 所确定的唯一插值多项式的次数为( 2次 )。 4、根据插值的定义,函数()x f x e -=在[0,1]上的近似一次多项式1()P x = ( 1(1)1e x --+ ),误差估计为( 18 )。 5、在做曲线拟合时,对于拟合函数x y ax b = +,引入变量变换y =( 1 y ),x =( 1 x )来线性化数据点后,做线性拟合y a bx =+。 6、在做曲线拟合时,对于拟合函数Ax y Ce =,引入变量变换( ln()Y y = )、 X x =和B C e =来线性化数据点后,做线性拟合Y AX B =+。 7、设3()1f x x x =+-,则差商[0,1,2,3]f =( 1 )。 8、在做曲线拟合时,对于拟合函数()A f x Cx =,可使用变量变换(ln Y y =)(ln X x = )和B C e =来线性化数据点后,做线性拟合Y AX B =+。 9、设(1)1,(0)0,(1)1,(2)5,()f f f f f x -====则的三次牛顿插值多项式为 ( 3211 66x x x +-),其误差估计式为( 4()(1)(1)(2),(1,2)24f x x x x ξξ+--∈-) 10、三次样条插值函数()S x 满足:()S x 在区间[,]a b 内二阶连续可导, (),,0,1,2,,,k k k k S x y x y k n ==(已知)且满足()S x 在每一个子区间1[,] k k x x +上是( 三次多项式 )。

matlab实现插值法和曲线拟合电子教案

m a t l a b实现插值法和 曲线拟合

插值法和曲线拟合 电子科技大学 摘要:理解拉格朗日多项式插值、分段线性插值、牛顿前插,曲线拟合,用matlab编程求解函数,用插值法和分段线性插值求解同一函数,比较插值余项;用牛顿前插公式计算函数,计算函数值;对于曲线拟 合,用不同曲线拟合数据。 关键字:拉格朗日插值多项式;分段线性插值;牛顿前插;曲线拟合 引言: 在数学物理方程中,当给定数据是不同散点时,无法确定函数表达式,求解函数就需要很大的计算量,我们有多种方法对给定的表格函数进行求解,我们这里,利用插值法和曲线拟合对函数进行求解,进一步了解函数性质,两种方法各有利弊,适合我们进行不同的散点函数求解。 正文: 一、插值法和分段线性插值 1拉格朗日多项式原理 对某个多项式函数,已知有给定的k + 1个取值点: 其中对应着自变量的位置,而对应着函数在这个位置的取值。 假设任意两个不同的x j都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为: 其中每个为拉格朗日基本多项式(或称插值基函数),其表达式为: [3] 拉格朗日基本多项式的特点是在上取值为1,在其它的点 上取值为0。 2分段线性插值原理 给定区间[a,b], 将其分割成a=x 0

插值和拟合区别

插值和拟合区别 运输1203黎文皓通过这个学期的《科学计算与数学建模》课程的学习,使我掌握了不少数学模型解决实际问题的方法,其中我对于插值与拟合算法这一章,谈一谈自己的看法可能不是很到位,讲得不好的地方也请老师见谅。 首先,举一个简单的例子说明一下这个问题。 如果有100个平面点,要求一条曲线近似经过这些点,可有两种方法:插值和拟合。 我们可能倾向于用一条(或者分段的多条)2次、3次或者说“低次”的多项式曲线而不是99次的曲线去做插值。也就是说这条插值曲线只经过其中的3个、4个(或者一组稀疏的数据点)点,这就涉及到“滤波”或者其他数学方法,也就是把不需要90多个点筛选掉。如果用拟合,以最小二乘拟合为例,可以求出一条(或者分段的多条)低次的曲线(次数自己规定),逼近这些数据点。具体方法参见《数值分析》中的“线性方程组的解法”中的“超定方程的求解法”。经过上面例子的分析,我们可以大致的得到这样一个结论。插值就是精确经过,拟合就是逼近。 插值和拟合都是函数逼近或者数值逼近的重要组成部分。他们的共同点都是通过已知一些离散点集M上的约束,求取一个定义在连续集合S(M包含于S)的未知连续函数,从而达到获取整体规律目的,即通过"窥几斑"来达到"知全豹"。 所谓拟合是指已知某函数的若干离散函数值{f1,f2,…,fn},通过调

整该函数中若干待定系数f(λ1, λ2,…,λ3), 使得该函数与已知点集的差别(最小二乘意义)最小。如果待定函数是线性,就叫线性拟合或者线性回归(主要在统计中),否则叫作非线性拟合或者非线性回归。表达式也可以是分段函数,这种情况下叫作样条拟合。 而插值是指已知某函数的在若干离散点上的函数值或者导数信息,通过求解该函数中待定形式的插值函数以及待定系数,使得该函数在给定离散点上满足约束。 从几何意义上将,拟合是给定了空间中的一些点,找到一个已知形式未知参数的连续曲面来最大限度地逼近这些点;而插值是找到一个(或几个分片光滑的)连续曲面来穿过这些点。 不过是插值还是拟合都是建立在一定的数学模型的基础上进行的。多项式插值虽然在一定程度上解决了由函数表求函数的近似表达式的问题,但是在逼近曲线上有明显的缺陷,很可能不能很好的表示函数的走向,存在偏差,在实际问题中我们往往通过函数近似表达式的拟合法来得到一个较为准却的表达式。

插值法与数据拟合法

第七讲插值方法与数据拟合 § 7.1 引言 在工程和科学实验中,常常需要从一组实验观测数据(x i , y i ) (i= 1, 2, …, n) 揭示自变量x与因变量y 之间的关系,一般可以用一个近似的函数关系式y = f (x) 来表示。函数f (x) 的产生办法因观测数据与要求的不同而异,通常可采用两种方法:插值与数据拟合。 § 7.1.1 插值方法 1.引例1 已经测得在北纬32.3?海洋不同深度处的温度如下表: 根据这些数据,我们希望能合理地估计出其它深度(如500米、600米、1000米…)处的水温。 解决这个问题,可以通过构造一个与给定数据相适应的函数来解决,这是一个被称为插值的问题。 2.插值问题的基本提法 对于给定的函数表 其中f (x) 在区间[a, b] 上连续,x0,x1,…,x n为[a, b] 上n + 1个互不相同的点,要求在一个性质优良、便于计算的函数类{P(x)} 中,选出一个使 P(x i ) = y i,i= 0, 1, …, n(7.1.1) 成立的函数P(x) 作为 f (x) 的近似,这就是最基本的插值问题(见图7.1.1)。 为便于叙述,通常称区间[a, b] 为插值区间,称点x0,x1,…,x n为插值节点,称函数类{P(x)} 为插值函数类,称式(7.1.1) 为插值条件,称函数P(x) 为插值函数,称f (x) 为被插函数。求插值函数P(x) 的方法称为插值法。 § 7.1.2 数据拟合 1.引例2 在某化学反应中,已知生成物的浓度与时间有关。今测得一组数据如下: 根据这些数据,我们希望寻找一个y = f (t) 的近似表达式(如建立浓度y与时间t之间的经验公式等)。从几何上看,就是希望根据给定的一组点(1, 4.00),…,(16, 10.60),求函数y = f (t) 的图象的一条拟合曲

MATLAB中的曲线拟合与插值

MATLAB中的曲线拟合和插值 在大量的使用领域中,人们经常面临用一个分析函数描述数据(通常是测量值)的任务。对这个问题有两种方法。在插值法里,数据假定是正确的,要求以某种方法描述数据点之间所发生的情况。这种方法在下一节讨论。这里讨论的方法是曲线拟合或回归。人们设法找出某条光滑曲线,它最佳地拟合数据,但不必要经过任何数据点。图11.1说明了这两种 方法。标有’0'的是数据点;连接数据点的实线描绘了线性内插,虚线是数据的最佳拟合。 11.1 曲线拟合 曲线拟合涉及回答两个基本问题:最佳拟合意味着什么?应该用什么样的曲线?可用 许多不同的方法定义最佳拟合,并存在无穷数目的曲线。所以,从这里开始,我们走向何方?正如它证实的那样, 当最佳拟合被解释为在数据点的最小误差平方和,且所用的曲线限定为多项式时,那么曲线拟合是相当简捷的。数 学上,称为多项式的最小二乘曲线拟合。如果这种描述使你混淆,再研究图11.1。虚线和标志的数据点之间的垂直距离是在该点的 误差。对各数据点距离求平方,并把平方距离全加起来,就是误差平方和。这条虚线是使 误差平方和尽可能小的曲线,即是最佳拟合。最小二乘这个术语仅仅是使误差平方和最小的省略说法。 图11.1 2阶曲线拟合 在MATLAB中,函数polyfit求解最小二乘曲线拟合问题。为了阐述这个函数的用法, 让我们以上面图11.1中的数据开始。 ? x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1]; ? y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2]; 为了用polyfit,我们必须给函数赋予上面的数据和我们希望最佳拟合数据的多项式的 阶次或度。如果我们选择n=1作为阶次,得到最简单的线性近似。通常称为线性回归。相反,如果我们选择n=2 作为阶次,得到一个2阶多项式。现在,我们选择一个2阶多项式。 ? n=2; % polyno mial order ? p=polyfit(x, y, n)

数据插值与数据拟合

数据插值与数据拟合 1、一维数据插值: y=interp1(x0,y0,x,’method’) ‘method’共有四种方法选择: ‘nearest’ 最近点插值法取较近点的值 ‘linear’线性插值法用直线连接数据点 ‘spline’样条插值法用三次样条曲线通过数据点 ‘cubic’立方插值法用三次曲线通过数据点 例:对,,用个节点(等分)作上述四种插值,用m=21个插值点(等分)作图比较结果; 练习: 根据程序washu.可得,x=0:3的193个数据,即对应得y值,现在将 x=0:56图形形状不变 从而得到x=1:56的对应的y值,并且比较分析,哪一种插值效果好 2、 数据拟合 P=polyfit(x,y,n)返回系数从高到低 polyval(p,x) 例、在化工生产中获得的氯气的等级随生产时间下降。假定在时,与之间有如下形式的非线性模型: 现收集了44组数据: 80.49160.43280.41 80.49180.46280.40 100.48180.45300.40 100.47200.42300.40 100.48200.42300.38 100.47200.43320.41 120.46200.41320.40 120.46220.41340.40 120.45220.40360.41 120.43240.42360.38

140.45240.40380.40 140.43240.40380.40 140.43260.41400.39 160.44260.40420.39 160.43260.41 要求利用该数据求的值,以确定模型。 练习 问题1: N P K 施肥量(kg/ha) 产量 (t/ha) 施肥量 (kg/ha) 产量 (t/ha) 施肥量 (kg/ha) 产量 (t/ha) 0 34 67 101 135 202 259 336 404 47115.18 21.36 25.72 32.29 34.03 39.45 43.15 43.46 40.83 30.75 24 49 73 98 147 196 245 294 342 33.46 32.47 36.06 37.96 41.04 40.09 41.26 42.17 40.36 42.73 47 93 140 186 279 372 465 558 651 18.98 27.35 34.86 38.52 38.44 37.73 38.43 43.87 42.77 46.22 (1)将上面第一个表中以施肥量为自变量n,产量为函数y,用最小二乘法拟合函数,输出a1,b1,c1的值,给出拟合误差R^2,并进行图形比较(2) 将上面第二个表中以施肥量为自变量p,产量为函数y,用最小二乘法拟合函数,输出a,b的值,给出拟合误差R^2,并进行图形比较 (3) 将上面第三个表中以施肥量为自变量k,产量为函数y,用最小二乘法拟合函数,输出a3,b3,c3的值,给出拟合误差R^2,并进行图形比较Quadratic: Compound:

数值分析实验插值与拟合

《数值分析》课程实验一:插值与拟合 一、实验目的 1. 理解插值的基本原理,掌握多项式插值的概念、存在唯一性; 2. 编写MATLAB 程序实现Lagrange 插值和Newton 插值,验证Runge 现象; 3. 通过比较不同次数的多项式拟合效果,理解多项式拟合的基本原理; 4. 编写MATLAB 程序实现最小二乘多项式曲线拟合。 二、实验内容 1. 用Lagrange 插值和Newton 插值找经过点(-3, -1), (0, 2), (3, -2), (6, 10)的三次插值公式,并编写MATLAB 程序绘制出三次插值公式的图形。 2. 设 ]5,5[,11 )(2 -∈+= x x x f 如果用等距节点x i = -5 + 10i /n (i = 0, 1, 2, …, n )上的Lagrange 插值多项式L n (x )去逼近它。不妨取n = 5和n = 10,编写MATLAB 程序绘制出L 5(x )和L 10(x )的图像。 3. 在某冶炼过程中,根据统计数据的含碳量与时间关系如下表,试求含碳量与时间t 的拟合曲线。

(1) 用最小二乘法进行曲线拟合; (2) 编写MATLAB 程序绘制出曲线拟合图。 三、实验步骤 1. (1) Lagrange 插值法:在线性空间P n 中找到满足条件: ?? ?≠===j i j i x l ij j i , 0,, 1)(δ 的一组基函数{}n i i x l 0)(=,l i (x )的表达式为 ∏ ≠==--= n i j j j i j i n i x x x x x l ,0),,1,0()( 有了基函数{}n i i x l 0)(=,n 次插值多项式就可表示为 ∑==n i i i n x l y x L 0)()( (2) Newton 插值法:设x 0, x 1, …, x n 是一组互异的节点,y i = f (x i ) (i = 0, 1, 2, …, n ),f (x )在处的n 阶差商定义为 1102110] ,,,[],,,[],,,[x x x x x f x x x f x x x f n n n n --= - 则n 次多项式 ) ())(](,,[) )(](,,[)](,[)()(11010102100100----++--+-+=n n n x x x x x x x x x f x x x x x x x f x x x x f x f x N 差商表的构造过程:

相关文档