Ch51基本应用领域.docx
《Ch51基本应用领域.docx》由会员分享，可在线阅读，更多相关《Ch51基本应用领域.docx（29页珍藏版）》请在冰点文库上搜索。
Ch51基本应用领域
第五章MATLAB基本应用领域
MATLAB的应用领域非常广泛，从最基本的线性代数、泛函分析，到应用广泛的信号处理、控制系统、通信系统，直到最新技术领域神经网络、小波理论等。
本章主要介绍MATLAB在基本应用领域：
线性代数、多项式与内插、数据分析与统计、泛函分析、常微分方程求解中的应用。
§5.1线性代数
5.1.1矩阵运算—参Ch1-基本操作
Ex.1设A=[13；24]判断A中元素是否大于2。
>>A=[13;24];
>>A>2
ans=
01
11
Ex.2找出A中小于0的元素的位置。
>>A=[1-23;45-4;5-67]
A=
1-23
45-4
5-67
>>B=find（A<0）
B=
4
6
8
5.1.2矢量范数和矩阵范数
1．矢量x的p-范数定义为
特别地，有
向量的1-范数：
向量的2-范数：
向量的-范数：
一般用norm（x,p）函数实现，p缺省时为p=2.如
>>v=[20-1];
>>n=[norm（v,1）,norm（v）,norm（v,inf）]
n=3.00002.23612.0000
2．矩阵A的p-范数定义为
特别地，有
矩阵的1-范数：
向量的2-范数：
向量的-范数：
一般用norm（A,p）函数实现，p缺省时为p=2.如
>>A=fix（10*rand（3,2））
A=
20
66
83
>>N=[norm（A,1）,norm（A）,norm（A,inf）]
N=16.000011.889512.0000
5.1.3线性代数方程求解
1．当A为非奇异时，解为
例1设矩阵
>>A=magic（3）
A=
816
357
492
举4种情形求解：
（1）Ax=b,b=[1:
3]’,x为未知列向量
>>x=A\b%给出方程组的解
x=
0.0500
0.3000
0.0500
>>A*x%验证解的正确性
ans=
1.0000
2.0000
3.0000
（2）xA=b,b=[1:
3]，x为未知行向量
>>b=[1:
3]%验证解的正确性
b=123
>>x=b/A%给出方程组的解
x=-0.03330.4667-0.0333
>>x*A
ans=1.00002.00003.0000
（3）AX=B
>>B=Pascal（3）
B=
111
123
136
>>X=A\B
X=
0.06670.05000.0972
0.06670.30000.6389
0.06670.0500-0.0694
>>A*X%验证解的正确性
ans=
1.00001.00001.0000
1.00002.00003.0000
1.00003.00006.0000
（4）XA=B%矩阵方程
>>X=B/A%给出方程组的解
X=
0.06670.06670.0667
-0.03330.4667-0.0333
-0.15281.0556-0.2361
>>X*A%验证解的正确性
ans=
1.00001.00001.0000
1.00002.00003.0000
1.00003.00006.0000
2．当A为mn（m>n）矩阵时,AX=B中方程个数多于未知量个数，应采用最小二乘法来求解。
例如，对一维测量数据，根据数据图形趋势（如下图中虚线上红圆点所示）选取延迟指数函数
来拟合这组数据（这将是6个方程两个未知量的线性方程组（也称超定方程组？
））：
>>t=[0.3.81.11.62.3]';
>>y=[.82.72.63.60.55.50]';
>>A=[ones（size（t）exp（-t））;
>>c=A\y
c=
0.4760
0.3413
>>T=[0:
.1:
2.5]';
Y=[ones（size（T））exp（-T）]*c;
plot（t,y,'--',T,Y,'-',t,y,'o'）
title（'最小二乘法曲线拟合'）
xlabel（'\itt'）,ylabel（'\ity'）
以上得到拟合函数
，其曲线与数据关系图为
3．当方程组有无穷多解情况（参Ch2,§2.9符号方程求解）
练习：
求解
.
解（按第2种情形，以最小二乘拟合求得：
）
>>A=[11;12;14;18;116];
>>b=[3.03,5.05,9.09,16.98,33.04]';
>>x=A\b
x=
1.0467
1.9986
5.1.4矩阵求逆
求解AX=B的方程组时，用X=A\B比用X=inv（A）*B所用内存、时间、误差小等特点更优。
pinv（A）用于计算非方阵的伪逆，例如：
>>c=[94;28;67]
c=
94
28
67
>>x=pinv（c）
x=
0.1159-0.07290.0171
-0.05340.11520.0418
则可使x*c为单位阵，即
>>q=x*c
q=
1.00000.0000
1.01.0000
但应注意c*x并非单位阵，即
>>p=c*x
p=
0.8293-0.19580.3213
-0.19580.77540.3685
0.32130.36850.3952
5.1.5LU、QR分解
若有A=LU，则可由[L，U]=lu（A）求得。
如：
>>A=[123;456;426]
A=
123
456
426
>>[L,U]=lu（A）
L=
0.2500-0.25001.0000
1.000000
1.00001.00000
U=
4.00005.00006.0000
0-3.00000
001.5000
通过LU分解，可很容易得到：
det（A）=det（L）*det（U）;
inv（A）=inv（U）*inv（L）.
求线性方程Ax=b时，可得到
x=U\（L\b）
这种方法运算速度更快。
正交矩阵或具有正交列的矩阵，其所有列的长度为1，且与其它列正交。
即如果Q为正交矩阵，则
Q’Q=I.
通过正交或QR分解，可将任意二维矩阵分解成一个正交阵和一个上三角阵的乘积，即A=QR。
如
>>A=[94;28;67]
A=
94
28
67
>>[Q,R]=qr（A）
Q=
-0.81820.3999-0.4131
-0.1818-0.8616-0.4739
-0.5455-0.31260.7777
R=
-11.0000-8.5455
0-7.4817
00
5.1.6矩阵求幂和矩阵指数
矩阵求幂如A2、B3可很容易求出：
A^2,B^3.元素对元素的求幂，可输入：
A.^2,B.^3.
函数dqrtm（A）可求出
。
函数expm（A）可计算出矩阵A的指数，即eA.这些函数对求解微分方程是很有用的。
例如，要求解微分方程：
其解为
因此可输入：
>>x0=[1;1;1];
>>x=[];
>>X=[];
>>fort=0:
.01:
1
X=[Xexpm（t*A）*x0];
end
plot3（X（1,:
）,X（2,:
）,X（3,:
）,'-o'）
gridon
5.1.7特征值
>>A=[0-6-1;62-16;-520-10];
>>lambda=eig（A）
lambda=
3.0710
-2.4645+17.6008i
-2.4645-17.6008i
>>[V,D]=eig（A）
V=
-0.83260.2003-0.1394i0.2003+0.1394i
-0.3553-0.2110-0.6447i-0.2110+0.6447i
-0.4248-0.6930-0.6930
D=
-3.071000
0-2.4645+17.6008i0
00-2.4645-17.6008i
5.1.8奇异值分解——P164
§5.2多项式与内插
多项式函数如roots,poly,polyval,polyvalm,residue,polyfit,polyder,conv,deconv等，它们使对多项式的操作变得便捷迅速。
5.2.1多项式表示
用行矢量表示，其元素按幂指数降序排列，例如
或表示成
p=[10–2–5];
5.2.2多项式的根
为求多项式的根，即p（x）=0的解，可利用roots函数：
>>r=roots（p）
r=
2.0946
-1.0473+1.1359i
-1.0473-1.1359i
利用poly函数可从多项式的根中恢复出多项式：
>>p2=poly（r）
p2=
1.00000-2.0000-5.0000
p=[10–2–5];
5.2.3特征多项式
poly函数还可以用于计算矩阵的特征多项式的系数：
>>A=[12-1;345;-192];
>>poly（A）
ans=1.0000-7.0000-38.000090.0000
5.2.4多项式计算
polyval函数可计算出多项式在指定点处的值，例如：
>>y1=polyval（p,4）
y1=51
5.2.5卷积和去卷积
多项式的乘和除就相当于卷积和去卷积操作，这可由函数conv和deconv实现。
例如：
设
>>a=[133];b=[456];
>>c=deconv（c,b）
c=
413282718
通过多项式除法，可以得到：
这可由deconv函数求出：
>>[q,r]=deconv（c,b）
q=
123
r=
00000
5.2.6多项式求导
polyder函数可用于求出单个多项式的导数，也可用于求两个多项式之积或之比的导数，这可由下列示例说明：
>>p=[10-2-5];
>>q=polyder（p）%求p的导数
q=30-2
>>a=[135];b=[246];
>>c=polyder（a,b）%求积a*b的导数
c=
8305638
>>[q,r]=polyder（a,b）%求商a/b的导数
q=
-2-8-2
r=
416404836
最后求得的q（s）/r（s）为a（s）/b（s）的导数。
5.2.7多项式曲线拟合
polyfit函数可在最小二乘意义下找出一多项式来拟合给定的一组数据。
其调用格式为
p=polyfit（x,y,n）
其中x,y为给定的数据,n为多项式的阶数。
例如，用三阶多项式来拟合下列数据：
>>x=[12345];y=[5.543.1128290.7498.4];
p=polyfit（x,y,3）;
x2=1:
.1:
5;
y2=polyval（p,x2）;
figure
（1）
subplot（2,1,1）
plot（x,y,'o',x2,y2）
gridon
title（'多项式曲线拟合'）
5.2.8部分分式展开
residue函数可将有理多项式进行部分分式展开：
其调用格式为
[r,p,k]=residue（b,a）
residue函数还可将部分分式形式的多项式还原成有理多项式，其调用格式为
[b,a]=residue（r,p,k）.
Ex.
>>a=[100-2];b=[1-1];
>>[q,r]=deconv（a,b）
q=
111
r=
000-1
得到
>>[r,p,k]=residue（a,b）
r=
-1
p=
1
k=
111
得到
Ex.
>>a=[1];b=[1-32];
>>[r,p,k]=residue（a,b）
r=
1
-1
p=
2
1
k=
[]
5.2.9一维内插
MATLAB中提供了两类一维内插：
多项式内插和其于FFT的内插。
1．多项式内插
interp1函数可完成一维内插，其调用格式为
yi=interp1（x,y,xi,method）
其中x,y为给定的数据对，xi为要内插的点矢量，method用于指定内插方法，可取：
nearest（最邻近内插）：
将内插点设置成最接近于已有数据点的值；
linear（线性内插）：
连接已有数据点作线性逼近。
这是interp1函数的缺省设置；
spline（三次样条内插）：
利用一系列样条函数获得内插数据点，从而确定已有数据点之间的函数；
cubic（三次曲线内插）：
通过y拟合三次曲线函数，从而确定内插点的值。
所有这四种方法都要求x中的数据为单调，每种方法不要求x为均匀间隔，但如果x已经为均匀间隔，则在method之前加上*，可使执行速度加快。
这四种方法对执行速度、内存要求及得到的平滑度是不同的。
总的说来，按nearst、linear、cubic、spline顺序，内存要求从小到大，执行速度由快到慢，平漫度由差到好。
2．基于FFT的内插
interpft完成基于FFT的一维内插，其调用格式为
y=interpft（x,n）
其中x中包含等间隔取样的周期函数值，n为要求得到的点数。
5.2.10二维内插——P168
§5.3数据分析与统计
MATLAB提供的数据分析与统计函数都是面向列的，即矩阵中的每一列代表一个变量的多项观测值，其列数相应于变量数，行数相应于测量点数。
Max和min函数可求出数据的最大值和最小值,mean和std函数可求出数据的均值和标准差,sum和prod函数可求出数据元素和与数据元素积.例如,对MATLAB内含的城市24小时的车流量数据count.dat,可作分析:
>>loadcount.dat
>>mx=max（count）
mx=
114145257
>>mu=mean（count）
mu=
32.000046.541765.5833
>>sigma=std（count）
sigma=
25.370341.405768.0281
对于有些函数还可给出位置,例如,在求出最小值的同时,可得到最小值所在的位置（行号）;
>>[mx,indx]=min（count）
mx=
797
indx=
22324
5.32.1协方差和协相关系数
cov函数可以求出单个变量的协方差,而corrcoef可求出两个变量之间的相关系数.例如
>>cv=cov（count）
cv=
1.0e+003*
0.64370.98021.6567
0.98021.71442.6908
1.65672.69084.6278
>>cr=corrcoef（count）
cr=
1.00000.93310.9599
0.93311.00000.9553
0.95990.95531.0000
5.3.2数据预处理
MATLAB中遇到超出范围的数据是均用NaN（非数值）表示,而且在任何运算中,只要包含NaN,就将它传递到结果中,因此在对数据进行分析前,应对数据中出现的NaN作剔除处理.例如:
>>a=[1,2,3;5NaN8;742];
>>sum（a）
ans=13NaN13
在矢量x中删除NaN元素,可有下列四种方法:
（1）i=find（~isnan（x））;x=x（i）;
（2）x=x（find（~isnan（x）））;
（3）x=x（~isnan（x））;
（4）x（isnan（x））=[];
在矩阵X中删除NaN所在的行,可输入:
X（any（isnan（x）’）,:
）=[]
经过这种预处理后的数据,可进行各种分析和统计操作.
5.3.3回归和曲线拟合
对给定数据进行拟合,可采用多项式回归,也可采用其它信号形式的回归,其基本原理是最小二乘法,这在MATLAB中显得轻而易举.
例1设通过测量得到一组时间t与变量y的数据:
t=[0.3.81.11.62.3]’;
y=[.5.821.141.251.351.40]’;
分别采用多项式
进行回归,可得到两种不同的结果.MATLAB程序如下:
X1=[ones（size（t））tt.^2];
a=X1\y;%解出y1的系数
X2=[ones（size（t））exp（-t）t.*exp（-t）];
b=X2\y;%解出y2的系数
T=[0:
.1:
2.5]';%为作出拟合曲线,将区间均分,步长为0.1
Y1=[ones（size（T））TT.^2]*a;%计算每个分点的函数值,即曲线上纵坐标
Y2=[ones（size（T））exp（-T）T.*exp（-T）]*b;
figure
（1）
subplot（2,2,1）
plot（T,Y1,'-',t,y,'o'）,gridon
title（'多项式回归'）
subplot（2,2,2）
plot（T,Y2,'-',t,y,'o'）,gridon
title（'指数函数回归'）
例2已知变量y与x1、x2有关，测得一组数据为
>>x1=[.2.5.6.81.01.1]';
>>x2=[.1.3.4.91.11.4]';
>>y=[.17.26.28.23.27.24]';
采用
来拟合，则有
>>a=X\y
a=
0.1018
0.4844
-0.2847
因此数据的拟合模型为
5.3.4滤波——P172
5.3.5傅里叶分析FFT（快速傅里叶变换）——P173
§5.4泛函分析
MATLABA提供了可对函数进行操作的函数，它们称之为泛函，如找出函数在区间上的极小值、求函数零值点、计算函数积分等，这些都属于泛函。
5.4.1函数在MATLAB中的表示
在MATLAB中，数学函数可用M文件表示，例如函数
在MATLAB中可表示成humps.m:
functiony=humps（x）
y=1./（x-0.3）.^2+0.01）+1./（（x-0.9）.^2+0.04）-6;
这样，humps可用作为某些函数的输入变量，从而构成泛函的计算。
5.4.2数学函数的绘图
fplot函数可绘制出指定函数的图形，它可以指定函数自变量和函数值的范围，还可在同一张图上绘制出多个图形。
例如，下面程序段执行后，可得到如下图线。
>>figure
（1）
subplot（2,2,1）,fplot（'humps',[-5,5]）,gridon
subplot（2,2,2）,fplot（'humps',[-55-1025]）,gridon
subplot（2,2,3）,fplot（'5*sin（x）',[-1,1]）,gridon
subplot（2,2,4）,fplot（'[5*sin（x）,humps（x）]',[-1,1]）,gridon
5.4.3函数极小值点和零值点
fmin函数用于求出指定单变量函数在特定区间上的极小值点，fmins函数用于求出多变量函数的极小值点。
Fzero函数可求出指定函数在特定区间上零值点。
例如，要求出函数humps在[0.3,1]区间上的极小值点，可输入：
>>x=fmin（'humps',0.3,1）
Warning:
FMINisobsoleteandhasbeenreplacedbyFMINBND.
FMINnowcallsFMINBNDwhichusesthefollowingsyntax:
[X,FVAL,EXITFLAG,OUTPUT]=FMINBND（FUN,x1,x2,OPTIONS,P1,P2,...）
UseOPTIMSETtodefineoptimizationoptions,ortype
'editfmin'toviewthecodeusedhere.FMINwillbe
removedinthefuture;pleaseuseFMINBNDwiththenewsyntax.
>InC:
\MATLAB6p1\toolbox\matlab\funfun\fmin.matline62
x=
1.6370
>>humps（x）%求出相应极小值
ans=
11.2528
从上图中可以看出，函数humps在[-1,0]上有一零值点，这时可用两种方式来得到：
>>a=fzero（'humps',-0.2）
a=
-0.1316
或
>>a=fzero（'humps',[-1,0]）
a=
-0.1316
5.4.4数值积分
quad和quad8函数可以求出指定函数在指定区间上的积分。
例如，q1=quad（‘humps’,0,1）可求出humps函数在[0,1]上的积分值，当然采用quad8函数可得到相同的结果，只是采用的积分算法有所不同。
1．如quad,quad8:
>>q=quad（'exp（-x）+x.^2',0,1）%注意x2应写成x.^2，即按元素计算函数对应值序列
q=
0.9655
或
>>q=quad8（'exp（-x）+x.^2',0,1）
Warning:
QUAD8isobsolete.QUADLisitsrecommendedreplacement.
>InC:
\MATLAB6p1\toolbox\matlab\funfun\quad8.matline35
q=
0.9655
2.dlquad函数可以计算函数的二重积分，例如要计算
，可输入：
>>q2=dblquad（'f（x,y）',x1,x2,y1,y2）
注：
书上的上面这种写法有问题。
dblquad函数的详细用法，可参看帮助文件：
>>helpdblquad。
从中找到应用的例子，再用来求解：
>>dblquad（'sin（x）+y-1'）,1,2,1,4）
ans=
7.3693
或先用
升级会员 