MATLAB与系统仿真14

1、MATLAB与系统仿真机电工程学院 勾燕洁 二一四年秋第九章 MATLAB在各学科中 的应用 从根本上讲，复变函数的运算是以前所讲所有从根本上讲，复变函数的运算是以前所讲所有实变函数运算的延伸，但由于其自身的一些特殊的性实变函数运算的延伸，但由于其自身的一些特殊的性质而有所不同，特别是当它引进了质而有所不同，特别是当它引进了“留数留数”的概念，的概念，且在引入了数学上的且在引入了数学上的TaylorTaylor级数展开和三大变换级数展开和三大变换（LaplaceLaplace变换、变换、FourierFourier变换和变换和ZetaZeta变换）之后而使变换）之后而使其显得更为重要了。其显得

2、更为重要了。9.1 MATLAB在复变函数中的应用9.1.1 9.1.1 复数和复矩阵的生成复数和复矩阵的生成复数的生成 复数可由z=a+b*i语句赋值生成，也可简写成z=a+bi；另一种生成复数的语句是z=r*exp(i*theta)，也可简写成z=r*exp(theta i)，其中theta为复数辐角的弧度值，r为复数的模。创建复数矩阵（1）C=1+5i 2+6i;3+7i 4+8i（2）C=1 2;3 4+1i*5 6;7 89.1.2 9.1.2 复数的运算复数的运算复数的实部和虚部共轭复数复数的模和辐角real(z)imag(z)conj(z)abs(z)angle(z)例例 求下列

3、复数的实部、虚部、共轭复数、模与辐角：求下列复数的实部、虚部、共轭复数、模与辐角：921113(34 )(25 )(1)(2)(3)(4)43212iiiiiiiiiia=1/(3+2i),1/i-3i/(1-i),(3+4i)*(2-5i)/2i,i9-4*i21+ireal(a)imag(a)conj(a)abs(a)angle(a) 在MATLAB工作空间中，复数的乘除法运算由“*”和“/”实现。复数的乘除法例例x=4*exp(pi/3i)y=3*exp(pi/5i)y1=3*exp(pi/5*i)x/yx/y1 复数的平方根运算由函数sqrt实现。复数的平方根 复数的幂运算的形式为xn

4、，结果返回复数x的n次幂。复数的幂运算sqrt(z)(-i)(1/6)例例 复数的指数和对数运算分别由函数exp和log实现。调用形式如下：复数的指数和对数运算exp(z)log(z)例例log(-i)log(-3+4i)复数方程求根 复数方程求根或实方程的复数根求解也由函数solve实现。例例 求方程求方程x3+8=0所有的根。所有的根。solve(x3+8=0) 复数的三角函数运算参见下表。复数的三角函数运算函数意义函数意义sin正弦asin反正弦cos余弦acos反余弦tan正切atan反正切cot余切acot反余切sec正割asec反正割csc余割acsc反余割sinh双曲正弦asin

5、h反双曲正弦cosh双曲余弦acosh反双曲余弦tanh双曲正切atanh反双曲正切MATLAB中，留数可由函数residue求出：9.1.3 9.1.3 留数留数R,P,K=residue(B,A)( )(1)(2)( )( )( )(1)(2)( )B sRRR nK sA ssPsPsP n 如果没有重根，则向量B和A为分子、分母以s降幂排列的多项式系数，留数返回为向量R、极点的位置在向量P，直接项返回到向量K。 极点的数目n=length(A)-1=length(R)=length(P)。如果length(B)d=74.001,74.005,74.003,74.001,74.000,7

6、3.998,74.006,74.002;mean(d) 对于有些情况下，各样本值还有与其相对应的权，这时只要先把样本值数组与权数组做点乘，再把结果做简单期望运算即可。例例 设随机变量设随机变量X的分布律见下表，求的分布律见下表，求E(X)和和E(3X2+5)的的值。值。x=-2 0 2;pk=0.4 0.3 0.3;sum(x.*pk)x-202pk0.40.30.3z=3*x.2+5;sum(z.*pk)数据比较 在给定的一组数据中，还常要对它们进行最大、最小、中值的查找或对它们排序等操作。 最大值和最小值最大值和最小值 MATLAB提供的求数据序列的最大值和最小值的函数分别为max和min

