蒙特卡洛方法解决威布尔密度分布函数

%P(X1.8)=1-P(X=1.8)%=1-P(0X=1.8)%画图观察一下图形a=2;k=3;x=0:0.1:3.6;%为了比较好的观察，设置到3.6f=(k/a).*((x/a).^(k-1)).*exp(-(x/a).^k);plot(x,f,'b')title('概率密度分布函数');gridon00.511.522.533.5400.10.20.30.40.50.60.7概率密度分布函数%函数积分的方法a=2;k=3;symsxfFx;f=(k/a).*((x/a).^(k-1)).*exp(-(x/a).^k);%当要求X1.8时，也就是da=int(f,1.8,inf)%最终答案0.4824da=1/exp(729/1000);%蒙特卡罗方法：随机试验的方法计算积分%方法1：%x范围(0,1.8),y的范围是(0,0.6)形成一个矩形%均匀布点N，计算落入曲线下面的数据点的个数acount%那么P(x=1.8)的(面积)概率也就是1.8*0.6*acount/N%当然，这个方法取决于布点的密度，也就是个数的多寡a=2;k=3;x=0:0.01:1.8;y=0:0.01:0.6;sx=size(x);sy=size(y);N=sx(1,2)*sy(1,2);%总共有N=11041个点acount=0;%计算落入曲线下方的点的个数fori=1:sx(1,2)forj=1:sy(1,2)t=(k/a).*((x(i)/a).^(k-1)).*exp(-(x(i)/a).^k);ify(j)t%在曲线下方acount=acount+1;endendendp=1.8*0.6*acount/N;PP=1-p;%跟实际答案的误差是：6.2681e-004wc=abs(PP-da);%方法2：采用随机布点10000个的方法a=2;k=3;N=10000;n=100;acount=0;%计算落入曲线下方的点的个数fori=1:nrx=rand*1.8;forj=1:nry=rand*0.6;t=(k/a).*((rx/a).^(k-1)).*exp(-(rx/a).^k);ifryt%在曲线下方acount=acount+1;endendendp=1.8*0.6*acount/N;PP=1-p;%跟实际答案的误差是：0.0160wc=abs(PP-da);%方法3：随机数是落在曲线上面函数的值%采用平均值估算最后的积分%试验的次数，a=2;k=3;x1=0;x2=1.8;l=x2-x1;%在（0,1.8）之间随机选取N个点，计算函数值，最后这些平均值再乘以1.8N=10000;x=unifrnd(0,1.8,N,1);p=1.8.*mean((k/a).*((x/a).^(k-1)).*exp(-(x/a).^k));PP=1-p;%跟实际答案的误差是：0.0028wc=abs(PP-da);点数累积概率概率误差方法1110410.481766.2681e-004方法2100000.42470.0160方法3100000.48520.0028