MATLAB编程与应用系列-第3章 矩阵运算（4）

矩阵的分解是矩阵相关运算中的重要内容，MATLAB提供了用于矩阵分解运算的多种函数。本节将集中介绍MATLAB所提供的矩阵分解运算函数的功能及使用。

矩阵的三角分解又称高斯消去法分解，它的目的是将一个矩阵分解成一个下三角矩阵L和一个上三角矩阵U的乘积，即A=LU。MATLAB提供了专门的函数lu来计算矩阵的LU分解。该函数的调用格式如下：

其中，返回矩阵U为上三角阵，矩阵L为下三角阵或其变换形式，且满足LU=X。返回矩阵P为单位矩阵的行变换矩阵，满足LU=PX。

奇异值分解在矩阵分析中占有极其重要的作用。MATLAB提供了用于矩阵奇异值分解的函数svd，该函数是利用LINPACK程序库中的ZSVDC编制而成的。在计算的过程中假如经过75步QR分解仍得不到一个奇异值，那么系统会给出“不收敛”的提示。奇异值分解函数svd的几种调用格式如下：

其中，命令①返回向量s包含矩阵X分解所得到的全部奇异值向量。命令② 返回一个与X同大小的对角矩阵S和两个酉矩阵U与V，且满足= U S V'。命令③ 得到一个“有效大小”的分解，如果m×n维矩阵X中m>n则只计算出矩阵U的前n列，矩阵S的大小为n×n。

MATLAB提供了eig函数来对矩阵进行特征值分解，该函数的几种调用格式如下：

其中，①计算矩阵A的特征值d，返回结果以向量形式存放。②计算方阵A和B的广义特征值d，返回结果以向量形式存放。③计算矩阵A的特征值对角阵D和特征向量阵V，使AV=VD成立。④计算矩阵A的特征值对角阵D和特征向量阵V，使AV=VD成立。当矩阵A中有与截断误差数量级相差不远的值时，该指令可能更精确。'nobalance'起误差调节作用。⑤计算矩阵A和B的广义特征值向量阵V和广义特征值阵D，满足AV=BVD。最后一条命令⑥由flag指定算法计算矩阵A和B的特征值D和特征向量V。其中，flag的可能值为：'chol' 和'qz' 。当flag值为'chol'时表示对B使用Cholesky分解算法，其中A为对称Hermitian矩阵，B为正定阵。当flag值为'qz'时表示使用QZ算法，其中A、B为非对称或非Hermitian矩阵。

MATLAB提供了chol函数来对矩阵进行Cholesky分解，该函数的调用格式为：

函数调用格式①如果X为n阶对称正定矩阵，则存在一个实的非奇异上三角阵R，满足R'*R = X；若X非正定，则产生错误信息。②不产生任何错误信息，若X为正定阵，则p=0，R与上相同；若X非正定，则p为正整数，R是有序的上三角阵。

正交矩阵是指矩阵的列向量相互正交，且各个列向量的长度相等。QR分解就是将矩阵A分解成一个正交矩阵与一个上三角矩阵的乘积。MATLAB提供了用于矩阵QR分解的函数，表3.7中介绍用于矩阵QR分解的函数调用格式和功能。

表3.7矩阵QR分解

Schur分解将使用schur函数，该函数的调用格式为：

命令行①-③返回正交矩阵U和schur矩阵T，满足A = U T U'。其中，若A有复特征根，则flag='complex'，否则flag='real'。

即使是实阵，在其特征值中也可能出现复数。实际使用中常需要把这一对对共轭复数特征值转化为一个（2x2）的实数块。函数调用格式为：

MATLAB提供了gsvd函数对矩阵进行广义奇异值分解，其具体调用格式为：

其中，函数调用格式①返回酉矩阵U和V、一个普通方阵X、非负对角矩阵C和S，满足A = U C X'，B = V S X'，C' C + S' S = I (I为单位矩阵)。A和B的列数必须相同，行数可以不同。函数调用格式②和①基本相同，而③则返回广义奇异值sigma值。