7、，两个函数的调用格式和操作过程类似。 y,I=max(X)：返回向量X的最大值存入y，最大值的序号存入I，如果X中包含复数元素，则按模取最大值。 y=max(X)：返回向量X的最大值存入y，如果X中包含复数元素，则按模取最大值。 求一个向量X的最大值的函数有两种调用格式，分别是：例 求向量x的最大值。x=-43,72,9,16,23,47;y=max(x) %求向量x中的最大值y,l=max(x) %求向量x中的最大值及其该元素的位置 求矩阵A的最大值的函数有3种调用格式，分别是： 求最小值的函数是min，其用法和max完全相同。 max(A,dim)：dim取1或2。dim取1时，该函数和m

8、ax(A)完全相同；dim取2时，该函数返回一个列向量，其第i个元素是A矩阵的第i行上的最大值。 Y,U=max(A)：返回行向量Y和U，Y向量记录A的每列的最大值，U向量记录每列最大值的行号。 max(A)：返回一个行向量，向量的第i个元素是矩阵A的第i列上的最大值。 平均值和中值平均值和中值 求数据序列平均值的函数是mean，求数据序列中值的函数是median。两个函数的调用格式为： median(A,dim)：当dim为1时，该函数等同于median(A)，返回一个行向量，其第i个元素是A的第i列的中值；当dim为2时，返回一个列向量，其第i个元素是A的第i行的中值。 mean(A,di

9、m)：当dim为1时，该函数等同于mean(A)，返回一个行向量，其第i个元素是A的第i列的算术平均值；当dim为2时，返回一个列向量，其第i个元素是A的第i行的算术平均值。 median(X)：返回向量X的中值。 mean(X)：返回向量X的算术平均值。nanmin求随机变量的忽略NAN的最小值元素nanmedian求随机变量的忽略NAN元素的中值mad求随机变量的绝对差分平均值(与数组平均值偏差的平均值)sort对随机变量由小到大排序sortrows对随机矩阵按首行进行行排序range求随机变量的极差，即最大值与最小值的差nanmax求随机变量的忽略NAN的最大值元素 MATLAB还提供了

10、其他几个求平均数的函数：nanmeangeomeanharmmeantrimmean求随机变量的算术平均求随机变量的几何平均求随机变量的和谐平均求随机变量的调整平均累积和累和 求向量或矩阵的元素累积或累和运算是比较常用的两种运算，在MATLAB中可由以下函数实现：sum求向量或矩阵元素的累和nansum忽略NAN求向量或矩阵的元素累和cumsum求此元素位置以前的元素和cumtrapz梯形累和函数cumprod求当前元素与所有前面位置的元素的积方差和标准差 对于一组采样数据，有时只从其累和或均值中还不能判断它们的质量好坏，为了表征随机变量与其均值的偏离程度，我们引入了方差或标准方差（标准差）。

11、方差表示为2( )( )( ) D xVar xExE x标准差表示为( )( )xD x 对于样本来说，样本方差为2211()1niisxxn标准差为211()1niisxxn 其中dim取1或2。当dim=1时，求各列元素的标准方差；当dim=2时，则求各行元素的标准方差。flag取0或1，当flag=0时，置前因子为1/（n1），当flag=1时，按置前因子为1/n计算标准方差。缺省flag=0，dim=1。std(A,flag,dim) 在MATLAB中，提供了计算数据序列的方差和标准方差函数。其一般调用格式为： 对于向量X，std(X)返回一个标准方差。对于矩阵A，std(A)返回一

