第5讲MATLAB多项式及插值.ppt

第第5讲讲 MATLAB数据分析与多项式计算数据分析与多项式计算5.1 数据统计处理数据统计处理5.2 多项式计算多项式计算5.3 插值运算插值运算5.4 多项式拟合多项式拟合15.1 数据统计处理数据统计处理5.1.1 最大值和最小值最大值和最小值 MATLAB提供的求数据序列的最大值和最小值的函数分别为max和min，两个函数的调用格式和操作过程类似1．求向量的最大值和最小值求一个向量X的最大值的函数有两种调用格式，分别是：(1) y=max(X)：返回向量X的最大值存入y，如果X中包含复数元素，则按模取最大值2) [y,I]=max(X)：返回向量X的最大值存入y，最大值的序号存入I，如果X中包含复数元素，则按模取最大值求向量X的最小值的函数是min(X)，用法和max(X)完全相同2例5-1 求向量x的最大值命令如下：x=[-43,72,9,16,23,47];y=max(x) %求向量x中的最大值[y,l]=max(x) %求向量x中的最大值及其该元素的位置32．求矩阵的最大值和最小值．求矩阵的最大值和最小值求矩阵A的最大值的函数有3种调用格式，分别是：(1) max(A)：返回一个行向量，向量的第i个元素是矩阵A的第i列上的最大值。

2) [Y,U]=max(A)：返回行向量Y和U，Y向量记录A的每列的最大值，U向量记录每列最大值的行号3) max(A,[],dim)：dim取1或2dim取1时，该函数和max(A)完全相同；dim取2时，该函数返回一个列向量，其第i个元素是A矩阵的第i行上的最大值求最小值的函数是求最小值的函数是min，其用法和，其用法和max完全相同完全相同4例5-2 分别求3×4矩阵x中各列和各行元素中的最大值，并求整个矩阵的最大值和最小值53．两个向量或矩阵对应元素的比较．两个向量或矩阵对应元素的比较 函数max和min还能对两个同型的向量或矩阵进行比较，调用格式为：(1) U=max(A,B)：A,B是两个同型的向量或矩阵，结果U是与A,B同型的向量或矩阵，U的每个元素等于A,B对应元素的较大者2) U=max(A,n)：n是一个标量，结果U是与A同型的向量或矩阵，U的每个元素等于A对应元素和n中的较大者min函数的用法和max完全相同6例5-3 求两个2×3矩阵x, y所有同一位置上的较大元素构成的新矩阵p75.1.2 求和与求积求和与求积 数据序列求和与求积的函数是sum和prod，其使用方法类似。

设X是一个向量，A是一个矩阵，函数的调用格式为：sum(X)：返回向量X各元素的和prod(X)：返回向量X各元素的乘积sum(A)：返回一个行向量，其第i个元素是A的第i列的元素和prod(A)：返回一个行向量，其第i个元素是A的第i列的元素乘积sum(A,dim)：当dim为1时，该函数等同于sum(A)；当dim为2时，返回一个列向量，其第i个元素是A的第i行的各元素之和prod(A,dim)：当dim为1时，该函数等同于prod(A)；当dim为2时，返回一个列向量，其第i个元素是A的第i行的各元素乘积8例5-4 求矩阵A的每行元素，每列元素的乘积和全部元素的乘积95.1.3 平均值和中值平均值和中值 求数据序列平均值的函数是mean，求数据序列中值的函数是median两个函数的调用格式为：mean(X)：返回向量X的算术平均值median(X)：返回向量X的中值mean(A)：返回一个行向量，其第i个元素是A的第i列的算术平均值median(A)：返回一个行向量，其第i个元素是A的第i列的中值mean(A,dim)：当dim为1时，该函数等同于mean(A)；当dim为2时，返回一个列向量，其第i个元素是A的第i行的算术平均值。

median(A,dim)：当dim为1时，该函数等同于median(A)；当dim为2时，返回一个列向量，其第i个元素是A的第i行的中值10例例5-5 5-5 别求别求4 4元素向量元素向量x x的平均值和中值的平均值和中值115.1.4 排序排序 MATLAB中对向量X是排序函数是sort(X)，函数返回一个对X中的元素按升序排列的新向量 sort函数也可以对矩阵A的各列或各行重新排序，其调用格式为： [Y,I]=sort(A,dim)其中dim指明对A的列还是行进行排序若dim=1，则按列排；若dim=2，则按行排Y是排序后的矩阵，而I记录Y中的元素在A中位置1213一个多项式的幂级数形式可表示为：一个多项式的幂级数形式可表示为：同时也可表为嵌套形式同时也可表为嵌套形式或因子形式或因子形式 N N阶多项式阶多项式n n个根，其中包含重根和复根若多项式所有个根，其中包含重根和复根。

