matlab在统计分析中的应用.ppt

MatlabMatlab在统计分析中的应用在统计分析中的应用MATLAB概述与运算基础概述与运算基础nMATLAB语言：广泛应用于工程计算及数值分析领域的新型高级语言n1984年由美国 MathWorks 公司推向市场n已成为国际公认的最优秀的工程应用开发环境n功能强大、简单易学、编程效率高n在欧美各高等院校，MATLAB已成为线性代数、自动控制理论、数字信号处理、时间序列分析、动态系统仿真、图像处理等课程的基本教学工具，成为大学生、硕士生以及博士生必须掌握的基本技能MATLAB：MATrix LABorotory一一. . MATLAB特点特点: :一一. MATLAB. MATLAB特点特点: :1. 数值计算和符号计算功能n数值计算功能包括：矩阵运算、多项式和有理分式运算、数据统计分析、数值积分、优化处理等n符号计算将得到问题的解析解2.MATLAB语言n除了命令行的交互式操作以外，还可以程序方式工作使用MATLAB很容易地实现C或FORTRAN语言的几乎全部功能3.图形功能n两个层次的图形命令：一种是对图形句柄进行的低级图形命令，另一种是建立在低级图形命令之上的高级图形命令。

利用MATLAB的高级图形命令可以轻而易举地绘制二维、三维乃至四维图形，并可进行图形和坐标的标识、视角和光照设计、色彩精细控制等等4.应用工具箱n包括：基本部分和各种可选的工具箱n基本部分中有数百个内部函数n工具箱分为两大类：功能性工具箱和学科性工具箱a. 功能性工具箱：主要用来扩充符号计算功能、可视建模仿真功能及文字处理功能等b. 学科性工具箱：专业性较强，如控制系统工具箱、信号处理工具箱、神经网络工具箱、最优化工具箱、金融工具箱等，用户可以直接利用这些工具箱进行相关领域的科学研究n命令窗口的命令编辑区 用于输入命令和显示计算结果n例：键入 A=[1 2 3;4 5 6;7 8 9] 输出 A = 1 2 3 　　　　　 4 5 6 　　　　 7 8 9二二. MATLAB. MATLAB集成环境集成环境矩阵运算矩阵运算n【例】 求解线性方程组AX=B 　　 1 1.5 2 9 7 　　 0 3.6 0.5 -4 4 其中A= 7 10 -3 22 33 　　 3 7 8.5 21 6 　　 3 8 0 90 -20 3 -4n B= 20 5 16 n 在MATLAB命令窗口输入命令： n a=[1,1.5,2,9,7; 0,3.6,0.5,-4,4; 7,10,-3,22,33; 3,7,8.5,21,6; 3,8,0,90,-20];n b=[3; -4; 20; 5; 16];n x=a\bn得到的结果是：n x =n 3.5653n -0.9255n -0.2695n 0.1435n 0.0101【例】 求方程 x^4+7x^3 +9x-20=0的全部根。

n命令窗口输入：n p=[1,7,0,9,-20]; %建立多项式系数向量n x=roots(p) %求根n得到的结果是：n x =n -7.2254n -0.4286 + 1.5405in -0.4286 - 1.5405in 1.0826二维图形二维图形n一、一、 plot函数函数n函数格式：plot(x,y) 其中x和y为坐标向量n函数功能：以向量x、y为轴，绘制曲线n【例例】 在区间0≤X≤2内，绘制正弦曲线Y=SIN（X）x=0:pi/100:2*pi;y=sin(x);plot(x,y)n【例例】同时绘制正、余弦两条曲线Y1=SIN(X)和Y2=COS(X)x=0:pi/100:2*pi;y1=sin(x);y2=cos(x);plot(x,y1,x,y2)n（一）线型与颜色（一）线型与颜色n格式：格式：plot(x,y1,’cs’,...)n其中c表示颜色， s表示线型一、一、 plot函数函数n【例】 用不同线型和颜色重新绘制上例图形 x=0:pi/100:2*pi;y1=sin(x);y2=cos(x);plot(x,y1,'go',x,y2,'b-.')其中：参数‘go’和‘b-.’表示图形的颜色和线型。