MATLAB提供了qz函数对矩阵进行特征值问题的QZ分解，该函数的调用格式为：

其中函数调用格式①中A、B为方阵，返回结果AA和BB为上三角阵，Q、Z为正交矩阵或其列变换形式，V为特征向量阵，且满足Q A Z= AA 和Q B Z = BB。命令行②产生由flag决定的分解结果，flag取值为'complex'表示复数分解（默认）；取值为'real'表示实数分解。

如果矩阵H的第一子对角线下元素都是0，则H为海森伯格(Hessenberg)矩阵。如果矩阵是对称矩阵，则它的海森伯格形式是对角三角阵。MATLAB可以通过相似变换将矩阵变换成这种形式，具体调用格式为：

本章将介绍如何使用MATLAB来解决一些基本的数学运算问题，主要包括多项式的相关计算，数据插值，曲线拟合以及数据统计处理等相关的内容。本章的主要内容如下：

在MATLAB中，多项式是以行向量的形式存放的，并且约定多项式以降幂的形式出现，如果多项式中缺少某幂次项，则该幂次项的系数为0。例如，多项式 可以表示为：p1=[1 21 20 0]，其中常数项为0。

本节将全面介绍与多项式有关的各种计算，包括多项式的四则运算、导函数运算、求值、求根以及分部展开。

多项式的加减运算并无特别，可以使用向量的加减运算实现。多项式的乘除运算比较复杂，为此MATLAB提供了专门的运算函数 conv 和 deconv 。

函数 conv 用于求多项式P1和P2的乘积，它的调用格式如下：

其中，P1、P2是两个多项式系数向量。

函数 deconv 用于对多项式P1和P2作除法运算，它的调用格式如下：

其中，Q返回多项式P1除以P2的商式，r返回P1除以P2的余式。返回的Q和r仍是多项式系数向量。

可以将除法运算deconv看作是乘法运算conv的逆运算，即有P1=conv(P2,Q)+r。

下面通过示例介绍多项式的表示和多项式的四则运算。

使用函数poly2str显示多项式p1和p2相乘后生成的新多项式，如下所示：

>>poly2str(y,'x') %以比较习惯的方式显示多项式

ans =

2 x^6 + 15 x^5 - 5 x^4 + 24 x^3 - 20 x^2 + 10 x - 30

MATLAB提供了polyder函数，用于求多项式的导函数。该函数的格式如下：

其中，参数P和Q是多项式的系数向量，返回结果p和q也是多项式的系数向量。

MATLAB提供了两种求多项式值的函数：polyval与polyvalm，它们的输入参数均为多项式系数向量P和自变量x，但是两者是有很大区别的，前者是按数组运算规则对多项式求值，而后者是按矩阵运算规则对多项式求值。具体的调用格式如下所示。

>>p3=[2 6 8 0 5 9 4] %生成多项式系数

>>A=rand(3) %生成随机矩阵

①使用函数polyval按数组运算规则求A中的每个元素对于多项式p3的值。在命令窗口中输入如下内容：

>>Y=polyval(p3,A)

运算结果如下：

Y =

7.2917 15.2885 5.4763

5.8986 7.2672 15.8387

5.9409 11.3612 8.3376

②使用函数polyvalm按矩阵运算规则求以方阵A为自变量的多项式p3的值。在命令窗口中输入如下内容：

>>Y1=polyvalm(p3,A)

运算结果如下：

Y1 =

13.6694 21.1448 16.7431

8.7641 22.5846 21.5403

8.4161 19.5396 22.2629

③如果函数polyval和polyvalm的第二个参数为数值，仍然可以按照数组和矩阵的运算规则计算求多项式在该参数下的结果。在命令窗口中输入如下内容：

>>A=3

分别使用函数polyval和polyvalm计算多项式的值，具体 *** 作及返回结如下：

>>Y=polyval(p3,A)

Y =

3640

>>Y1=polyvalm(p3,A)

Y1 =

3640

