[MATLAB基础] MATLAB应用大全 书连载

[lili123]

第5章 MATLAB数据分析
针对数据分析和处理，MATLAB提供了大量的函数，非常方便和灵活。本章将详细介绍如何利用MATLAB进行一些基本的数据分析，主要包括多项式及其函数，插值，以及函数的极限。MATLAB能够很好地解决多项式运算问题，这些函数用于多项式的建立、多项式求值、多项式乘法和除法、多项式求导和多项式展开和拟合等。插值函数，主要包括一维插值、二维插值、样条插值和高维插值等内容。MATLAB提供了非常强大的函数，可以非常方便和灵活地求得函数的极限。
5.1 多项式及其函数
MATLAB提供了一些处理多项式的专用函数，用户可以很方便地进行多项式的建立、多项式求值、乘法和除法运算，以及求多项式的导数和微分、多项式的根、多项式的展开和拟合等。
5.1.1 多项式的建立
MATLAB语言中，对于多项式 ，用多项式的系数按照降幂次序存放在向量 中。顺序必须是从高到低进行排列。例如，多项式 可以用系数向量 来表示。对多项式的操作就转换为对多项式系数向量的操作问题， 次多项式用一个 维的行向量表示，在多项式中缺少的幂次要用“0”来补齐。在MATLAB中，可以采用直接输入多项式系数建立多项式，也可以采用多项式的根来建立多项式，下面分别进行介绍。
1．直接输入多项式系数法
MATLAB中多项式是以向量的形式存储的，输入向量后，MATLAB将按照降幂顺序自动把向量的元素分配给多项式各项的系数。该向量可以是行向量，也可以是列向量。
【例5-1】 使用向量来创建多项式 ，并进行显示。
首先创建系数向量，然后使用函数poly2sym()将该向量转换为多项式，直接进行输出，也可以采用函数disp()进行多项式的显示，代码如下：

>> clear all;
p1=[4 3 2 1];
y=poly2sym(p1) %由向量创建多项式
disp(y) %显示多项式

运行程序，输出结果如下：

y =
4*x^3+3*x^2+2*x+1
4*x^3+3*x^2+2*x+1

在MATLAB中，多项式的输出采用降幂形式。如果多项式的系数为0，则在输出的时候不显示，可采用函数disp()来显示多项式。
2．由多项式的根来建立多项式
如果多项式的根已知，可以采用函数poly()建立多项式。其调用格式为p=poly(r)，其中r为由多项式的根组成的向量，p为输出的多项式的系数向量。
【例5-2】 已知某多项式的根为2、3、4，求该多项式。代码如下：

>> r=[2 3 4];
p=poly(r); %由多项式的根创建多项式
y=poly2sym(p) %显示多项式

运行程序后，输出结果如下：

y =
x^3-9*x^2+26*x-24

在程序中，函数poly()通过多项式的根创建多项式，多项式的根向量为r。通过函数poly2sym()对创建的多项式进行显示。

[lili123]

5.1.2 多项式的求值与求根
在MATLAB中，通过函数polyval()和polyvalm()可以对多项式进行求值，两者的区别为前者是代数多项式求值，后者是矩阵多项式求值。通过函数roots()对多项式进行求根，如果已经知道多项式的根，也可以通过函数poly()建立相对应的多项式。下面分别进行 介绍。
1．多项式求值
在MATLAB中，提供了两个函数对多项式进行求值，函数polyval()和polyvalm()。函数polyval()以数组或矩阵中的元素为计算单位，函数polyvalm()以矩阵为计算单位。
函数polyval()的调用方式为y=polyval(p, x)，其中参数p为行向量形式的多项式，参数x为带入多项式的值。参数x可以是标量，也可以是向量和矩阵。如果参数x为向量或矩阵，则该函数对向量或矩阵中的每一个元素计算多项式的值，其返回值y也分别为向量或矩阵。
【例5-3】 求多项式 在点2、3、4、5的值。代码如下：

>> clear all;
p=[1 -1 -6]; %多项式的系数向量
x=2:5;
y=polyval(p,x) %求多项式在x处的值