g：绿色 o：线型为圆圈b：蓝色 -.：表示图形线型为点划线n（二）图形标记（二）图形标记ntitle(‘加图形标题'); nxlabel('加X轴标记'); nylabel('加Y轴标记'); ntext(X,Y,'添加文本'); 一、一、 plot函数函数三维图形三维图形n一、一、 plot3函数函数函数格式：plot3(x1,y1,z1,c1,x2,y2,z2,c2,……) 其中x1,y1,z1…表示三维坐标向量；c1,c2…表示线形或颜色n函数功能：以向量x，y，z为坐标，绘制三维曲线n【例例】 绘制三维螺旋曲线nt=0:pi/50:10*pi;nx=sin(t),y=cos(t);nplot3(x,y,t);ntitle('helix'),text(0,0,0,'origin');nxlabel('sin(t)'),ylabel('cos(t)'),zlabel('t');ngrid;n二、二、mesh函数函数n绘制三维网格图三维网格图函数格式：mesh(x,y,z,c)n其中：x，y控制X和Y轴坐标 矩阵z是由(x，y)求得Z轴坐标 (x,y,z)组成三维空间的网格点 c用于控制网格点颜色n【例】 绘制三维网格曲面图nx=[0:0.15:2*pi];ny=[0:0.15:2*pi];nz=sin(y')*cos(x); %矩阵相乘nmesh(x,y,z);【例例】画出由函数形成的立体网状图：x=linspace(-2,￿2,￿25);￿%￿在x轴上取25点 y=linspace(-2,￿2,￿25);￿%￿在y轴上取25点 [xx,yy]=meshgrid(x,y);￿%￿xx和yy都是21x21的矩阵 zz=xx.*exp(-xx.^2-yy.^2);￿%￿计算函数值，zz也是21x21的矩阵 mesh(xx,￿yy,￿zz);￿%￿画出立体网状图￿n三、三、surf函数函数n绘制三三维维曲曲面面图图，各线条之间的补面用颜色填充。