④如果函数ployval和ployvalm的第二个参数为一向量，前者按照数组运算规则仍然可以计算求多项式在该参数下的结果，但是后者按矩阵运算规则计算则会提示错误信息。在命令窗口中输入如下内容：

>>p4=[4 8 0 0 0 3 6]

分别使用函数polyval和polyvalm计算多项式的值，具体 *** 作及返回结果如下：

>>Y=polyval(p3,p4) %第二个参数为向量

Y =

16504 754060 4 4 4 3640 150574

>>Y1=polyvalm(p3,p4) %第二个参数必须为方阵或数值

??? Error using ==>polyvalm

Matrix must be square.

n次多项式具有n个根，这些根可能是实根，也可能含有若干对共轭复根。MATLAB提供了roots函数用于求多项式的全部根，该函数的调用格式为：

其中，P为多项式的系数向量，返回向量x为多项式的根，即x(1),x(2),…,x(n)分别代表多项式的n个根。

另外，如果已知多项式的全部根，MATLAB还提供了函数poly用来建立该多项式，该函数的调用格式为：

其中，x为多项式的根，返回向量P为多项式的系数向量。

对于一个方阵s，可以用函数poly来计算矩阵的特征多项式的系数。特征多项式的根即为特征值，可以用roots函数来计算。

MATLAB提供函数 residue 可以实现将分式表达式进行多项式的部分分式展开。

对于 ,函数的调用格式如下：

其中，b和a分别是分子和分母多项式系数行向量；返回值r是[r1 r2 …rn]留数行向量，p为[p1 p2 …pn]极点行向量，k为直项行向量。下面通过示例来讲述该函数的使用。

多项式的微分MATLAB提供了函数 polyder 来实现，前面介绍多项式的导函数时已经介绍了该函数的具体使用。但是对于多项式的积分运算MATLAB没有提供专门的函数，但可以用 [p./length(p):-1:1,k] 的方法来完成积分，其中k为常数。下面通过示例讲解如何进行多项式的积分运算。

插值运算是根据数据点的规律，首先找到一个多项式连接这些已知的数据点，然后根据该多项式计算出要得到的与已知数据点相邻的点对应的数值。数据的插值运算在信号和图象处理等领域使用比较广泛。MATLAB提供了专用的函数来处理数据的插值问题，下面将详细的介绍使用这些插值函数的方法。

一维插值是指对一个自变量的插值，实现一维数据插值的函数是interp1，该函数的调用格式为：

interp1函数的功能是根据X，Y的值，计算出在X1处的值，并返回给Y1。其中，X和Y是两个等长的已知向量，分别描述采样点和样本值；X1是一个向量或标量，描述欲插值的点，返回值Y1是与X1等长的插值结果；method是插值函数的类型，允许的取值有“linear”（线性插值）、“nearest”（用最接近的相邻点插值）、“cubic”（三次插值）和“spline”（三次样条插值），linear为默认值。

除此之外，MATLAB还提供了一个专门的用于3次样条插值的函数spline，功能与函数 Y1=interp1(X,Y,X1,‘spline’) 完全相同，使用方法也类似。该函数的调用格式如下：

除前面介绍的一维数据的插值，MATLAB还提供用于解决二维插值问题的函数 interp2 ，该函数的调用格式为：

其中，X和Y是两个向量，分别描述两个参数的采样点，Z是与参数采样点对应的函数值，X1，Y1是两个向量或标量，描述欲插值的点。返回值Z1是根据相应的插值方法得到的插值结果。method的取值与一维插值函数相同。X，Y，Z也可以是矩阵形式。

多项式曲线拟合是用一个多项式来逼近一组给定的数据，拟合的准则是最小二乘法，即找出使 的 .

在MATLAB中，用 polyfit 函数来求得最小二乘拟合多项式的系数，计算得到多项式后可以用 polyval 函数计算所给出点的近似值。polyfit函数的调用格式为：

