第5章 matlab 多项式、插值与数据拟合
Matlab中数据处理和多项式插值与曲线拟合
Matlab中数据处理和多项式插值与曲线拟合⼀、基本统计处理1、查取最⼤值MAX函数的命令格式有：[Y,I]= max (X)：将max(X)返回矩阵X的各列中的最⼤元素值及其该元素的位置赋予⾏向量Y与I；当X为向量时，则Y与I为单变量。
[Y,I]=max(X,[],DIM)：当DIM=1时按数组X的各列查取其最⼤的元素值及其该元素的位置赋予向量Y与I;当DIM=2时按数组X的各⾏查取其最⼤的元素值及其该元素的位置赋予向量Y与I.max(A,B)：返回⼀个与A,B同维的数组，其每⼀个元素是由A,B同位置上的元素的最⼤值组成。
【例1】查找下⾯数列x的最⼤值。
x=[3 5 9 6 1 8] % 产⽣数列xx = 3 5 9 6 1 8y=max(x) % 查出数列x中的最⼤值赋予yy = 9[y,l]=max(x) % 查出数列x中的最⼤值及其该元素的位置赋予y,ly = 9l = 3【例2】分别查找下⾯3×4的⼆维数组x中各列和各⾏元素中的最⼤值。
x=[1 8 4 2;9 6 2 5;3 6 7 1] % 产⽣⼆维数组xx = 1 8 4 29 6 2 53 6 7 1y=max(x) % 查出⼆维数组x中各列元素的最⼤值产⽣赋予⾏向量yy = 9 8 7 5[y,l]=max(x) % 查出⼆维数组x中各列元素的最⼤值及其这些% 元素的⾏下标赋予y,ly = 9 8 7 5l = 2 1 3 2[y,l]=max(x,[ ],1) % 本命令的执⾏结果与上⾯命令完全相同y = 9 8 7 5l = 2 1 3 2[y,l]=max(x,[ ],2) % 由于本命令中DIM=2,故查找操作在各⾏中进⾏y = 897l = 213[y,l]=max(x) % 查出⼆维数组x中各列元素的最⼤值及其这些% 元素的⾏下标赋予y,ly = 9 8 7 5l = 2 1 3 2[y,l]=max(x,[ ],1) % 本命令的执⾏结果与上⾯命令完全相同y = 9 8 7 5l = 2 1 3 2[y,l]=max(x,[ ],2) % 由于本命令中DIM=2,故查找操作在各⾏中进⾏y = 897l = 2132、查取最⼩值MIN函数⽤来查取数据序列的最⼩值。
matlab的数据拟合与插值
matlab的数据拟合与插值Matlab 的数据的分析处理-拟合与插值在数学建模过程中，常常需要确定⼀个变量依存于另⼀个或更多的变量的关系，即确定这些变量之间的函数关系。
但在实际中确定这些变量之间函数函数关系时往往没有先验的依据，只能在收集的实际数据的基础上对若⼲合乎理论的形式进⾏试验，从中选择⼀个最有可能反映实际的函数形式，这就是统计学中的拟合和回归⽅程问题。
本节我们主要介绍如何分析处理实际中得到的数据。
下⾯先看⼀个例⼦。
例1 “⼈⼝问题”是我国最⼤社会问题之⼀，估计⼈⼝数量和发展趋势是我们制定⼀系列相关政策的基础。
有⼈⼝统计年鉴，可查到我国从1949年⾄1994⼀般地，我们采⽤下⾯的分析处理⽅法：⾸先，在直⾓坐标系上作出⼈⼝数与年份的散点图象。
观察随着年份的增加⼈⼝数与年份变化关系，初步估计出他们之间的关系可近似地可看做⼀条直线。
那么我们如何把这条直线⽅程确定出来呢？并⽤他来估计1999年我国的⼈⼝数。
⽅法⼀：先选择能反映直线变化的两个点，如（1949，541.67），（1984，1034.75）⼆点确定⼀条直线，⽅程为 N = 14.088 t – 26915.842 ，代⼊t =1999,得N ≈12.46亿⽅法⼆：可以多取⼏组点对，确定⼏条直线⽅程，将t = 1999代⼊，分别求出⼈⼝数，在取其算数平值。
⽅法三：可采⽤“最⼩⼆乘法”求出直线⽅程。
这就是曲线拟合的问题。
⽅法⼀与⽅法⼆都具有⼀定的局限性，下⾯我们重点介绍数据的曲线拟合。
所谓曲线拟合是指给定平⾯上的n 个点(x i ,y i ),i=1,2,….,n,找出⼀条曲线使之与这些点相当吻合，这个过程称之为曲线拟合。
最常见的曲线拟合是使⽤多项式来作拟合曲线。
曲线拟合最常⽤的⽅法是最⼩⼆乘法。
其原理是求f(x),使21])([i ni i y x f -=∑=δ达到最⼩。
matlab 提供了基本的多项式曲线拟合函数命令polyfit格式::polyfit(x,y,n)说明:polyfit(x,y,n)是找n 次多项式p(x)的系数,这些系数满⾜在最⼩⼆乘法意义下p(x(i)) ~= y(i).已知⼀组数据，⽤什么样的曲线拟合最好呢?可以根据散点图进⾏直观观察，在此基础上，选择⼏种曲线分别拟合，然后⽐较，观察那条曲线的最⼩⼆乘指标最⼩。
MATLAB数据拟合与插值PPT课件
% linear interpolation » zcubic=interp2(width, depth, temps, wi,d, ' cubic ') ;
% cubic interpolation »plot(wi, zlinear, ' - ' , wi, zcubic)
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)
p =-9.8108 20.1293 -0.0317 polyfit 的输出是一个多项式系数的行向量。其解是
y = －9.8108x^2 ＋20.1293x－0.0317 为了将曲线拟合解与数据点比较，让我们把二者都绘成图。
» xi=linspace(0, 1, 100); % x-axis data for plotting » z=polyval(p, xi); » plot(x, y, ' o ' , x, y, xi, z, ' : ' ) » xlabel(' x '), ylabel(' y=f(x) '), title(' Second Order Curve Fitting ')
(完整版)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’——三次插值；默认为分段线性插值。
例1 从1点12点的11小时内，每隔1小时测量一次温度，测得的温度的数值依次为：5，8，9，15，25，29，31，30，22，25，27，24．试估计每隔1/10小时的温度值。
Matlab数据插值与拟合
end
end
end
第16页，共49页。
例4-3 根据下表的数据点求出其拉格朗日 插值多项式，并计算当x=1.6时y的值。
x
1
y 0.8415
1.2
0.9320
1.8
2. 5
0.9738 0.5985
4
-0.7568
解：
>> x=[1 1.2 1.8 2.5 4]; >> y=[0.8415 0.9320 0.9738 0.5985 -0.7568]; >> f=language(x,y)
同‘pchip’，三次Hermite多项式插值
第5页，共49页。
1.Linear（分段线性插值）
它 在的区算间法[xi是,xi在+1]每上个的小子区插间值[多xi,x项i+式1]上为采：用简单的线性插值。
Fi
x xi1 xi xi1
f
(xi )
x xi xi1 xi
f (xi1)
由此整个区间[xi,xi+1]上的插值函数为：
邻近的已知点的线性函数插值计算该区间内插值点上的函数
值。
第11页，共49页。
例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];
end;
%计算拉格朗日基函数
f = f + l; simplify(f);
%计算拉格朗日插值函数 %化简
if(i==n)
matlab插值与拟合
matlab插值与拟合
在MATLAB中，插值和拟合都是通过函数来实现的。
插值是通过创建新的数据点来填充在已知数据点之间的空白。
MATLAB提供了几种不同的插值方法，例如分段线性插值、三次样条插值、立方插值等。
具体使用哪种插值方法取决于数据的特性和所需的精度。
插值函数的一般形式是`interp1(x, y, xi, 'method')`，其中`x`和`y`是已知的数据点，`xi`是待插值点的横坐标向量，`method`是插值方法，例如最近邻点插值、线性插值、三次样条插值、立方插值等。
拟合是通过调整一个数学模型来使得该模型尽可能地接近给定的数据点。
在MATLAB中，可以使用`polyfit`函数进行多项式拟合。
该函数的一般形式是`p = polyfit(x, y, n)`，其中`x`和`y`是已知的数据点，`n`是多项式的阶数。
该函数返回一个向量`p`，表示多项式的系数。
可以使用`polyval`函数来评估这个多项式模型在给定数据点上的值。
需要注意的是，插值和拟合都是数学上的近似方法，它们只能尽可能地逼近真实的情况，而不能完全准确地描述数据的变化。
因此，选择合适的插值和拟合方法是非常重要的。
matlab插值与拟合(命令与示例)
目录【一维插值】interp1 (1)yi = interp1(x,y,xi,method) (1)例1 (1)例2 (2)【二维插值】interp2 (4)ZI = interp2(X,Y,Z,XI,YI,method) (4)插值方式比较示例 (4)例3 (8)例4 (9)【三角测量和分散数据插值】 (13)【数据拟合】 (17)例5 (17)例6 (18)【一维插值】interp1yi = interp1(x,y,xi,method)例1在1-12 的11 小时内，每隔1 小时测量一次温度，测得的温度依次为：5，8，9，15，25，29，31，30，22，25，27，24。
试估计每隔1/10 小时的温度值。
建立M文件temp.mhours=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,temps,'kp',h,t,'b');35302520151050 2 4 6 8 10 12例2已知飞机下轮廓线上数据如下，求x每改变0.1时的y值。
X0357911 12131415 Y0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6建立M文件plane.mx0=[0 3 5 7 9 11 12 13 14 15 ];y0=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6 ];x=0:0.1:15;y1=interp1(x0,y0,x,'nearest');y2=interp1(x0,y0,x);y3=interp1(x0,y0,x,'spline');plot(x0,y0,'kp',x,y1,'r')2.521.510.50 5 10 15 plot(x0,y0,'kp',x,y2,'r')2.521.510.50 5 10 15 plot(x0,y0,'kp',x,y3,'r')2.521.510.50 5 10 15 【二维插值】interp2ZI = interp2(X,Y,Z,XI,YI,method)插值方式比较示例用较大间隔产生peaks 函数数据点[x,y] = meshgrid(-3:1:3);z = peaks(x,y);surf(x,y,z)642-2-4-642 40 2-2 -2-4 -4●产生一个较好的网格[xi,yi] = meshgrid(-3:0.25:3);●利用最近邻方式插值zi1 = interp2(x,y,z,xi,yi,'nearest');surf(xi,yi,zi1)●双线性插值方式zi2 = interp2(x,y,z,xi,yi,'bilinear');surf(xi,yi,zi2)●双立方插值方式zi3 = interp2(x,y,z,xi,yi,'bicubic');surf(xi,yi,zi3)●不同插值方式构造的等高线图对比contour(xi,yi,zi1)321-1-2-3-3 -2 -1 0 1 2 3 contour(xi,yi,zi2)321-1-2-3-3 -2 -1 0 1 2 3 contour(xi,yi,zi3)321-1-2-3-3 -2 -1 0 1 2 3例3测得平板表面3*5 网格点处的温度分别为：82 81 80 82 8479 63 61 65 8184 84 82 85 86试作出平板表面的温度分布曲面z=f(x,y)的图形。
matlab中插值与拟合
计算可视化
1 插值与数据拟合 1.1 一维数据的插值问题 1.1.1 一维插值问题的求解
例1-1：已知的数据点来自函数25f ()(35)sin x x x x e x -=-+，根据生成的数据进行插值处理，得出较平滑的曲线直接生成数据。
例1-2：编写一段程序，允许利用插值方法手工绘制一条光滑的曲线。
1.1.2 Lagrange插值算法及应用
例1-3：对2()1/(125),11f x x x =+-≤≤进行Lagrange 插值。
1.2 已知样本点的定积分计算
例1-4：利用样条插值算法求解3/2
cos(15)x dx π⎰。
例1-5：已知其中的150个数据点，用quadspln()计算出该定积分的值
2
()x
t erf x e dt -=
⎰。
1.3 二维网格数据的插值问题
例1-6：由2
22(,)(2)x
y xy
z f x y x x e ---==-可计算出一些较稀疏的网格数据，对整
1.4 二维一般分布数据的插值问题1.5 高维插值问题。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
两个多项式的和与差:
ya = a1 x + a2 x
m
m 1 n 1
+ + am x + am +1 + + bn x + bn +1
yb = b1 x + b2 x
n
命令poly_add:求两个多项式的和,其调用格式为: c= poly_add(a,b) 多项式a减去b,可表示为: c= poly_add(a,-b)
功能:求多项式积分 调用格式:py=poly_itg(p) p:被积多项式的系数 py:求积后多项式的系数 poly_itg.m function py=poly_itg(p) n=length(p); py=[p.*[n:-1:1].^(-1),0] 不包括最后一项积分常数
多项式微分:
y = c1 x + c2 x
例: y = (x 1)6 = x6 6x5 +15x4 20x3 +15x2 6x +1
>> r=roots([1 -6 15 -20 15 -6 1]) r= 1.0042 + 0.0025i 1.0042 - 0.0025i 1.0000 + 0.0049i 1.0000 - 0.0049i 0.9958 + 0.0024i 0.9958 - 0.0024i 舍入误差的影响,与计算精度有关.
polyval: 可用命令polyval计算多项式的值. 例: y = 3x4 7x3 +2x2 + x +1 计算y(2.5) >> c=[3,-7,2,1,1]; xi=2.5; yi=polyval(c,xi) yi = 23.8125 如果xi是含有多个横坐标值的数组,则yi也 为与xi长度相同的向量. >> c=[3,-7,2,1,1]; xi=[2.5,3]; >> yi=polyval(c,xi) yi = 23.8125 76.0000
y = 0.2015x3 +1.4385x2 2.7477x + 5.4370
多项式积分:
y = c1 x n + c 2 x n 1 + + c n x + c n + 1
Y =
∫
cn 2 c1 n +1 c 2 n ydx = x + x + + x + c n +1 x + c n + 2 n +1 n 2
m阶多项式与n阶多项式的乘积是d=m+n阶的多项式:
y a = a1 x + a 2 x
m
m 1
+ + a m x + a m +1
y b = b1 x n + b 2 x n 1 + + b n x + b n + 1
yc = ya yb = c1 x + c2 x
d
d 1
+ + cd x + cd +1
为解决Rung问题,引入分段插值.
5.2.4 分段插值
算法分析:所谓分段插值就是通过插值点用折 线或低次曲线连接起来逼近原曲线. MATLAB实现 可调用内部函数.
– 命令1 interp1
功能 : 一维数据插值(表格查找).该命令对数据点之 间计算内插值.它找出一元函数f(x)在中间点的数值.其 中函数f(x)由所给数据决定. 格式1 yi = interp1(x,Y,xi) %返回插值向量yi,每一元素对应于参量xi,同时由向 量x与Y的内插值决定.参量x指定数据Y的点.若Y为一矩阵, 则按Y的每列计算. 算例 对于t,beta ,alpha分别有两组数据与之对应,用分段线 性插值法计算当t=321, 440, 571时beta ,alpha的值.
polyfit:给定n+1个点将可以唯一确定一个n阶多项式.利 用命令polyfit可容易确定多项式的系数. 例: >> x=[1.1,2.3,3.9,5.1]; >> y=[3.887,4.276,4.651,2.117]; >> a=polyfit(x,y,length(x)-1) a= -0.2015 1.4385 -2.7477 5.4370 >> poly2sym(a) ans = -403/2000*x^3+2877/2000*x^2-27477/10000*x+5437/1000 多项式为 Polyfit的第三个参数是多项式的阶数.
)
j =1 j≠k
MATLAB实现
function y=lagrange(x0,y0,x) n n x xj ) ii=1:length(x0); y=zeros(size(x)); y ( x ) = ∑ y k (∏ k =1 j =1 x k x j for i=ii j≠k ij=find(ii~=i); y1=1; for j=1:length(ij), y1=y1.*(x-x0(ij(j))); end y=y+y1*y0(i)/prod(x0(i)-x0(ij)); end 算例:给出f(x)=ln(x)的数值表,用Lagrange计算 ln(0.54)的近似值. >> x=[0.4:0.1:0.8]; >> y=[-0.916291,-0.693147,-0.510826,-0.356675,-0.223144]; >> lagrange(x,y,[0.54,0.55,0.78]) ans = -0.6161 -0.5978 -0.2484 ( 精确解-0.616143)
例: >> r=roots(p) 得到 r= 0.2500 + 1.5612i 0.2500 - 1.5612i -1.0000 所有零点由一个列向量给出.
Poly: 由零点可得原始多项式的各系数,但可能相差 一个常数倍. 例: >> poly(r)
ans = 1.0000 0.5000 2.0000 2.5000 注意:若存在重根,这种转换可能会降低精度.
幂系数:在MATLAB里,多项式用行向量表示,其 元素为多项式的系数,并从左至右按降幂排列.
例:
y = 2x + x + 4x + 5
3 2
被表示为 >> p=[2 1 4 5] >> poly2sym(p) ans = 2*x^3+x^2+4*x+5
Roots: 多项式的零点可用命令roots求的.
计算 yc 系数的MATLAB命令是:c=conv(a,b) 多项式 yb 除多项式 ya 的除法满足:
y a = y q yb + y r
其中 yq 是商, yr 是除法的余数.多项式 yq 和 yr 可由命令deconv算出. 例:[q, r]=deconv(a,b)
例 >> a=[2,-5,6,-1,9]; b=[3,-90,-18]; >> c=conv(a,b) c= 6 -195 432 -453 9 -792 -162 >> [q,r]=deconv(c,b) q= 2 -5 6 -1 9 r= 0 0 0 0 0 0 0 >> poly2sym(c) ans = 6*x^6-195*x^5+432*x^4-453*x^3+9*x^2-792*x-162
>> x=[-5:1:5]; y=1./(1+x.^2); x0=[-5:0.1:5]; >> y0=lagrange(x,y,x0); >> y1=1./(1+x0.^2); %绘制图形 >> plot(x0,y0,'--r') %插值曲线 >> hold on >> plot(x0,y1,'-b') %原曲线
5.2.2 Hermite插值
方法介绍 不少实际问题不但要求在节点上函数值相等,而且 要求导数值也相等,甚至要求高阶导数值也相等,满足 这一要求的插值多项式就是Hermite插值多项式.下面 只讨论函数值与一阶导数值个数相等且已知的情况. 已知n个插值点 x1, x2 ,, xn 及对应的函数值 y1, y2 ,, yn 和一阶导数值 y1' , y2' ,, y'n .则对插值区间 内任意x的函数值y的Hermite插值公式:
n
n 1
+ + cn x + cn +1
n2
y = nc1 x
'
n 1
+ ( n 1)c2 x
+ + cn
Polyder: 求多项式一阶导数的系数. 调用格式为: b=polyder(c ) c为多项式y的系数,b是微分后的系数, 其值为:
[ nc1 , ( n 1) c 2 , , c n ]
功能:两个多项式相加 调用格式:b=poly_add(p1,p2) b:求和后的系数数组
poly_add.m function p3=poly_add(p1,p2) n1=length(p1); n2=length(p2); if n1==n2 p3=p1+p2;end if n1>n2 p3=p1+[zeros(1,n1-n2),p2];end if n1<n2 p3=[zeros(1,n2-n1),p1]+p2;end
一个多项式的幂级数形式可表示为:
y = c1 x + c2 x
n
n 1
+ + cn x + cn+1