运行程序后输出结果如下：

y =
-4 0 6 14

函数polyvalm()的调用方式为y=polyvalm(p, x)，其中参数p为行向量形式的多项式，参数x必须为方阵，输出结果仍然为方阵。
【例5-4】 求 时，多项式 的值。
利用函数polyvalm()时，输出结果计算公式为 ，常数值变为常数乘以相同阶次的单位矩阵。如果采用函数polyval()，计算矩阵中每个元素对应的多项式的值。下面分别采用函数polyvalm()和polyval()对该矩阵进行计算，代码如下：

>> x=[1 2 ;3 4];
p=[2 3 4];
y1=polyvalm(p,x) %采用polyvalm()函数，以矩阵为计算单位
x=[1 2 ;3 4];
p=[2 3 4];
y2=polyval(p,x) %采用polyval()函数，以矩阵的元素为计算单位

运行程序后，输出结果如下：

y1 =
21 26
39 60
y2 =
9 18
31 48

当采用函数polyval()时，虽然输入参数是矩阵，但在计算时将矩阵的每个元素分别输入多项式，计算对应的输出。
2．多项式求根
在MATLAB中，利用roots()函数来求多项式的根，其调用格式为x=roots(p)，其中参数p为多项式系数，输出参数x为多项式的根。对于n次多项式具有n个根，这些根可能是实根，也可能是共轭复根。在MATLAB中，如果已经知道多项式的根，可以利用函数poly()求多项式的系数，其调用格式为y=poly(x)，输入参数x为根，输出参数y为得到的多项式系数向量。
【例5-5】 求多项式 的根，以及以4和5为根的多项式。代码如下：

>> clear all;
p=[1 0 0 -1 -6];
x1=roots(p) %对多项式p求根
x2=[4 5];
y=poly(x2); %求以x2为根的多项式
y=poly2sym(y)

运行程序后，输出结果如下：

x1 =
1.6638
-0.1021 + 1.5684i
-0.1021 - 1.5684i
-1.4597
y =
x^2 - 9*x + 20

利用函数roots()计算多项式的根，非常方便，函数的返回值x是一个向量，其长度等于多项式的根的个数。

[lili123]

5.1.3 多项式乘法和除法
在MATLAB中，使用函数conv()对多项式进行乘法运算。其调用格式为c=conv(a, b)，a和b为多项式的系数向量，该函数实现向量a和b的卷积，在代数上相当于多项式a乘以多项式b，其中c为相乘所产生的多项式的系数向量。
【例5-6】 求多项式 和 的乘积。采用函数conv()实现，其MATLAB程序如下：

p1=[4 2 0 5]; %缺少的幂次用0补齐
p2=[5 8 1];
y1=poly2sym(p1)
y2=poly2sym(p2)
p3=conv(p1,p2); %多项式相乘
y=poly2sym(p3)

运行程序后，输出结果如下：

y1 =
4*x^3+2*x^2+5
y2 =
5*x^2+8*x+1
y =
20*x^5+42*x^4+20*x^3+27*x^2+40*x+5

在MATLAB中，采用poly2sym()函数将向量作为多项式的系数进行输出，和其相对应的函数是sym2poly()，该函数将输入多项式的系数提取出来，作为向量进行输出。对于例5-6程序也可以用poly2sym()函数实现（和例5-6的计算结果相同），其MATLAB程序如下：

>> syms x
p1=sym2poly(4*x^3+2*x^2+5)
p2=sym2poly(5*x^2+8*x+1)
p3=conv(p1,p2); %多项式相乘
y=poly2sym(p3)

运行程序后，输出结果如下：

p1 =
4 2 0 5
p2 =
5 8 1
y =
20*x^5+42*x^4+20*x^3+27*x^2+40*x+5

在MATLAB中，使用函数deconv()对多项式进行除法运算。其调用格式为[q, r] = deconv(a, b)，实现解卷积运算。其中a和b为多项式的系数向量，在代数上相当于多项式a除以b，得到的商为q和余多项式r，它们之间的关系为a = conv(b, q) + r。
【例5-7】 求多项式 除以多项式 的商和余数，代码如下：