若多项式所有系数均为实数，则全部复根都将以共轭对的形式出现系数均为实数，则全部复根都将以共轭对的形式出现 5.2 多项式运算多项式运算1、多项式表示142 2、相关函数、相关函数、相关函数、相关函数幂系数:在MATLAB里，多项式用行向量表示，其元素为多项式的系数，并从左至右按降幂排列 例： 被表示为 >> p=[2 1 4 5] >> poly2sym(p) ans = 2*x^3+x^2+4*x+5Roots: 多项式的零点可用命令roots求的。

例： >> r=roots(p) 得到 r = 0.2500 + 1.5612i 0.2500 - 1.5612i -1.0000 所有零点由一个列向量给出15Poly: 由零点可得原始多项式的各系数，但可能相差一个常数倍例：根据上例： >> poly(r) ans = 1.0000 0.5000 2.0000 2.5000注意：若存在重根，这种转换可能会降低精度例： 舍入误差的影响，与计算精度有关16polyval: polyval: 可用命令可用命令polyvalpolyval计算多项式的值计算多项式的值例：例： 计算计算y(2.5)y(2.5) >> c=[3,-7,2,1,1]>> c=[3,-7,2,1,1]；； xi=2.5; yi=polyval(c,xi) xi=2.5; yi=polyval(c,xi) yi = yi = 23.8125 23.8125 如果xi是含有多个横坐标值的数组，则yi也为与xi长度相同的向量。

>> c=[3,-7,2,1,1]; xi=[2.5,3];>> yi=polyval(c,xi)yi = 23.8125 76.0000175.35.3插值运算插值运算1、Lagrange插值方法介绍 对给定的n个插值点 及对应的函数值 ，利用n次Lagrange插值多项式，则对插值区间内任意x的函数值y可通过下式求的： MATLAB中没有直接实现拉格朗日算法的函数，我们已经介绍过该函数的书写：18function y=lagrange (a,b,x)y=0;for i=1:length(a) l=1; for j=1:length(b) if j==i l=l; else l=l.*(x-a(j))/(a(i)-a(j)); end end y=y+l*b(i);end19算例：给出算例：给出f(x)=ln(x)f(x)=ln(x)的数值表，用的数值表，用LagrangeLagrange计算计算ln(0.54)ln(0.54)的近似值。

的近似值>> x=[0.4:0.1:0.8]>> x=[0.4:0.1:0.8]；；>> y=[-0.916291,-0.693147,-0.510826,-0.356675,-0.223144];>> y=[-0.916291,-0.693147,-0.510826,-0.356675,-0.223144];>> lagrange(x,y,0.54)>> lagrange(x,y,0.54)ans =ans = -0.6161 -0.6161 （精确解（精确解-0.616143)-0.616143)202 2、、、、RungeRunge现象和分段插值现象和分段插值现象和分段插值现象和分段插值n n问题的提出：根据区间问题的提出：根据区间[a,b][a,b]上给出的节点做插值多项式上给出的节点做插值多项式p(x)p(x)的近似值，一般总认为的近似值，一般总认为p(x)p(x)的次数越高则逼近的次数越高则逼近f(x)f(x)的精度就越的精度就越好，但事实并非如此。

好，但事实并非如此n n反例：反例： 在区间在区间[-5,5][-5,5]上的各阶导数存在，但在此区间上取上的各阶导数存在，但在此区间上取n n个节个节点所构成的点所构成的LagrangeLagrange插值多项式在全区间内并非都收敛插值多项式在全区间内并非都收敛 取取n=10,n=10,用用LagrangeLagrange插值法进行插值计算插值法进行插值计算21>> x=[-5:1:5]; y=1./(1+x.^2); x0=[-5:0.1:5];>> x=[-5:1:5]; y=1./(1+x.^2); x0=[-5:0.1:5];>> y0=lagrange(x,y,x0);>> y0=lagrange(x,y,x0);>> y1=1./(1+x0.^2);>> y1=1./(1+x0.^2);％绘制图形％绘制图形>> plot(x0,y0,'--r')>> plot(x0,y0,'--r')％插值曲线％插值曲线>> hold on>> hold on>> plot(x0,y1,‘-b')>> plot(x0,y1,‘-b')％原曲线％原曲线n n为解决Rung问题，引入分段插值。

22n n算法分析：所谓分段插值就是通过插值点用折线或低次曲线算法分析：所谓分段插值就是通过插值点用折线或低次曲线连接起来逼近原曲线连接起来逼近原曲线n nMATLABMATLAB实现实现 可调用内部函数可调用内部函数ØØ命令命令 interp1interp1功能 ： 一维数据插值（表格查找）该命令对数据点之间计算内插值它找出一元函数f(x)在中间点的数值其中函数f(x)由所给数据决定 格式：yi = interp1(x,Y,xi,method) %用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值 ’cubic’： 分段三次Hermite插值 23t = 1900:10:1990;p = [75.995 91.972 105.711 123.203 131.669... 150.697 179.323 203.212 226.505 249.633];对应于美国从1900年到1990年的每10年的人口数，求1975年的人口。

