Matlab经济应用.ppt_三一办公31ppt.com

资源描述

《Matlab经济应用.ppt》由会员分享，可在线阅读，更多相关《Matlab经济应用.ppt（74页珍藏版）》请在三一办公上搜索。

1、数组操作与代数运算,MATLAB应用2,内容,数组（向量、矩阵）的建立数组元素的查找常见的代数运算初等变换的实现线性方程组特征值、特征向量,3,向量、矩阵(数组)的建立,在MATLAB中最常用的基本数据结构就是一、二维数组，而标量反到是数组的特例(11数组)。下面我们先来讨论一下常规数组的特例。向量(一维数组)的建立矩阵(二维数组)的建立,4,向量、矩阵的直接输入,输入矩阵的最基本方法直接输入法直接输入矩阵的元素,用方括号表示矩阵,同一行元素间用空格或逗号分隔,不同行间用分号或回车符分隔。例如：A=1,2,3;4,5,6;7,8,9B=0 2 4,5,向量的快捷输入,利用冒号“:”生成等差数列

2、a=i:j%初值：终值a=i:k:j%初值：步长：终值利用linspace函数生成线性等分行向量a=linspace(i,j)%生成100个在i,j之间等间隔分布的元素a=linspace(i,j,n)利用logspace函数生成对数等分行向量a=logspace(i,j)%生成50个从10 i 10 j 对数等间隔元素a=logspace(i,j,n),6,矩阵的快捷生成,用函数建立矩阵ones(m,n);%生成mn的1矩阵zeros(m,n);%生成mn的零矩阵eye(n);%生成n阶单位矩阵rand(m,n);%生成mn的0,1均匀分布随机矩阵randn(m,n);%生成mn的标准正态分

3、布随机矩阵用小矩阵合成大矩阵（分块矩阵方法）A=A1，A2；A3，A4,7,标准数组生成函数,存取访问EXCEL文件中的数据,读出数据：N=xlsread(filename,sheet,range)写入数据：xlswrite(filename,M,sheet,range)示例：A=ones(2,4)xlswrite(d:tmp.xls,A,sheet1,A1:D2)B=xlsread(d:tmp.xls,sheet1,A1:D3),9,矩阵的调用、访问元素,MATLAB可以用矩阵的名称调用全矩阵,也可用下标访问矩阵的某个或某些元素(注意下标可以是向量)。例如：B=A;a=A(2,3)A(:,3

4、),A(2,:)A(1:4,1:3),A(1 3 5,2,4)可以用删除矩阵的元素，如：A(:,2)=,A=,10,数组元素的标识与访问,MATLAB数组元素通过数组下标来标识，起始下标为1；访问数组元素可通过数组名后跟圆括号内的下标(对)来实现。如：A(2),B(1,3),C(1,3,5,2),D(:,2)等。注意：数组的索引下标可以是普通向量和逻辑向量(分量为逻辑值,必须由关系运算产生)，此时表示的是原数组的一个子数组。例：x=-3:5;L=abs(x)2;y=x(L)I=1,3,4;x=x(I),11,元素“全下标”与“单下标”表示法,二维以上数组在标识元素时，可采用“全下标”或“单下标

