随机过程matlab程序Word格式文档下载.docx

1、plot（x,y2,k-）;legend（sinx,cosx）;2*pi;y=sin（x）;figure（1）plot（x,y,r-gridon以二元函数图z=xexp（-x2-y2）为例讲解基本操作，首先需要利用meshgrid函数生成X-Y平面的网格数据，如下所示：xa=-2:0.2:2;ya=xa;x,y=meshgrid（xa,ya）;z=x.*exp（-x.2-y.2）;mesh（x,y,z）;建立M文件functionfenshu（grade）ifgrade95.0disp（ThegradeisA.else86.0ThegradeisB.76.0ThegradeisC.66.0Th

2、egradeisD.ThegradeisF.endfunctiony=func（x）ifabs（x）1y=sqrt（1-x2）;elsey=x2-1;functionsumm（n）i=1;sum=0;while（i=n）sum=sum+i;i=i+1;str=?1?ao,num2str（sum）;disp（str）求极限symsxlimit（1+x）（1/x）,x,0,right求导数symsx;f=（sin（x）/x）;diff（f）diff（log（sin（x）求积分int（x2*log（x）int（abs（x-1）,0,2）常微分方程求解dsolve（Dy+2*x*y=x*exp（-x2

3、）,x计算偏导数x/（x2+y2+z2）（1/2）diff（x2+y2+z2）（1/2）,x,2）重积分int（int（x*y,y,2*x,x2+1）,x,0,1）级数symsn;symsum（1/2n,1,inf）Taylor展开式求y=exp（x）在x=0处的5阶Taylor展开式taylor（exp（x）,0,6）矩阵求逆A=0-6-1;62-16;-520-10det（A）inv（A）特征值、特征向量和特征多项式-520-10;lambda=eig（A）v,d=eig（A）poly（A）多项式的根与计算p=10-2-5;r=roots（p）p2=poly（r）y1=polyval（p,

4、4）例子：x=-3:3y=3.03,3.90,4.35,4.50,4.40,4.02,3.26;A=2*x,2*y,ones（size（x）;B=x.2+y.2;c=inv（A*A）*A*B;r=sqrt（c（3）+c（1）2+c（2）2）例子ezplot（-2/3*exp（-t）+5/3*exp（2*t）-2/3*exp（-t）+2/3*exp（2*t）,0,1）gridon;axis（0,12,0,5）密度函数和概率分布设xb（20,0.1）,binopdf（2,20,0.1）分布函数设xN（1100,502），yN（1150,802），则有normcdf（1000,1100,50）=0.

5、0228，1-0.0228=0.9772normcdf（1000,1150,80）=0.0304,1-0.0304=0.9696统计量数字特征x=29.827.628.3mean（x）max（x）min（x）std（x）symspk;Ex=symsum（k*p*（1-p）（k-1）,k,1,inf）symsxy;f=x+y;Ex=int（int（x*y*f,y,0,1）,0,1）参数估计例：对某型号的20辆汽车记录其5L汽油的行驶里程（公里），观测数据如下：29.827.628.327.930.128.729.928.027.928.728.427.229.528.528.030.029.12

6、9.829.626.9设行驶里程服从正态分布，试用最大似然估计法求总体的均值和方差。x1=29.827.628.327.930.128.729.928.027.928.7;x2=28.427.229.528.528.030.029.129.829.626.9;x=x1x2p=mle（norm,x）;muhatmle=p（1）,sigma2hatmle=p（2）2m,s,mci,sci=normfit（x,0.5）假设检验下面列出的是某工厂随机选取的20只零部件的装配时间（分）：9.810.410.69.69.79.910.911.19.610.210.39.69.911.210.69.810.

7、510.110.59.7设装配时间总体服从正态分布，标准差为0.4，是否认定装配时间的均值在0.05的水平下不小于10。解：在正态总体的方差已知时MATLAB均值检验程序：x1=9.810.410.69.69.79.910.911.19.610.2;x2=10.39.69.911.210.69.810.510.110.59.7;m=10;sigma=0.4;a=0.05;h,sig,muci=ztest（x,m,sigma,a,1）得到：h=1，%PPT例2一维正态密度与二维正态密度s=1;t=2;mu1=0;mu2=0;sigma1=sqrt（s2）;sigma2=sqrt（t2）;x=-6

8、:0.1:6;f1=1/sqrt（2*pi*sigma1）*exp（-（x-mu1）.2/（2*sigma12）;f2=1/sqrt（2*pi*sigma2）*exp（-（x-mu2）.2/（2*sigma22）;plot（x,f1,x,f2,k-.rho=（1+s*t）/（sigma1*sigma2）;f=1/（2*pi*sigma1*sigma2*sqrt（1-rho2）*exp（-1/（2*（1-rho2）*（x-mu1）2/sigma12-2*rho*（x-mu1）*（y-mu2）/（sigma1*sigma2）+（y-mu2）2/sigma22）;ezsurf（f）%P34例3.1.

9、1p1=poisscdf（5,10）p2=poisspdf（0,10）p1,p2%输出p1=0.0671p2=4.5400e-005ans=0.06710.0000%P40例3.2.1p3=poisspdf（9,12）p3=0.0874%P40例3.2.2p4=poisspdf（0,12）p4=6.1442e-006%P35-37（Th3.1.1）解微分方程%输入：symsp0p1p2;S=dsolve（Dp0=-lamda*p0Dp1=-lamda*p1+lamda*p0Dp2=-lamda*p2+lamda*p1,p0（0）=1p1（0）=0p2（0）=0S.p0,S.p1,S.p2%输出