polyfit 函数根据采样点X和采样点函数值Y，返回一个m次多项式P及供polyval使用的结构数组S，S有三个域：S.R给出QR分解后满足Q·R=V的矩阵R，S.df给出相应χ2量的自由度，S.normr给出拟合残数的2—范数。其中X,Y是两个等长的向量，P是一个长度为m+1的向量，P的元素为多项式系数。

>>Y=polyval(P,x) %根据多项式系数向量计算对应点x处的拟合函数值

即可计算得到拟合多项式在给定点的函数值。

本节介绍数据统计处理方法，包括最大（小）值运算、求和（积）运算、平均值（中值）运算、累加（乘）运算、标准方差、相关系数以及排序运算。

MATLAB提供的求数据序列的最大值和最小值的函数分别为max和min，两个函数的调用格式和 *** 作过程类似，可以分别用来求向量或矩阵的最大值和最小值。

（1）求向量的最大值和最小值

求向量的最大值和最小值的函数调用格式见表5.1。

表5.1 求向量最大值、最小值函数

（2）求矩阵的最大值和最小值

求矩阵的最大值和最小值的函数调用格式见表5.2。

表5.2 求矩阵的最大值、最小值函数

（3）两个向量或矩阵对应元素的比较

函数max和min还能对两个同型的向量或矩阵进行比较，函数调用格式见表5.3。

表5.3 最大值、最小值函数

MATLAB提供的数据序列求和与求积的函数分别是sum和prod，这两个函数的使用方法类似，分别可以用来对向量和矩阵求和与求积。函数调用格式及功能见表5.4。

表5.4求和与求积函数

在命令窗口中输入：

>>prod(B) %返回各列元素的积

计算得到的各列元素的积如下：

ans =

0.0648 0.0057 0.1780 0.0487

③可以采用下列的方式返回矩阵B各列元素的和与矩阵B各列元素的乘积，具体输入内容和计算返回结果如下：

>>sum(B,1) %返回各列元素的和

ans =

2.2741 2.1284 2.6735 2.2420

>>prod(B,1) %返回各列元素的积

ans =

0.0648 0.0057 0.1780 0.0487

④返回矩阵B各行所有元素的和与矩阵B各行所有元素的积，在命令窗口中输入：

>>sum(B,2) %返回各行元素的和

>>prod(B,2) %返回各行元素的积

MATLAB提供了求数据序列平均值的函数mean与数据序列中值的函数median，函数调用格式及功能见表5.5。

表5.5求平均值与中值函数

④求矩阵A的各行的算术平均值与中值，在命令窗口中输入：

>>mean(A,2) %计算得到矩阵A各行的算术平均值

>>median(A,2) %计算得到矩阵A各行的中值

在MATLAB中，使用cumsum和cumprod函数能方便地求得向量和矩阵元素的累加和与累乘积向量，函数调用格式及功能见表5.6。

表5.6累加和与累加积函数

在MATLAB中，提供了计算数据序列的标准方差的函数std。该函数对于向量X返回一个标准方差；对于矩阵A返回一个行向量，它的各个元素便是矩阵A各列或各行的标准方差。调用格式为：

Y=std(A,flag,dim)

其中，dim可以取1或2。当dim=1时，求各列元素的标准方差；当dim=2时，则求各行元素的标准方差。flag可以取0或1，当flag=0时，置前因子为 ；否则置前因子为 。缺省flag=0和dim=1。

MATLAB提供了corrcoef函数，可以求出数据的相关系数矩阵。调用格式为：

corrcoef函数返回从矩阵X形成的一个相关系数矩阵。此相关系数矩阵的大小与矩阵X一样。它把矩阵X的每列作为一个变量，然后求它们的相关系数。其中X，Y是向量，与corrcoef([X,Y])的作用一样。

>>corrcoef(A) % 求解矩阵A形成的一个相关系数矩阵

返回相关系数矩阵如下：

ans =

1.0000 -0.2608 0.5478 -0.7232

-0.2608 1.0000 -0.9397 0.2996

0.5478 -0.9397 1.0000 -0.3984

-0.7232 0.2996 -0.3984 1.0000