5、”方式。全下标表示法：就是完整的标出元素的的各维下标(逗号分隔)，即指出是“第几行,第几列”的元素。单下标表示法：将数组的所有列按先左后右顺序，从上到下排成一列，进行编号，单下标就是按此顺序索引元素。(公式：l=(c-1)m+r)转换函数：ind=sub2ind(size(A),i,j)全下标=单下标I,J=ind2sub(size(A),ind)单下标=全下标,12,数组及元素的寻访与赋值总结,13,数组寻访、赋值示例,rand(state,0)%初始化,置0.x=rand(1,5)x(3:end)x(5:-1:1)%倒序x(find(x0.5)x(1 2 3 4 4 3 2 1)x(1 4

6、)=1 1A=zeros(2,5);A(:)=-4:5L=abs(A)3islogical(L)X=A(L),A=zeros(2,4)A(:)=1:8s=2 3 5;A(s)Sa=10 20 30A(s)=SaA(:,2 3)=ones(2),14,数组运算与矩阵运算,MATLAB提供的运算可分数组运算与矩阵运算.从结构形状看,二维数组与矩阵无区别.但,矩阵的运算有明确的线性代数规则;而数组运算是MATLAB自定义的,针对数组每个元素的分别操作,目的是便于数据管理、使用简便、自然、有效。数组运算也可引入到函数定义中,即 X=xij mn,则 f(X)=f(xij)mn。,15,数组运算（点运算

7、）,为了能实现对应元素之间的四则运算、乘方等,MATLAB提供了一种“数组运算”,识别特征是在原有运算符前加一个点号“.”。主要有：.*、./、.、.、.等.例如：A=magic(4);B=2*ones(4);k=3;A.*B A./B,A.BA.B,A.k,k.AA.,16,数组中元素的查找,用find()函数可以查找数组中满足条件的元素下标。如：i=find(A=3)%返回数值为3元素的单下标。s,t=find(A=3)%返回数值不超过3元素的全下标。A(s,t)用比较运算获取逻辑下标来访问数组也可实现数组元素查找。I=(A=3)A(I)B=A(abs(A-2)=1),17,数组元素求和,

8、用sum(A,n)函数可以对数组元素沿第n维(增加)方向求和：sum(A)表示对一维数组求和;若A是多维数组时与sum(A,1)相同,即沿第一维方向求和(按列求和)。sum(A,2)表示对多维数组沿第二维方向求和(按行求和)。思考:若A为一矩阵,sum(sum(A)表示何意?,18,数组元素的最值,min()、max()函数可以求数组元素的最值:min(A)-若A是一维数组则求出最小值;若A是多维数组则分别求第一维方向上的最小值.min(X,Y)-比较同结构数组,求出对应元素的最小值.min(A,dim)-沿第dim维分别求出最小值.Y,I=min(.)-求出最小值外,还要该维上的相应下标.m

9、ax()用法与min类似,19,执行数组运算的常用函数,取整、求余函数以及坐标变换函数。,20,数组操作函数,关系运算,关系运算符主要用来对数与矩阵、矩阵与矩阵进行比较，并返回反映二者之间大小关系的由数0和1组成的矩阵。基本的关系运算符主要有:、=、=、=这6个。A=1,2;3,4,B=0,2;4,5C=A=BD=A=BA=B,逻辑运算,逻辑运算在计算机语言中是普遍存在的.在Matlab中包含与、或、非、异或4种基本的逻辑运算。逻辑表达式和逻辑函数的值应该为一个逻辑量“真”或“假”。MatLat系统在给出逻辑运算的结果时，以逻辑值“1”代表逻辑“真”，以“0”代表“假”，但在判断一个量是否为“

10、真”时，以0代表“假”，以任意的非零值代表“真”。Matlab的逻辑运算也是以矩阵为基本运算单元的。,23,关系、逻辑函数(1),24,对矩阵、数组进行操作、运算示例,【例】标准数组产生的演示。ones(1,2)randn(state,0)randn(2,3)E=eye(3)diag(E)diag(diag(E)A=diag(-5:5)B=diag(ones(10,1),1)C=diag(ones(10,1),-1)D=A+B+C,【例 2.8-2】reshape的使用演示。a=-4:4A=reshape(a,3,3)【例】数组转置、对称交换和旋转操作后果的对照比较。A.flipud(A)fl

12、视化中，裁剪图形。“空”数组：某(几)维长度为零的数组。经常用于改变数组大小。如：a=,b=ones(2,2,0),A(:,1)=,28,字符串数组（字符型变量）,字符串数组的内容是字符或字符串，它与数值型数组是不同的类型。可以用函数class(Arrary)判别数组Arrary的类型。可用赋值法直接建立字符串数组,只要将所赋值内容用单引号对括住，MATLAB 就识别为字符串。内存中存放一个字符需要2个字节，而存放一个数值标量需要8个字节。例：字符串建立与倒排序。a=This is a VIP test!,size(a)b=a(end:-1:1),size(b),a(12:14),常用字符串操

13、作,删除字符串结尾处的空格：S=deblank(C)输入空格符：blanks(n)产生n个空格将字符串进行大小写转换：upper(S)lower(S)将字符串作为命令执行：a=eval(字符串表达式)字符串比较操作：strcmp(stringl,string)返回0,1字符串中查找子串：findstr(str1,str2)返回短串的开始位置字符串替换操作：strrep(s1,s2,s3)用s3替换s1中的s2字符串与数值转换：int2str(A)、num2str(A)、str2num(S)abs(str)把串翻译成ASCII码;setstr(asc)把ASCII码翻译为串示例：a=Hello“

14、World”s1=upper(a)s2=lower(a),字符串操作示例,abs(a)a=This is a good man!s=a,blanks(2)deblank(s)strrep(s,good,bad)i=findstr(s,is)b=eval(1,2;3,4),应用举例,fid=fopen(data.txt,r)%打开数据文件返回文件句柄line=fgetl(fid)%读取文件句柄对应文件的一行字符k=findstr(line,cgg)%返回字符串中字串cgg位置s=length(k)line=fgetl(fid)k=findstr(line,cgg)s=length(k)fclos

15、e(fid)%关闭文件句柄对应的文件,32,单元数组(Cell arrary),单元数组的基本组分(Element)是单元(cell)。每个单元本身在数组中地位平等,它们只能以下标区分。单元中可以存放任何类型、任何大小的数组。同一单元数组的各单元内容允许存放类型、大小都不同数据。单元数组的维数定义不受限制，数组对各单元的编址方法也有“全下标”与“单下标”两种。对单元数组而言,单元与单元的内容是不同范畴的东西。访问单元靠“单元索引(cell index)”,即带圆括号的下标;而访问单元的内容靠“单元内容编址(content addressing)”,即用带花括号的下标访问单元内容。,33,直接赋

16、值法建立单元数组,一、按单元索引进行赋值如：A(1,1)=1,2;3,4,A(1,2)=abc二、按内容编址进行赋值如：B1,1=1,2;3,4,B1,2=abc注意：花括号用于访问单元的内容，圆括号（）用于访问单元。如：e=A(1,1),f=A1,1,34,cell专用函数法建立单元数组,使用cell函数可为单元数组先定制外形，然后再赋值成为一个具体的单元数组。如：c=cell(2,3)c(1,1)=This does workc2,4=This works two注：cellfun可将一指定的函数应用到单元数组的所有单元上。如：cellfun(length,A)cellfun(isclas

17、s,A,char),35,矩阵的四则运算、乘方,矩阵相加减：A+B；A-B；矩阵相乘：A*B矩阵左除：AB;相当与求AX=B的解矩阵右除：B/A;相当与求XA=B的解方阵的乘方：Ap(p为0或正负整数)(注：p非整数时方阵的乘方涉及矩阵的特征值、特征向量*),36,其它常用矩阵运算,det(A)方阵的行列式A共扼转置 A.转置inv(A)逆矩阵orth(A)正交阵rank(A)秩trace(A)迹d=eig(A),V,D=eig(A)特征值、标准正交化的特征向量组、特征值poly(A)特征多项式,代数应用-行列式的试算,n阶行列式-4、5阶行列式-估计n阶行列式syms aA4=a,0,0,1

18、;0,a,0,0;0,0,a,0;1,0,0,aA5=a,0,0,0,1;0,a,0,0,0;0,0,a,0,0;0,0,0,a,0;1,0,0,0,adA4=det(A4),dA5=det(A5)syms x a0 a1 a2 a3 a4B4=a0,-1,0,0;a1,x,-1,0;a2,0,x,-1;a3,0,0,xB5=a0,-1,0,0,0;a1,x,-1,0,0;a2,0,x,-1,0;a3,0,0,x,-1;a4,0,0,0,xdB4=det(B4),dB5=det(B5)syms a b c dA=a b c d;-b a-d c;-c d a-b;-d-c b asimplif

19、y(det(A),38,矩阵的秩与逆矩阵、伪逆矩阵,用rank(A)可求矩阵A的秩。用inv(A)函数可以求出方阵A的逆矩阵。用pinv(A)函数可以求一般矩阵A的伪逆矩阵。(在求Ax=b的解时，若A不是可逆阵，则可用pinv(A)*b表示其一种最小二乘解。),39,行列式、矩阵计算示例,用det函数可以求矩阵的行列式.例如:A=1,3,5;2,4,6;1,1,1a=det(A),40,矩阵的迹与对角元素的提取,用trace(A)可求矩阵A的迹(主对角元素之和)用diag(A)可提取矩阵的对角元素，构成列向量。用diag(A,n)可提取与主对角线平行的右上第n条线上的矩阵元素(n为负数则为左下

20、第n条)。用diag(x)可将向量x作为对角元素构造对角矩阵;diag(x,n)将向量作为平行主对角第n条线上的元素构造矩阵。,41,矩阵的初等变换,A(i,j,:)=A(j,i,:)A(i,:)=k*A(i,:)A(i,:)=A(i,:)+k*A(j,:)A(:,i,j)=A(:,j,i)A(:,i:)=k*A(:,i)A(:,i)=A(:,i)+k*A(:,j)辅助指令“format rat”可将数据用分数格式显示。练习用inv(A)指令和初等变换分别求A的逆矩阵,其中,42,化矩阵为最简阶梯形矩阵,化最简阶梯形矩阵是求解线性方程组和判断向量线性相关/无关的基础,MATLAB除了用初等变换

21、方法化最简阶梯形矩阵外，还有更直接方便的指令rref来实现。rref(A)例如 A=1,-2,-1,0,2;-2,4,2,6,-6;2,-1,0,2,3;3,3,3,3,4;format rat rref(A),43,向量组的相关性和极大无关组,要判断向量组1,2,n的相关性,并求其极大无关组,可先将向量组作为矩阵A的列向量,再将A化为最简阶梯形;最简阶梯形中基本单位向量所在列对应的向量就是极大无关组,由初等变换保持向量间的线性关系,易知其余向量如何用极大无关组表示。练习求1,2,n的极大无关组,并指出其余向量如何用其表示。其中,44,化矩阵为最简阶梯形矩阵,化最简阶梯形矩阵是求解线性方程组和

22、判断向量线性相关/无关的基础,MATLAB除了用初等变换方法化最简阶梯形矩阵外，还有更直接方便的指令rref来实现。rref(A)例如 A=1,-2,-1,0,2;-2,4,2,6,-6;2,-1,0,2,3;3,3,3,3,4;format rat rref(A),45,向量组的相关性和极大无关组,要判断向量组1,2,n的相关性,并求其极大无关组,可先将向量组作为矩阵A的列向量,再将A化为最简阶梯形;最简阶梯形中基本单位向量所在列对应的向量就是极大无关组,由初等变换保持向量间的线性关系,易知其余向量如何用极大无关组表示。练习求1,2,n的极大无关组,并指出其余向量如何用其表示。其中,46,线

23、性方程组求特解,指令x=Ab可求线性方程组Ax=b的特解。(假如A不是可逆阵,即方程组有可能是欠定或超定的,则x是众多解中(含0最多)的一个或最小二乘解。)用x=lsqnonneg(A,b)可求出方程组Ax=b的一个非负最小二乘解。例：A=2 1 3;1 3 2;b=5;10;x=lsqnonneg(A,b),47,求极大无关组,要求向量组1,2,n的极大无关组,可将它们作为列向量组构造矩阵A,然后用指令rref(A)将A化为行简化矩阵,再由其中的单位向量容易判断出那些是极大无关向量。例如：A=2 1 3;1 3 2;rref(A),48,求齐次线性方程组基础解系,用null(A)指令可以求出

24、齐次线性方程组Ax=0的基础解系。例如：A=1,-2,3,-4;0,1,-1,1;-1,0,-1,2;1,-3,4,-5;B=null(A)B中各列就是Ax=0的基础解系。,49,线性方程组Ax=b求通解,方法一：用null(A)求出相应齐次线性方程组基础解系;用x0=Ab得到原方程组的一个特解;最后合成为原方程组的通解。方法二：构造增广矩阵=A,b;用指令rref()得到行简化矩阵,直接利用行简化阵写出通解形式。,50,线性无关向量组的标准正交化,将线性无关向量组作为列向量组构成矩阵A,再用orth(A)指令可得标准正交列向量组。例如：A=1 1 1;0 1 1;0 0 1;1 0 0ort

25、h(A),51,特征值与矩阵对角化、标准型,用指令V,D=eig(A)可以对方阵A生成其特征值为对角元素的对角阵D以及相应于特征值的特征向量为列向量的矩阵V,它们间满足关系AV=VD。假如A是实对称矩阵,则V必为正交矩阵,即V的列向量组是标准正交基。实二次型经正交替换化标准型的过程可以用上述MATLAB指令来实现。简化指令E=eig(A)将返回由方阵A特征值构成的向量E。,52,若当标准型与若当化*,用指令V,J=jordan(A)可求出方阵A的若当矩阵J并满足VA*V=J。例如：A=gallery(5);V,J=jordan(A),53,矩阵的LU分解*,指令LL,U,P=lu(A)可得到一

26、个上三角阵U、一个下三角阵LL和一个置换矩阵P,使得LL*U=P*A。指令L,U=lu(A)可得到一个上三角阵U和一个矩阵L使得A=L*U,其中L=inv(P)*LL。,54,矩阵的QR分解*,QR分解即矩阵的正交三角分解,它将矩阵分解为一个正交矩阵和一个上三角矩阵的乘积。指令Q,R=qr(A)将得到一个与A同阶的上三角阵R和一个酉矩阵Q,使得A=QR。指令Q,R,S=qr(A)将得到一个置换矩阵S、一个对角元素递减的上三角阵R和一个酉矩阵Q,使得AS=QR。,55,矩阵的奇异值分解*,指令s=svd(A)将得到矩阵A的奇异值组成的向量s。指令U,S,V=svd(A)将得到一个与A具有相同维数

27、的矩阵S,其对角元素为递减非负值;同时得到酉矩阵U和V,使得A=USV。,56,多项式的运算*,多项式与其系数行向量1-1对应,故MATLAB可用行向量来表示降幂排列的多项式。polyval(p,x);polyvalm(p,A)求值roots(p)求根conv(a,b)多项式相乘g,r=deconv(c,b)%相除：c/b=gb+rpolyder(p)求导,57,方阵的特征多项式*,p=poly(A)指令将求出方阵的特征多项式。p=poly(r)指令将求出以向量r中元素为根的多项式p(p为按降幂方式排列的多项式系数构成的向量,MATLAB中常采用向量来表示多项式！如用5维向量向量表示4次多项式

28、。)。相关指令：r=roots(p)表示求出多项式p的根r。相关指令：polyval(p,X)表示计算数组X关于多项式的值。(此处执行数组运算.若要进行矩阵运算则用polyvalm指令.),迁徙模型(矩阵乘法、特征值应用),设在一个大城市中的总人口是固定的。人口的分布则因居民在市区和郊区之间迁徙而变化。每年有6%的市区居民搬到郊区去住，而有2%的郊区居民搬到市区。假如开始时有30%的居民住在市区，70%的居民住在郊区，问10年后市区和郊区的居民人口比例是多少？30年、50年后又如何？这个问题可以用矩阵乘法来描述。把人口变量用市区和郊区两个分量表示，一年以后，市区人口为xc1(10.06)xc0

29、0.02xs0，郊区人口xs1 0.06xc0(10.02)xs0，,问题的矩阵描述,用矩阵乘法来描述，可写成：从初始到k年，此关系保持不变，因此上述算式可扩展为，故可用Matlab程序进行计算：A0.94,0.02;0.06,0.98,x00.3;0.7x1A*x0,x10A10*x0,x30A30*x0,x50A50*x0得到：,本题特征值和特征向量的意义,无限增加时间k，市区和郊区人口之比将趋向一组常数0.25/0.75。为了弄清为什么这个过程趋向于一个稳态值，我们改变一下坐标系统。在这个坐标系统中可以更清楚地看到乘以矩阵A的效果，先求A的特征值和特征向量，得到令它是特征向量的整数化，得

30、到,迁徙模型(动物繁殖),某农场饲养的某种动物所能达到的最大年龄15 岁，将其分成3个年龄组，第一组:05 岁;第二组:610 岁;第二组 1115 岁.动物从第二年龄组起开始繁殖后代，经过长期统计，第二年龄组的动物在其年龄段平均繁殖4个后代，第三年龄组的动物在其年龄段平均繁殖3个后代，第一年龄组和第二年龄组的动物能顺利进入下一个年龄组的存活率分别为 1/2、1/4，假设农场现有3个年龄段的动物各1000 头，问 15 年后农场3个年龄段的动物各有多少头?因为年龄组分为5岁一段，故将时间周期也取为5年，15 年后经过了3个时间周期.设 xi(k)表示k个时间周期第i组年龄阶段动物的数量(k=1

31、,2,3).这类应用问题实际上是所谓马尔可夫过程的一个类型。所得到的向量序列x1,x2,.,xk称为马尔可夫链。马尔可夫过程的特点是k时刻的系统状态xk完全可由其前一个时刻的状态xk-1所决定，与k-1时刻之前的系统状态无关。,中草药配伍-极大无关组应用,问题：某中药厂用9种中草药（A-I），根据不同的比例配制成了7种特效药，各用量成分见表1（单位：克）,1现3号药与6号药脱销，能否用其他5种特效药配置出来？2新药v1、v2、v3含9种中草药量如下表，能否用上面7中特效药配制？,解决方案,把每一种特效药ui看成一个九维列向量，分析7个列向量构成向量组的线性相关性。若向量组线性无关，则无法配制脱

32、销的特效药；若向量组线性相关，并且能找到不含u3,u6的一个最大线性无关组，则可以配制3号和6号药品。Matlab程序：u1=10;12;5;7;0;25;9;6;8;u2=2;0;3;9;1;5;4;5;2;u3=14;12;11;25;2;35;17;16;12;u4=12;25;0;5;25;5;25;10;0;u5=20;35;5;15;5;35;2;10;0;u6=38;60;14;47;33;55;39;35;6;u7=100;55;0;35;6;50;25;10;20;U=u1,u2,u3,u4,u5,u6,u7,U0,r=rref(U)易知u1、u2、u4、u5、u7是极大无关

33、组，u3、u6能配制！v1=40,67,14,44,53,50,71,41,14;3=88,67,8,51,7,80,38,21,30;v2=162,141,27,102,60,155,118,68,52;.(用解方程组方法怎么做？),应用矩阵编制Hill密码,在密码学中将信息代码称为密码，没有转换成密码的文字信息称为明文，把密码表示的信息称为密文。从明文转换为密文的过程叫加密，反之则为解密。1929年，希尔（Hill）通过矩阵理论对传输信息进行加密处理，提出了在密码学史上有重要地位的希尔加密算法。,希尔加密算法(用Matlab做个实例),假设我们要发出“attack”这个消息。首先把每个字母

34、a，b，c，dx，y，z映射到数1，2，3，424，25，26。例如1表示a，3表示c，20表示t，11表示k，另外用0表示空格，用27表示句号等。于是可以用以下数集来表示消息“attack”：,把这个消息按列写成矩阵的形式：,第一步：“加密”工作。现在任选一个三阶的可逆矩阵，如:,把矩阵M经过乘以A变成“密码”（B）后发出。,第二步：“解密”。解密是加密的逆过程：M=A-1B,矩阵在文献管理中的应用,搜索与多个指定关键词匹配度高的书名,=A,现要搜索书名内与三个关键词“应用，线性，代数”匹配度高的书籍。,每本书与指定三个关键词的匹配度为：P=Ax,用Matlab实现此例!,67,多项式乘法与

35、向量卷积*,用指令w=conv(u,v)可以求出多项式u和v的乘积w，它也是向量u和v的卷积。多项式乘积与向量卷积算法同构。,68,多项式除法与向量的解卷运算*,用指令q,r=deconv(v,u)可求出多项式v除以u所得的商多项式q和余数多项式r，这也可以看作对向量v,u进行解卷运算。例如：p1=3,4,6;p2=5,2,-4,7;p=conv(p1,p2)q,r=deconv(p,p1),69,多项式求微分*,用指令d=polyder(p)可求出多项式p的导数。用指令d=polyder(a,b)可求出多项式a与b乘积的导数多项式d。用指令q,d=polyder(v,u)可求出多项式v与u的商v/u的导数分式q/d。,70,部分分式展开*,r,p,k=residue(v,u)指令可将分式v(x)/u(x)化为：k(x)+r1/(x-p1)+r2/(x-p2)+注意:若pj是u(x)的m重根,则有相应m项,即 rj1/(x-pj)+rj2/(x-pj)2+rjm/(x-pj)m例如：u=1,1,-1,-1;v=1,0,-2,2,1;r,p,k=residue(v,u),练习二,参考：for i=1:4 eval(a,int2str(i),=A(:,int2str(i),);end,续1,续2,续3,

展开阅读全文