10、：ans=exp（-lamda*t）,exp（-lamda*t）*t*lamda,1/2*exp（-lamda*t）*t2*lamda2%P40泊松过程仿真%simulate10timesclear;lamda=1;x=;fori=1:ms=exprnd（lamda,seed,1）;%seed是用来控制生成随机数的种子,使得生成随机数的个数是一样的.x=x,exprnd（lamda）;t1=cumsum（x）;x,t10.65090.65092.40613.05700.10023.15720.12293.28000.82334.10330.24634.34961.90746.25700.478

11、36.73531.34478.08000.80828.8882N=;fort=0:（t1（m）+1）iftt1（1）N=N,0;elseift=t1（i）&tt1（i+1）N=N,i;t1（m）N=N,m;%P48非齐次泊松过程仿真%exprnd（lamda,setseedsT=;T=T,t.3;%timeisadjusted,cumulativeintensityfunctionist3.plot（T,N,%output0.42200.42203.33233.75430.16353.91780.06833.98610.38754.37360.27744.65100.29694.94790.9

12、3595.88380.42246.30621.76508.071210timessimulation100timessimulation%P50复合泊松过程仿真niter=100;%iteratenumber%arrivingratet=input（Inputatime:sniterrand（state,sum（clock）;x=exprnd（lamda）;%intervaltimet1=x;whilet1tt1=sum（x）;%arrivingtimey=trnd（4,1,length（t1）;%rand（1,length（t1）;gamrnd（1,1/2,1,length（t1）;frnd

13、（2,10,1,length（t1）;t2=cumsum（y）;,y,t2X=;m=length（t1）;X=X,0;X=X,t2（i）;X=X,t2（m）;（t1（m）+1）,X,跳跃度服从0,1均匀分布情形跳跃度服从分布情形跳跃度服从t（10）分布情形%Simulatetheprobabilitythatsalesrevenuefallsinsomeinterval.（e.g.example3.3.6inteachingmaterial）niter=1.0E4;%numberofiterationslamda=6;%arrivingrate（unit:minute）t=720;%12hou

14、rs=720minutesabove=repmat（0,1,niter）;%setupstoragen=1;whilex=tn=n;n=n+1;z=binornd（200,0.5,1,n）;%generatensalesy=sum（z）;above（i）=sum（y432000）;pro=mean（above）Output:pro=0.3192%Simulatethelosspro.ForaCompoundPoissonprocessniter=1.0E3;below=repmat（0,1,niter）;x=x+exprnd（lamda）;r=normrnd（0.05/253,0.23/sqr

15、t（253）,1,n）;%generatenrandomjumpsy=log（1.0E6）+cumsum（r）;minX=min（y）;%minmumreturnovernextnjumpsbelow（i）=sum（minXlog（950000）;pro=mean（below）t=50,pro=0.45%P75（Example5.1.5）马氏链chushivec0=001000P=0,1/2,1/2,0,0,0;1/2,0,1/2,0,0,0;1/4,1/4,0,1/4,1/4,0;0,0,1,0,0,0,;0,0,1/2,0,0,1/2;0,0,0,0,1,0jueduivec1=chush

16、ivec0*Pjueduivec2=chushivec0*（P2）%计算1到n步后的分布chushivec0=001000;0,0,0,0,1,0;n=10t=1/6*ones（16）;jueduivec=repmat（t,n1）;fork=1:njueduiveck=chushivec0*（Pk）;jueduivec（k,1:6）=jueduiveck%比较相邻的两行n=70;jueduivecn=chushivec0*（Pn）n=71;%Replacethefirstdistribution,Comparingtwoneighbourabsolute-vectorsoncemorechus

17、hivec0=1/61/61/61/61/61/6;%赌博问题模拟（带吸收壁的随机游走：结束1次游走所花的时间及终止状态）a=5;p=1/2;m=0;whilem0&a10ifa=0t1=t1+m;m1=m1+1;t2=t2+m;m2=m2+1;fprintf（Theaveragetimesofarriving0and10respectivelyare%d,%d.n,t1/m1,t2/m2）;Thefrequenciesofarriving0and10respectivelyare%d,%d.n,m1/N,m2/N）;%verify:Theprobabilityofarriving0andit

18、sapproximaterespectivelyare%d,%d.n,5/10,m1/N）;Theexpectationofarriving0or10anditsapproximaterespectivelyare%d,%d.n,5*（10-5）/（2*p）,（t1+t2）/N）;%p=qp=1/4;t1=zeros（1,N）;t2=zeros（1,N）;15t1（1,n）=m;t2（1,n）=m;,sum（t1,2）/m1,sum（t2,2）/m2）;,（p10*（1-p）5-p5*（1-p）10）/（p5*（p10-（1-p）10）,m1/N）;,5/（1-2*p）-10/（1-2*p）*

19、（1-（1-p）5/p5）/（1-（1-p）10/p10）,（sum（t1,2）+sum（t2,2）/N）;%连续时间马尔可夫链通过Kolmogorov微分方程求转移概率输入：symsp00p01p10p11lamdamu;P=p00,p01;p10,p11;Q=-lamda,lamda;mu,-muP*Q输出：-p00*lamda+p01*mu,p00*lamda-p01*mu-p10*lamda+p11*mu,p10*lamda-p11*mup00,p01,p10,p11=dsolve（Dp00=-p00*lamda+p01*muDp01=p00*lamda-p01*muDp10=-p10*lamda+p11*muDp11=p10*lamda-p11*mup00（0）=1,p01（0）=0,p10（0）=0,p11（0）=1