②可以求向量B形成的一个相关系数矩阵。在命令窗口中输入以下内容：

>>corrcoef(B) %求取向量B形成的一个相关系数矩阵

返回相关系数矩阵如下：

ans =

1

MATLAB提供了sort函数来实现排序功能，调用格式如下：

函数返回一个对X中的元素按升序排列的新向量，Y是排序后的矩阵，而I记录Y中的元素在A中位置。其中，dim指明对矩阵A的列还是行进行排序。若dim=1，则按列排；若dim=2，则按行排。

Y =

0.5226 0.1730 0.0118 0.1991

0.7948 0.2523 0.1365 0.2987

0.8801 0.2714 0.7373 0.6614

0.9568 0.9797 0.8757 0.8939

I =

3 1 4 2

1 4 3 3

本章重点介绍了MATLAB提供的基本的数学运算功能，主要包括多项式的相关运算、数据的插值与拟合运算、数据统计处理运算以及傅立叶变换等内容。

多项式部分通过实例重点介绍了多项式运算相关的内容，主要包括多项式的表示方法、多项式的四则运算、多项式的求导运算、多项式的求值与求根运算、多项式的展开以及多项式的积分运算等内容。掌握多项式运算内容是深入学习MATLAB其他内容的基础。

数据插值部分和拟合部分通过实例介绍了MATLAB提供的数据处理时经常使用的数据插值和拟合运算函数，其中插值部分包括一维和二维数据的插值运算。曲线拟合运算时要正确的选择所要拟合的多项式的阶，并不是拟合多项式的阶越高精度越好，一般拟合多项式的阶不超过5阶。

数据统计处理部分通过实例详细介绍一些常用的数据统计处理方法，主要包括数据的最大值与最小值运算、求和与求积运算、平均值与中值运算、累加和与累乘积运算、标准方差、相关系数以及排序等运算。

结构数组的基本组成是结构(Structure)，每个结构包含多个域(Fields)。数据只能存储在域中，与元胞数组一样，结构数组也可以存储任何类型的数据，但是从一定意义上讲，结构数组组织数据的能力比元胞数组要强一些。本节讲述结构数组的一些基本 *** 作。

结构数组可以通过直接对域进行赋值来创建，也可以通过MATLAB提供的函数struct来创建。下面通过示例来讲解结构数组的创建方法。

结构数组可以通过结构数组标识直接获取整个结构数组，或通过域名来获取结构数组中某个域中存储的具体内容，或通过MATLAB提供的函数getfield来获取某个域中存储的具体内容。要获取域中的具体内容，需要预先知道结构数组的各域的名称，为此MATLAB还提供了查询结构数组中各域名的函数fieldnames方便使用。

MATLAB提供了修改结构数组中域存储内容的函数setfield。下面通过示例演示该函数的具体使用。

本章重点介绍了MATLAB提供的3种重要的数据类型的创建及相关 *** 作，这3种数据类型分别是：字符数组、元胞数组和结构数组。理解和掌握这3种数据类型是进一步学习MATLAB的基础。

字符数组部分使用大量的实例重点介绍了简单和复杂字符数组的创建方法、字符数组的比较 *** 作、字符串的查找和替换 *** 作以及字符数组与其他数据类型之间的相互转化 *** 作等内容。大量的实例使得读者可以轻松的掌握字符数组的相关 *** 作。

元胞数组部分首先介绍了元胞数组中元胞和元胞中所包含的内容为两个不同的范畴，并且对应的是两种截然不同的 *** 作，分别为元胞外标识和元胞内编址。在此基础上通过实例介绍了元胞数组的多种创建方法、元胞数组内容的访问以及元胞数组内容的显示等内容。

结构数组部分首先介绍了结构数组的定义，然后通过实例重点介绍了结构数组的创建以及对元胞数组中元素进行访问和修改等 *** 作相关的内容。与元胞数组一样，结构数组也可以存储任何类型的数据，但是从一定意义上讲，结构数组组织数据的能力比元胞数组要强一些。

---