surf函数和mesh函数的调用格式一致n函数格式: surf (x,y,z)n其中x，y控制X和Y轴坐标，矩阵z是由x，y求得的曲面上Z轴坐标n【例例】 绘制三维曲面图nx=[0:0.15:2*pi];ny=[0:0.15:2*pi];nz=sin(y')*cos(x); %矩阵相乘nsurf(x,y,z);nxlabel('x-axis'),ylabel('y-axis'),zlabel('z-axis');ntitle('3-D surf');【例例】剔透玲珑球剔透玲珑球[X0,Y0,Z0]=sphere(30); %产生单位球面的三维坐标X=2*X0;Y=2*Y0;Z=2*Z0; %产生半径为2的球面的三维坐标surf(X0,Y0,Z0); %画单位球面shading interp %采用插补明暗处理hold on; mesh(X,Y,Z);hold off %画外球面hidden off %产生透视效果axis off %不显示坐标轴【例例】卫星返回地球的运动轨线示意。

R0=1;%以地球半径为一个单位以地球半径为一个单位a=12*R0;b=9*R0;T0=2*pi; %T0是轨道周期是轨道周期T=5*T0;dt=pi/100;t=[0:dt:T]';f=sqrt(a^2-b^2); %地球与另一焦点的距离地球与另一焦点的距离th=12.5*pi/180;%卫星轨道与卫星轨道与x-y平面的倾角平面的倾角E=exp(-t/20);%轨道收缩率轨道收缩率x=E.*(a*cos(t)-f);y=E.*(b*cos(th)*sin(t));z=E.*(b*sin(th)*sin(t));plot3(x,y,z,'g')%画全程轨线画全程轨线[X,Y,Z]=sphere(30);X=R0*X;Y=R0*Y;Z=R0*Z; %获得单位球坐标获得单位球坐标grid on,hold on,surf(X,Y,Z),shading interp %画地球画地球x1=-18*R0;x2=6*R0;y1=-12*R0;y2=12*R0;z1=-6*R0;z2=6*R0;axis([x1 x2 y1 y2 z1 z2]) %确定坐标范围确定坐标范围view([117 37]),comet3(x,y,z,0.02),hold off %设视角、画运动轨线设视角、画运动轨线 n五、等高线图五、等高线图n函数contour3n【例例】 多峰函数peaks的等高线图n[x,y,z]=peaks(30);ncontour3(x,y,z,16);nxlabel('x-axis'),ylabel('y-axis'),zlabel('z-axis');ntitle('contour3 of peaks')4.6 动画设计动画设计动画功能函数：getframe、moviein和movien【例例】 播放一个不断变化的眼球程序。

nm=moviein(20); %建立一个20个列向量组成的矩阵nfor j=1:20n plot(fft(eye(j+10))) %绘制出每一幅眼球图并保存到m矩阵中n m(:,j)=getframe;nendnmovie(m,10);%以每秒10幅的速度播放画面【例例 】求下列三阶线性代数方程组的近似解MATLAB程序为：A=[2 -5 4;1 5 -2;-1 2 4]；b=[5;6;5]；x=A\bn解解 法法 1：： 分 别 解 方 程 组 (1)Ax=b1；(2)Ay=b2nA=[1 -1 1;5 -4 3;2 1 1];nb1=[2;-3;1];nb2=[3;4;-5];nx=A\b1nx =n -3.8000n 1.4000n 7.2000y=A\b2 -3.6000 -2.2000 4.4000得两个线性代数方程组的解： (1) x1= -3.8， x2= 1.4， x3= 7.2； (2) y1= -3.6， y2=－2.2， y3= 4.4n解解法法2：：将两个方程组连在一起求解：Az=bnb=[2 3;-3 4;1 -5]nz=A\bnz =n -3.8000 -3.6000n 1.4000 -2.2000n 7.2000 4.4000一、 基本统计处理n1、查取最大值、查取最大值nMAX函数的命令格式：n[Y,I]= max (X)：将max(X)返回矩阵X的各列中的最大元素值及其该元素的位置赋予行向量Y与I；当X为向量时，则Y与I为单变量。

n[Y,I]=max(X,[ ],DIM)：按数组X的第DIM维的方向查取其最大的元素值及其该元素的位置赋予向量Y与In【例例】查找下面数列x的最大值nx=[3 5 9 6 1 8] % 产生数列xnx = 3 5 9 6 1 8ny=max(x) % 查出数列x中的最大值赋予yny = 9n[y,l]=max(x) % 查出数列x中的最大值及其该元素的位置赋予y,lny = 9nl = 3n【例例】分别查找下面3×4的二维数组x中各列和各行元素中的最大值nx=[1 8 4 2;9 6 2 5;3 6 7 1] % 产生二维数组xnx = 1 8 4 2n 9 6 2 5n 3 6 7 1ny=max(x) % 查出二维数组x中各列元素的最大值产生赋予行向量yny = 9 8 7 5n2、查取最小值、查取最小值 MIN函数n3、求中值、求中值nY=median(X)：将median(X)返回矩阵X各列元素的中值赋予行向量Y。

若X为向量，则Y为单变量nY=median(X,DIM)：按数组X的第DIM维方向的元素求其中值赋予向量Y若DIM=1，为按列操作；若DIM=2，为按行操作若X为二维数组，Y为一个向量；若X为一维数组，则Y为单变量n【例例】分别求下面数列x1与x2的中值nx1=[9 -2 5 7 12]; % 奇数个元素ny1=median(x)ny1 =n 7nx2=[9 -2 5 6 7 12]; % 偶数个元素ny2=median(x)ny2 =n 6.5000n4、求和、求和n命令格式：nY=sum(X)：将sum(X)返回矩阵X各列元素之和赋予行向量Y；若X为向量，则Y为单变量Y=sum(X,DIM)：按数组X的第DIM维的方向的元素求其和赋予Y若DIM=1，为按列操作；若DIM=2，为按行操作若X为二维数组，Y为一个向量；若X为一维数组，则Y为单变量n例：nx=[4 5 6;1 4 8]nx =n 4 5 6n 1 4 8ny=sum(x,1)ny =n 5 9 14ny=sum(x,2)ny =n 15n 13n5、求平均值、求平均值nMEAN函数调用的命令格式：nY= mean(X)：将mean (X)返回矩阵X各列元素之的平均值赋予行向量Y。

若X为向量，则Y为单变量nY= mean(X,DIM)：按数组X的第DIM维的方向的元素求其平均值赋予向量Y若DIM=1，为按列操作；若DIM=2，为按行操作若X为二维数组，Y为一个向量；若X为一维数组，则Y为单变量n例：nx=[4 5 6;1 4 8];ny1= mean(x,1)ny1 =n 2.5000 4.5000 7.0000ny2= mean(x,2)ny2 =n 5.0000n 4.3333n6、求积、求积n命令格式：nY= prod(X)：将prod(X)返回矩阵X各列元素之积赋予行向量Y若X为向量，则Y为单变量nY= prod(X,DIM)：按数组X的第DIM维的方向的元素求其积赋予向量Y若DIM=1，为按列操作；若DIM=2，为按行操作若X为二维数组，Y为一个向量；若X为一维数组，则Y为单变量n例：nx=[4 5 6;1 4 8];ny1= prod(x,1)ny1 =n 4 20 48ny2= prod(x,2)ny2 =n 120n 327、、 求累计和、累积积、标准方差与升序排序求累计和、累积积、标准方差与升序排序n累计和函数CUMSUMn累积积CUMPROD n标准方差STD微积分n求极限nsyms x; %定义符号变量nf=(x*(exp(sin(x))+1)-2*(exp(tan(x))-1))/sin(x)^3; %确定符号表达式nw=limit(f) %求函数的极限解方程n解方程函数的格式：solve(expr1,expr2,...,exprN,var1,var2,...varN) n或 solve(expr1,expr2,...,exprN) n功能：求解代数方程组expr1,expr2,...,exprN的根，未知数为var1,var2,...varN。

n说明：若不指明符号表达式expr1,expr2,...,exprN的值，系统默认为0如给出一个表达式x^2-3*x-8,则系统将按x^2-3*x-8=0进行运算；n【例】解代数方程：a*x2-b*x-6=0nsyms a b xnsolve(a*x^2-b*x-6) nans =n [ 1/2/a*(b+(b^2+24*a)^(1/2))]n[ 1/2/a*(b-(b^2+24*a)^(1/2))]n即该方程有两个根： x1=1/2/a*(b+(b^2+24*a)^(1/2))；x2=1/2/a*(b-(b^2+24*a)^(1/2)) 例 研究矮壮素使玉米矮化的效果，在抽穗期测定喷矮壮素小区8株、对照区玉米9株，其观察值如下表：y y1 1( (喷施矮壮素喷施矮壮素) )160160160160200200160160200200170170150150210210y y2 2( (对照对照) )170170270270180180250250270270290290270270230230170170Matlab实现实现t检验检验 从理论上判断，喷施矮壮素只可能矮化无从理论上判断，喷施矮壮素只可能矮化无效而不可能促进植物长高，因此假设效而不可能促进植物长高，因此假设H H0 0：喷施：喷施矮壮素的株高与未喷的相同或更高，，即矮壮素的株高与未喷的相同或更高，，即H H0 0：： μμ1 1≥≥μμ2 2对对H HA A：： μμ1 1＜＜μμ2 2，即喷施矮壮素的株，即喷施矮壮素的株高较未喷的为矮。

显著水平高较未喷的为矮显著水平αα=0.05=0.05测验计算：测验计算： 按按νν=7+8=15=7+8=15，查，查t t 表得一尾表得一尾t t0.050.05=1.753(=1.753(一尾测验一尾测验t t0.050.05等于两尾测验的等于两尾测验的t t0.100.10) )，现实得，现实得t t=-3.05=-3.05＜＜- - t t0.050.05=-1.753=-1.753，，故故P P＜＜0.050.05推断：否定推断：否定H H0 0：： μμ1 1≥≥μμ2 2，接受，接受H HA A：： μμ1 1＜＜μμ2 2，即认为玉米喷施矮壮素后，其株，即认为玉米喷施矮壮素后，其株高显著地矮于对照高显著地矮于对照x=[160 160 200 160 200 170 150 210];y=[170 270 180 250 270 290 270 230 170];[h,p,ci,stats]=ttest2(x,y,0.05,‘left')n例例 选生长期、发育进度、植株大小和其它方面皆比较一致的两株番茄构成一组，共得7组，每组中一株接种A处理病毒，另一株接种B处理病毒，以研究不同处理方法的纯化的病毒效果。

组别组别y y1 1(A(A法法) )y y2 2(B(B法法) )d d1 110102525-15-152 2131312121 13 38 81414-6-64 43 31515-12-125 520202727-7-76 620202020 0 07 76 61818-12-12假设：两种处理对纯化病毒无不同效果，即：假设：两种处理对纯化病毒无不同效果，即：H H0 0：：μμd d=0=0 ；对；对H HA A：：μμd d≠0≠0显著水平显著水平αα=0.01=0.01测验计算：测验计算： 查附表查附表4 4，， νν=7-1=6=7-1=6时，时，t t0.010.01=3.707=3.707实得| |t t| |＞＞ t t0.01,0.01,故故P P＜＜0.010.01推断：否定推断：否定H H0 0：：μμd d=0=0，接受，接受H HA A：：μμd d≠0≠0，即，即A A、、B B两法对纯化病毒的效应有极显著差异两法对纯化病毒的效应有极显著差异y1=[10 13 8 3 20 20 6];y2=[25 12 14 15 27 20 18];[h,p,ci,stat]=ttest(y1,y2,0.05,'both')例例 以以A A、、B B、、C C、、D4D4种药剂处理水稻种子，种药剂处理水稻种子，其中其中A A为对照，每处理各得为对照，每处理各得4 4个苗高观察值个苗高观察值(cm)(cm)，试做方差分析。

试做方差分析药剂药剂苗高观察值苗高观察值总和总和T Ti i平均数平均数 A A181821212020131372721818B B202024242626222292922323C C101015151717141456561414D D28282727292932321161162929T T=336=336=21=21变异来源变异来源DFDFSSSSMSMSF F显著显著F F值值药剂处理间药剂处理间3 3504504168.00168.0020.520.56 6F F0.050.05=3.49=3.49药剂处理内药剂处理内( (误差误差) )121298988.178.17F F0.010.01=5.95=5.95总变异总变异1515602602水稻药剂处理苗高方差分析表水稻药剂处理苗高方差分析表x=[18 20 10 28; 21 24 15 27; 20 26 17 29; 13 22 14 32][p,anovatab,stats]=anova1(x)例例 采用采用5 5种生长素处理豌豆，未处理为对照，待种子种生长素处理豌豆，未处理为对照，待种子发芽后，分别每盆中移植发芽后，分别每盆中移植4 4株，每组株，每组6 6盆，每盆一个处盆，每盆一个处理，试验共有理，试验共有4 4组组2424盆，并按组排列于温室中，使同盆，并按组排列于温室中，使同组各盆的环境条件一致。

当各盆见第一朵花时记录组各盆的环境条件一致当各盆见第一朵花时记录4 4株豌豆的总节间数，结果见下表，试作方差分析株豌豆的总节间数，结果见下表，试作方差分析处理（处理（A A））组（组（B B））总和总和T Ti.i.平均平均I IIIIIIIIIIIIVIV未处理未处理(CK)(CK)606062626161606024324360.860.8赤霉素赤霉素656565656868656526326365.865.8动力精动力精636361616161606024524561.361.3吲哚乙酸吲哚乙酸硫酸腺嘌呤硫酸腺嘌呤马来酸马来酸64646262616167676565626263636262626261616464656525525525325325025063.863.863.363.362.562.5总和　总和　T.T.j j375375382382377377375375T T＝＝15091509变异来源变异来源DFDFSSSSMSMSF FF F0.050.05F F0.010.01组间组间3 35.455.451.821.82＜＜1 1处理间处理间5 565.8765.8713.1713.174.564.562.902.904.564.56误差误差151543.3043.302.892.89总变异总变异2323114.62114.62方差分析表方差分析表x=[60 62 61 60; 65 65 68 65; 63 61 61 60; 64 67 63 61; 62 65 62 64; 61 62 62 65][p,table,stats]=anova2(x)例例 一些夏季害虫盛发期的早迟和春季温度高低有关。

江苏武进连续9年测定3月下旬至4月中旬旬平均温度累积值(x)和水稻一代三化螟盛发期(y，以5月10日为0)的关系，得结果于表1试计算其直线回归方程x累积温累积温35.534.131.740.336.840.231.739.244.2y盛发期盛发期12169273139-1由观察值计算一级数据因而有：因而有：从而得到回归方程：从而得到回归方程：x=[35.5 34.1 31.7 40.3 36.8 40.2 31.7 39.2 44.2]’y=[12 16 9 2 7 3 13 9 -1]’[p,s]=polyfit(x,y,1)相关分析x=[35.5 34.1 31.7 40.3 36.8 40.2 31.7 39.2 44.2]’y=[12 16 9 2 7 3 13 9 -1]’[R,P,RLO,RUP]=CORRCOEF(x,y)。