12、个行向量，它的各个元素便是矩阵A各列或各行的标准方差。 var(x,1) 若X为向量，var(X)返回向量的样本方差值。若X为矩阵，则返回一个行向量，它的各个元素便是矩阵X各列的方差。有参数1时，返回X的简单方差（即置前因子为1/n）。例例 随机地取随机地取8只活塞环，测得它们的直径为（以只活塞环，测得它们的直径为（以mm计）：计）：74.001，74.005，74.003，74.001，74.000，73.998，74.006，74.002。对样本值对样本值d，求其方差值、样本方差值、，求其方差值、样本方差值、标准方差、样本标准差的值。标准方差、样本标准差的值。d=74.001,74.005

13、,74.003,74.001,74.000,73.998,74.006,74.002;var(d,1)var(d)std(d,1)std(d)偏斜度和峰度 为描述随机变量分布的形状与对称形式或正态分布型的偏离程度，引入了特征量的偏斜度和峰度。 偏斜度偏斜度偏斜度定义为21( )() ( )xE xvED x 此函数表征分布形状偏斜对称的程度，若v10则可以认为分布是对称的，若v10则称为右偏态，此时位于均值右边的值比位于左边的值多一些；反之，则称为左偏态即位于均值左边的值比位于右边的值多一些。 此函数在MATLAB中由功能函数skewness实现。其调用形式为： 峰度峰度峰度的定义为42( )

14、() ( )xE xvED xskewness(X) 若X为向量则返回向量的偏斜度，若X为矩阵，则返回矩阵列向量的偏斜度行向量。 若v20，表示分布有沉重的“尾巴”，即数据中含有较多偏离均值的数据，对于正态分布，v2=0，故v2的值也可以看成是数据偏离正态分布的尺度。 在MATLAB中此函数由kurtosis实现。kurtosis(X) 若X为向量则返回向量的峰度，若为矩阵则返回矩阵列向量的峰度行向量。w=75.0,64.0,47.4,66.9,62.2,62.2,58.7,63.5,66.6,64.0,57.0,61.0,56.9,50.0,72.0;M=mean(w)D=std(w)S=s

15、kewness(w)K=kurtosis(w)例例 有有15名学生的体重（单位为名学生的体重（单位为kg）为）为75.0，64.0，47.4，66.9，62.2，62.2，58.7，63.5，66.6，64.0，57.0，61.0，56.9，50.0，72.0，计算此，计算此15名学生体重的均值、标准方名学生体重的均值、标准方差、偏斜度和峰度。差、偏斜度和峰度。协方差和相关系数 对于二维随机变量（X，Y），除了讨论X与Y的数学期望和方差以外，还需讨论描述X与Y之间相互关系的数学特征，我们引入了协方差和相关系数。协方差的定义为( , )( )( )cov x yExE xyE y相关系数的定义为

16、cov( , )( , )( )( )x ycof x yD xD y 在MATLAB中这两个量分别由功能函数cov和corrcoef来实现。corrcoef(X,Y)：返回列向量X,Y的相关系数。corrcoef(X)：返回从矩阵X形成的一个相关系数矩阵。此相关系数矩阵的大小与矩阵X一样。它把矩阵X的每列作为一个变量，然后求它们的相关系数。cov(X)：X为向量时，函数返回此向量的协方差，X为矩阵时，函数返回此矩阵的协方差矩阵，其中协方差矩阵的对角元素是X矩阵的列向量的方差值。cov(X，Y)：函数返回向量X，Y的协方差矩阵，且X和Y的维数必须相同。cov(X)：返回向量X的样本协方差（矩阵