由此推断美国1900年到2000年每一年的人口数，并画出图形与1阶拉格朗日算法比较：lagrand([1970 1980],[203.212 226.505],1975)ans = 214.858524推广到多个点计算：>> t = 1900:10:1990;>> p = [75.995 91.972 105.711 123.203 131.669...150.697 179.323 203.212 226.505 249.633];>> x=[1965 1975]x = 1965 1975>> y=interp1(t,p,x)y = 191.2675 214.8585lagrand([1960 1970],[179.323 203.212],1965)ans = 191.2675lagrand([1970 1980],[203.212 226.505],1975)ans = 214.858525 例例例例>> x0=-1+2*[0:10]/10; >> y0=1./(1+25*x0.^2);>> x=-1:.01:1; >> y=lagrange(x0,y0,x); % Lagrange 插值>> ya=1./(1+25*x.^2); >> plot(x,ya,x,y,':')26>> y1=interp1(x0,y0,x);>> y1=interp1(x0,y0,x);>> plot(x,ya,x,y1);>> plot(x,ya,x,y1);27 用插值的方法对一函数进行近似用插值的方法对一函数进行近似, ,要求所得到的插值要求所得到的插值多项式经过已知插值节点多项式经过已知插值节点; ;在在n n比较大的情况下比较大的情况下, ,插值多项插值多项式往往是高次多项式式往往是高次多项式, ,这也就容易出现振荡现象（龙格现这也就容易出现振荡现象（龙格现象），即虽然在插值节点上没有误差象），即虽然在插值节点上没有误差, ,但在插值节点之外但在插值节点之外插值误差变得很大插值误差变得很大, ,从从“ “整体整体” ”上看上看, ,插值逼近效果将变得插值逼近效果将变得“ “很差很差” ”。

所谓数据拟合是求一个简单的函数,例如是一个低次多项式,不要求通过已知的这些点,而是要求在整体上“尽量好”的逼近原函数这时,在每个已知点上就会有误差,数据拟合就是从整体上使误差,尽量的小一些5.4 5.4 数据最小二乘拟合数据最小二乘拟合数据最小二乘拟合数据最小二乘拟合28多项式最小二乘拟合多项式最小二乘拟合多项式最小二乘拟合多项式最小二乘拟合n n次多项式：次多项式：下面我们练习一个具体例题 某种铝合金的含铝量为x%，其溶解温度为y摄氏度，由试验测得的x与y的数据表如下，试用最小二乘算法建立x,y的经验公式i ix xi iy yi i1 136.936.91811812 246.746.71971973 363.763.72352354 477.877.82702705 584842832836 6878729229229（1）根据所给出数据，画出图形，观察数据关系（2）通过图形，我们可以确定该曲线可以通过一次方程y=ax+b的形式来进行拟合可以确定，n=1，L=630（3）建立法方程组带入具体数据>> A=[6 396.9; 396.9 28365.28];>> B=[1458;101176.3];>> X=inv(A)*BX = 94.7479 2.2412 则 a=94.7479，b=2.2414拟合方程为y=94.7479 + 2.2414x31画出比较图形：x=[36.9 46.7 63.7 77.8 84 87.5];y=[181 197 235 270 283 292];plot(x,y,x,y,'o') hold on x1=[35:5:90]; y1=2.2412*x1+ 94.7479;plot(x1,y1,'r')32而该拟合公式在而该拟合公式在而该拟合公式在而该拟合公式在matlabmatlabmatlabmatlab中有内部函数：中有内部函数：中有内部函数：中有内部函数：多项式拟合多项式拟合多项式拟合多项式拟合MATLABMATLABMATLABMATLAB命令：命令：命令：命令：polyfitpolyfitpolyfitpolyfit格式：格式：格式：格式：p=polyfit(x,y,n)p=polyfit(x,y,n)p=polyfit(x,y,n)p=polyfit(x,y,n)其中，xy为原始样本点构成的向量 n为选定的多项式阶数 P为多项式系数（降幂排列）则上例采用命令求解：>> x=[36.9 46.7 63.7 77.8 84 87.5];>> y=[181 197 235 270 283 292];>> a=polyfit(x,y,1)a = 2.2337 95.352433>> x0=0:.1:1; >> y0=(x0.^2-3*x0+5).*exp(-5*x0).*sin(x0);>> p3=polyfit(x0,y0,3); %先进行三次拟合p3 = 2.8400 -4.7898 1.9432 0.0598多项式如下：2.8400*x^3-4.7898*x^2+1.9432*x+0.0598 已知数据点来自 用多项式拟合方法在不同阶次下进行拟合。

拟合该数据的程序34绘制拟合曲线：>> x=0:.01:1; >> ya=(x.^2-3*x+5).*exp(-5*x).*sin(x); %原函数数据>> y1=polyval(p3,x); %获得新拟合曲线数据>> plot(x,y1,x,ya,x0,y0,'o')35就不同的次数进行拟合：>> p4=polyfit(x0,y0,4); y2=polyval(p4,x);>> p5=polyfit(x0,y0,5); y3=polyval(p5,x);>> p8=polyfit(x0,y0,8); y4=polyval(p8,x);>> plot(x,ya,x0,y0,'o',x,y2,x,y3,x,y4)36。