>> p1=[4 3 8 1 4];
p2=[2 3 1];
[q,r]=deconv(p1,p2); %多项式p1除以p2
y1=poly2sym(q) %商
y2=poly2sym(r) %余数

运行程序后，输出结果如下：

y1 =
2*x^2-3/2*x+21/4
y2 =
-53/4*x-5/4

5.1.4 多项式的导数和积分
在MATLAB中，通过函数polyder()和polyint()分别对多项式进行求导和积分。求导和积分互为逆运算，如果先对多项式进行积分，然后再求导，结果仍然为原来的多项式。下面对多项式的求导和积分分别进行讲解。
1．多项式的导数
在MATLAB中，采用函数polyder()进行多项式的求导，调用方式如下。
 y=polyder(p)：对以向量p为系数的多项式求导。
 y=polyder(a, b)：对以a和b为系数的多项式乘积进行求导。
 [q,d]=polyder(b, a)：返回以b为系数的多项式除以以a为系数的多项式的商的导数，并以q/d格式表示。
【例5-8】 对多项式求导，其MATLAB程序如下：

>> p1=[4 3 2];
p2=[2 2 1];
y1=polyder(p1); %对多项式p1求导
y1=poly2sym(y1)
y2=polyder(p1,p2); %对多项式p1和p2的乘积求导
y2=poly2sym(y2)
[q,d]=polyder(p1,p2); %对多项式p1除以p2的商求导
q=poly2sym(q)
d=poly2sym(d)

运行程序后，输出结果如下：

y1 =
8*x + 3
y2 =
32*x^3 + 42*x^2 + 28*x + 7
q =
2*x^2 - 1
d =
4*x^4 + 8*x^3 + 8*x^2 + 4*x + 1

在MATLAB中，通过函数polyder()对多项式进行求导，通过对输入参数和输出参数个数的不同，对相对应的多项式进行求导计算。对于函数[q,d]=polyder(b, a)，相当于对多项式 求导，结果为 。
2．多项式的积分
在MATLAB中，使用函数polyint()对多项式进行积分运算，其调用方式如下。
 polyint(p, k)：返回以向量p为系数的多项式的积分，积分的常数项为k。
 polyint(p)：返回以向量p为系数的多项式的积分，积分的常数项为默认值0。
【例5-9】 对多项式 进行积分运算，其常数项分别为3和0，其实现的MATLAB程序代码如下：

p1=[3 2 2];
y1=polyint(p1,3); %对多项式p1进行积分，常数项为3
y1=poly2sym(y1)
y2=polyint(p1); %对多项式p1进行积分，常数项为0
y2=poly2sym(y2)

运行程序后，输出结果如下：

y1 =
x^3 + x^2 + 2*x + 3
y2 =
x^3 + x^2 + 2*x

通过polyint()函数对多项式进行积分运算，积分的常数项通过参数k进行设置。如果不对参数k进行设置，则k取默认值0。

[lili123]

5.1.5 多项式展开
在MATLAB中，有理多项式用它们的分子多项式和分母多项式进行表示，函数residue()可以将多项式之比用部分分式展开，也可以将一个部分分式用多项式之比进行表示。函数residue()的调用方式如下。
 [r, p, k]=residue(b, a)：求多项式之比b/a的部分分式展开，函数的返回值r是余数，p是部分分式的极点，k是常数项。如果多项式a没有重根，部分分式展开的形式如下：

其中向量r、p的长度和向量a、b的长度有如下关系：

当向量b的长度小于a时，向量k中没有元素，否则应满足：

 [b, a]=residue(r, p, k)：通过部分分式得到多项式，该多项式的形式为b/a。
【例5-10】 将多项式 和 展开成几个简单多项式的和。其实现的MATLAB代码如下：