17、），即置前因子为1/(n-1)，也可写成cov(X，0) 。cov(X，1)：返回向量的协方差（矩阵），即置前因子为1/n。cov(X，Y)与cov（X，Y，1）的区别同上。例例 协方差矩阵函数和相关系数函数应用示例。协方差矩阵函数和相关系数函数应用示例。a=1,2,1,2,2,1;var(a)cov(a) 对向量对向量d=rand(2,3)cov1=cov(d) 对矩阵对矩阵例例 生成满足正态分布的生成满足正态分布的100005随机矩阵，然后求各列随机矩阵，然后求各列元素的均值和标准方差，再求这元素的均值和标准方差，再求这5列随机数据的相关系数列随机数据的相关系数矩阵。矩阵。X=randn(

18、10000,5);M=mean(X)D=std(X)R=corrcoef(X)协方差矩阵 对于二维随机变量（X，Y），定义协方差矩阵其中2111221222( ) ( )( )( )( )( ) CExE xCExE xyE yCEyE yxE xCEyE y11122122cov(, )CCX YCC 在MATLAB中实现协方差矩阵的函数与实现协方差的函数cov一样。例如例如 X=rand(3,3)cov(X)corrcoef(X)9.2.2 9.2.2 常用的统计分布量常用的统计分布量 常用的统计分布量包括各种统计分布的期望和方差、各种分布的概率密度函数、概率累积函数和分值点函数及常用的随

19、机数等等。函数名数学期望方差调用格式二项分布E=NPD=NP(1-P)M,V=binostat(N,P)泊松分布E=lD=lM,V=poisstat( l)几何分布M,V=geostat(P)超几何分布MN,V=hygestat(N,M,K)b 分布M,V=betastat(A,B)1PEP21PDP2()()(1)KM NMNKDNN2(1)()ABDABABAEABNMKEN期望和方差函数名数学期望方差调用格式均匀分布M,V=unifstat(A,B)指数分布M,V=expstat( m )正态分布M,V=normstat(m,s)分布M,V=gamstat(A,B)t-分布M,V=tst

20、at(N)2BAM12)(2ABVM2VM2VABM 2ABV 1, 0NM2,2NNNV例例 求参数为求参数为0.12和和0.34的的b b分布的期望和方差。分布的期望和方差。 m,v=betastat(0.12,0.34)例例 求参数为求参数为6的泊松分布的期望和方差。的泊松分布的期望和方差。 m,v=poisstat(6) 由此可见泊松分布参数的值与它的期望和方差是相同的。 m,v=binostat(20,0.2)例例 按规定，某型号的电子元件的使用寿命超过按规定，某型号的电子元件的使用寿命超过1500小小时为一级品，已知一样品时为一级品，已知一样品20只，一级品率为只，一级品率为0.2

21、，问样，问样品中一级品元件的期望和方差为多少？品中一级品元件的期望和方差为多少？ 由上面的例子可以看出，求期望和方差的函数都是以stat结尾的，几乎所有的分布都可以找到对应的函数，输入输出参数直观，调用方便。 样品数目20为试验次数，分布概率为0.2。选用二次分布函数binostat。概率密度函数 通用函数介绍通用函数介绍 在MATLAB中求概率密度的函数中有一个通用的函数pdf，通过此函数可以求任何输入分布的概率密度函数。Y=pdf(name,X,A,B,C) 返回以name为分布，在X值的参数为A，B，C的概率密度。name的值可以如下表所示), 1 , 0(1NxNy), 1 , 0()

22、1 (NxppCyxNxxN), 2 , 1 , 0()1 (1xppCyrxxxrl), 1 , 0(! nxexyxll函数名调用格式数学意义离散均匀分布Y=unipdf(X,N)二项分布Y=binopdf(X,N,P)负二项分布Y=Nbinpdf(X,R,P)泊松分布Y=Poisspdf(X, )几何分布Y=Geopdf(X,P), 3 , 2 , 1()1 (1xppyx概率密度函数表 专用函数列表专用函数列表)0(1/xeyxbb,)(21222/)(xeyx)0()2/(212/12/2/xexnyxnn)0()(1/1xexabybxaab) 10()1 (),(111xxxba

23、Byba函数名调用格式数学意义指数分布Y=Exppdf(X, )正态分布Y=Normpdf(X, )c2分布Y=chi2pdf(X,n)分布Y=Gampdf(X,A,B)b分布Y=betapdf(X,A,B)说明：（1）其他参数与“期望和方差”表中的相同。 （2）输入对应向量参数的维数必须相等。 pdf(norm,0.7733,0,1) normpdf(0.7733,0,1)例例 计算正态分布计算正态分布N（0，1）在点）在点0.7733下的值。下的值。clfx=0:0.1:30y1=chi2pdf(x,1);plot(x,y1, :)hold on;y2=chi2pdf(x,5)plot(x

24、,y2, +)y3=chi2pdf(x,15)plot(x,y3, o)axis(0,30,0,0.2)例例 绘制卡方分布密度函数在绘制卡方分布密度函数在n分别等于分别等于1，5，15的图。的图。概率值函数 在概率统计计算中，常常要做这样的计算，即给出分布类型、分布参数及计算点，需要得到在计算点的概率值。在MATLAB中，这种计算由概率累积函数cdf来实现。函数名函数名调用形式调用形式函数注释函数注释cdfP=cdf(name,A1,A2,A3)估值选定的累积分布函数估值选定的累积分布函数nctcdfP=nctcdf(X,nu,delta)非中心非中心T T累积分布函数累积分布函数gamcdf

25、 P=gamcdf(X,A,B)g g累积分布函数累积分布函数nxtcdfP=nxtcdf(X,v,delta)非中心卡方累积分布函数非中心卡方累积分布函数函数名函数名调用形式调用形式函数注释函数注释BetacdfP=betacdf(X,A,B)b累积累积分布函数分布函数binocdfY=binocdf(X,N,P)二项式累积函数二项式累积函数chi2cdfP=chi2cdf(X,v)卡方累积分布函数卡方累积分布函数expcdfP=expcdf(X,v)指数累积分布函数指数累积分布函数fcdfP=fcdf(X,v)F累积分布函数累积分布函数gamcpdY=gamcdf(X,P)几何累积分布函数

26、几何累积分布函数geocdfP=geopdf(X,M,K,N)超几何累积分布函数超几何累积分布函数hygecdfP=hygepdf(X,mu,sigma)对数累积分布函数对数累积分布函数logecdfY=logepdf(X,R,P)负二项式累积分布函数负二项式累积分布函数ncfcdfP=ncfcdf(X,nu1,nu2,delta)非中心非中心F累积分布函数累积分布函数 poisscdf(0,1.5) poisscdf(0,2.5)例例 某一公安局在长度为某一公安局在长度为t的时间间隔内收到的紧急呼救的时间间隔内收到的紧急呼救次数服从参数为次数服从参数为t/2的泊松分布，而与时间间隔的起点的泊

27、松分布，而与时间间隔的起点无关（时间以小时计）。求：（无关（时间以小时计）。求：（1）在某天中午）在某天中午12时至时至下午下午3时没有收到呼救的概率；（时没有收到呼救的概率；（2）在某天中午）在某天中午12时时至下午至下午5时至少收到一次紧急呼救的概率。时至少收到一次紧急呼救的概率。 中午12时至下午5时至少收到一次呼救的概率为（1-0.0821）=0.9179。分值点函数（逆概率累积函数） 已知分布及分布中的一点，求此点处的概率值要用到概率累积函数cdf；当已知概率值而要求对应该概率的分布点时，就要用到逆概率累积函数。函数名调用形式函数注释betainvX= betainv(P,A,B)b

28、 b分布逆累积分布函数分布逆累积分布函数binoinvX=binoinv(Y,N,P)二项式二项式逆累积分布函数逆累积分布函数chi2invX=chi2inv(P,v)卡方卡方逆累积分布函数逆累积分布函数expinvX=expinv(P,mu)指数指数逆累积分布函数逆累积分布函数FinvX=finv(P,v1,v2)F逆累积分布函数逆累积分布函数gaminvX=gaminv(P,A,B)g g逆累积分布函数逆累积分布函数geoinvX=geoinv(Y,P)几何逆累积分布函数几何逆累积分布函数函数名调用形式函数注释hygeinvX=hygeinv(P,M,K,N)超几何超几何逆累积分布函数逆累

29、积分布函数icdfX=icdf(name,P,A,B,C)选定的选定的逆累积分布函数逆累积分布函数logeinvX=logeinv(P,mu,sigma)对数对数逆累积分布函数逆累积分布函数nbininvX=nbininv(Y,R,P)负二项式负二项式逆累积分布函数逆累积分布函数ncfinvX=ncfinv(P,nu1,nu2,delta)非中心非中心F逆累积分布函数逆累积分布函数nctinvX=nctinv(P,nu,delta)非中心非中心T逆累积分布函数逆累积分布函数ncx2invX=ncx2inv(P,v,delta)非中心卡方非中心卡方逆累积分布函数逆累积分布函数续表续表例：由某商店

30、过去的销售记录知道，某种商品每月的销例：由某商店过去的销售记录知道，某种商品每月的销售数可以用参数售数可以用参数l l=25的泊松分布来描述，为了有的泊松分布来描述，为了有95%以以上的把握不使商品脱销，问商店在每月底应该进该种商上的把握不使商品脱销，问商店在每月底应该进该种商品多少件？品多少件？ 由计算结果知，这家商店只要在月底进货该种商品33件，就可以有95%以上的把握保证该种商品在下月不脱销。 poissinv(0.95,25) ans = 33例：生产线生产某种产品每个批次的数量为例：生产线生产某种产品每个批次的数量为1000件，质件，质检部门对每个批次的产品都要从中抽取检部门对每个批

31、次的产品都要从中抽取50件产品进行检件产品进行检验，看其是否存在缺陷。如果每个批次中存在缺陷的产验，看其是否存在缺陷。如果每个批次中存在缺陷的产品数不超过品数不超过10件，问从该批所抽的件，问从该批所抽的50件产品中有缺陷的件产品中有缺陷的产品的数量最多是多少？产品的数量最多是多少？ 设从该批所抽的50件产品中所含的有缺陷的产品数量为X,则X服从超几何分布Hyge(x,1000,10,50),以99%的概率估算有Hygeinv(0.99,1000,10,50) ans = 3 结果表明从该批所抽的50件产品中所含的有缺陷的产品最多是3件，可信度为99%。随机数生成函数函数名调用形式函数注释be

32、tarndR=betarnd(A,B)b分布的随机矩阵分布的随机矩阵binorndR=binornd(N,P,MM,NN)二项式分布的随机矩阵二项式分布的随机矩阵chi2rndR=chi2rnd(v)卡方分布的随机矩阵卡方分布的随机矩阵frndX=frnd(v1,v2)F分布的随机矩阵分布的随机矩阵georndR=geornd(P)几何分布的随机矩阵几何分布的随机矩阵hygerndR=hygernd(M,K,N)超几何分布的随机矩阵超几何分布的随机矩阵mvnrndR=mvnrnd(mu,sigma,CASES)多元正态分布的随机矩阵多元正态分布的随机矩阵normordR=normrnd(mu,

33、sigma)正态分布的随机矩阵正态分布的随机矩阵unifrndR=unifrnd(A,B)连续均匀分布的随机矩阵连续均匀分布的随机矩阵trndR=trnd(v)T（学生学生）分布的随机矩阵分布的随机矩阵9.2.3 9.2.3 参数估计参数估计 利用样本对总体进行统计推断主要有两大类，一类是参数估计，另一类是检验函数估计。估计的问题是多种多样的，如总体参数的估计、分布函数和密度函数的估计等。估计的方法也很多，如最大似然方法、矩估计法、最小二乘法等等。在MATLAB中用于参数估计的函数见下表。函数名函数名调用形式调用形式函数注释函数注释betafitbetafit(X)PRAT,PCI=betafit(X,ALPHA)返回返回b b分布的最大似然估计；分布的最大似然估计；返回最大似然估计值和返回最大似然估计值和a a水平水平的置信区间的置信区间binofitbinofit(X,N)PHAT,PCI=binofit(X,N,ALPHA)二项分布的最大似然估计；二项分布的最大似然估计；a a水