>> clear all;
clear all;
b=[1 -1 -7 -1]; %分子多项式
a=poly([1;5;6]); %分母多项式
[r,p,k]=residue(b,a) %进行多项式b/a展开
[b1,a1]=residue(r,p,k); %通过余数、极点和常数项来求多项式b1/a1
b1=poly2sym(b1)
a1=poly2sym(a1)
b=[1 -1 -7 -1]; %多项式a有三重根
a=poly([1;1;1]); %分母多项式
[r,p,k]=residue(b,a) %展开多项式b/a

运行程序后，输出结果如下：

r =
27.4000
-16.0000
-0.4000
p =
6.0000
5.0000
1.0000
k =
1
b1 =
x^3 - x^2 - 7*x - 1
a1 =
x^3 - 12*x^2 + 41*x - 30
r =
2.0000
-6.0000
-8.0000
p =
1.0000
1.0000
1.0000
k =
1

利用函数[r, p, k]=residue(b, a)将多项式b/a进行展开，结果为余数、极点和常数项。对该多项式进行展开后的结果如下：


将余数、极点和常数项带入函数[b, a]=residue(r, p, k)中，可以求得对应的多项式，并通过b/a的形式给出。
当多项式a有三重根1时，对多项式进行展开后的结果如下：

[lili123]

5.1.6 多项式拟合
在MATLAB中，函数polyfit()采用最小二乘法对给定的数据进行多项式拟合，得到该多项式的系数。该函数的调用方式为p = polyfit(x, y, n)，采用n次多项式来拟合数据x和y，得到以p为系数的多项式。该函数使得p(x)与y最小均方误差最小。
【例5-11】 某数据的横坐标为x=[0.2 0.3 0.5 0.6 0.8 0.9 1.2 1.3 1.5 1.8]，纵坐标为y=[1 2 3 5 6 7 6 5 4 1]，对该数据进行多项式拟合，代码如下：

clear all;
x=[0.2 0.3 0.5 0.6 0.8 0.9 1.2 1.3 1.5 1.8];
y=[1 2 3 5 6 7 6 5 4 1];
p5=polyfit(x,y,5); %5阶多项式拟合
y5=polyval(p5,x);
p5=vpa(poly2sym(p5),5) %显示5阶多项式
p9=polyfit(x,y,9); %9阶多项式拟合
y9=polyval(p9,x);
figure; %画图显示
plot(x,y,'bo');
hold on;
plot(x,y5,'r:');
plot(x,y9,'g--');
legend('原始数据','5阶多项式拟合','9阶多项式拟合');
xlabel('x');
ylabel('y');

运行程序后，得到的5阶多项式如下：

p5 =
- 10.041*x^5 + 58.244*x^4 - 124.54*x^3 + 110.79*x^2 - 31.838*x + 4.0393

运行程序后，得到的输出结果如图5.1所示。由图5.1可以看出，使用5次多项式拟合时，得到的结果比较差。当采用9次多项式拟合时，得到的结果与原始数据符合的比较好。当使用函数polyfit()进行拟合时，多项式的阶次最大不超过length(x)-1。

[lili123]

5.1.7 曲线拟合图形用户接口
为了方便用户的使用，在MATLAB中提供了曲线拟合的图形用户接口。它位于MATLAB图形窗口Tools菜单下的Basic Fitting菜单中。在使用该工具时，首先将需要拟合的数据采用函数plot()画图，其MATLAB代码如下：

>> clear all;
x=[0.2 0.3 0.5 0.6 0.8 0.9 1.2 1.3 1.5 1.8];
y=[1 2 3 5 6 7 6 5 4 1];
figure;
plot(x,y,'bo');

该程序运行后，得到Figure窗口，如图5.2所示。然后选择Tools | Basic Fitting命令，弹出Basic Fitting对话框。单击该窗口右下角的 按钮，将会全部展开Basic Fitting对话框，如图5.3所示。

[lili123]

在Basic Fitting对话框的Plot fits选项区域中，勾选5th degree polynomial复选框；在Numerical results选项区域中，会自动列出曲线拟合的多项式系数和残留误差，如图5.4所示。同时，在Figure窗口中会把拟合曲线绘制出来，如图5.5所示。