《用蒙特卡洛方法估计积分方法及matlab编程实现》由会员分享，可在线阅读，更多相关《用蒙特卡洛方法估计积分方法及matlab编程实现（36页珍藏版）》请在金锄头文库上搜索。
1、用蒙特卡洛方法估计积分方法及mat lab编程实现专业班级：材料43学生姓名：王宏辉学 号：2140201060指导教师:李耀武 完成时间：2016年6月8日用蒙特卡洛方法估计积分方法及 matlab 编程实现实验内容：1用蒙特卡洛方法估计积分f xsinxdx， f e”dx和口 ex2+y2dxdy的值，00x2 + y 2 1并将估计值与真值进行比较。2用蒙特卡洛方法估计积分/ ex2dx和ff .1 dxdy的值，并对误1 + x 4 + y 40x2 + y2 0(a, b)，h(x)=疸)；令 Y=h(x)，则积分 S=E(Y)；利用 f(x)matlab 软件，编程产生随机变量
2、x 的随机数，在由y = h(x)I (x),I(x) T:爲，得到随机变量Y的随机数，求出样本均值，以此估计积分值。积分S =ff g(x,y)dxdy的求法与上述方法类似，在此不赘述。A概率密度函数的选取：一重积分，由于要求f(x)的支持域xl f (x) 0_(a, b)，为使方法普遍适用，考虑到标准正态分布概率密度函数f (x) = -e弋支持域为J iC 耳2兀R，故选用f (x)=三x2-e - 2。类似的，二重积分选用f (x,y)=丄e-亍，支持域为R2。 2兀估计评价：进行重复试验，通过计算样本均值以评价估计的无偏性；通过计 算均方误（针对第1类题，积得出）或样本方差（针对第
3、2类题，积 不出）以评价估计结果的精度。程序设计：依据问题分四类：第一类一重积分；第一类二重积分；第二类一 重积分，第二类二重积分，相应程序设计成四类。为了使程序具有一般性以及方便以后使用：一重积分，程序保存 为一个.m文本，被积函数，积分区间均采用键盘输入；二重积分， 程序主体保存为一个 .m 文本，被积函数键盘输入，示性函数用 function 语句构造，求不同区域二重积分，只需改变 function 函 数内容。编程完整解决用蒙特卡洛方法估计一重、二重积分值问题。程序代码及运行结果：第一类一重积分程序代码：%构造示性函数function I=I1（x,a,b）if x=a&x=bI=1;
4、elseI=0;%保存为Il.m%第一类一重积分,程序主体：%保存为fll.mfunction outf11=f11()g1=input(输入一元被积函数如x.*sin(x):,s)%输入被积 函数g1=inline(g1);a=inpu t(输入积分下界a:);%输入积分上下限b=input(输入积分上界b:);Real二input(积分真值:);%输入积分真值fprintf(输入样本容量 101102:丫) V=zeros(1,2);V二inpu t(Vl:);%输入样本容量V(2)=input(V2:);for m=V(1):V(2)% 样本容量 10F110F2n=10八mfor j=
5、1:10x=randn(1,n);for i=1:ntl(i)=Il(x(i),a,b);%示性及求和向量yl二gl(x)*(pi*2厂0.5).*exp(x.八2/2);Y1(j)=y1*t1/n; %单次实验样本均值 end t=ones(l,l0);EY=Yl*t/l0; %十次均值 D=abs(EY-Real); %绝对误差 RD=D/Real;%绝对误差d=0;for i=l:l0d二d+(Y1(i)-Real厂2;endd=d/(l0-l);EYl(m-V(l)+l)=EY;%样本容量为10 m时的样本均值Dl(m-V(l)+l)=D;%绝对误差RDl(m-V(l)+l)=RD;%
6、绝对误差MSEl(m-V(l)+l)=d;%方差endReal,EYl,Dl,RDl,MSEloutfll=EYl;Dl;RDl;MSEl; %存放样本数字特征%保存为fll.m运行结果：匹%估计积分fxsin xdx，积分真值为10m=f11输入一元被积函数如 x.*sin(x):x.*sin(x)g1 =x.*sin(x)输入积分下界 a:0输入积分上界 b:pi/2积分真值:1输入样本容量10AV1-10AV2:V1:1V2:5n =100100010000100000Real =EY1 =1.26351.00881.00661.01091.0018D1 =0.26350.00880.0
7、0660.01090.0018RD1 =0.26350.00880.00660.01090.0018MSE1 =0.64390.02050.00280.00060.0001m=1.26351.00881.00661.01091.00180.26350.00880.00660.01090.00180.26350.00880.00660.01090.00180.64390.02050.00280.00060.0001%估计积分Te x2dx真值为0.88620M=f11输入一元被积函数如x.*sin(x):exp(-x42)g1 =exp(-x.A2)输入积分下界 a:0输入积分上界 b:+inf
8、积分真值:piA0.5/2%0.8862输入样本容量 10AV1-10AV2:V1:1V2:410100100010000Real =0.8862EY1 =0.93330.90770.88730.8871D1 =0.04700.02150.00100.0009RD1 =0.05310.02430.00120.0010MSE1 =0.93330.90770.88730.88710.04700.02150.00100.00090.05310.02430.00120.00100.19270.01120.00160.0000第一类二重积分程序代码：%构造示性函数，求不同区域上积分只需更改示性函数 fu
9、nction I=I2(x,y)ifx八2+y八2=1I=1;elseI=0;end%保存为I2.m%第一类二重积分程序主体%保存为fl2.mfunction outf12=f12()g2=input(输入二元被积函数如 exp(x2+y.八2):,s)%输 入被积函数g2=inline(g2,x,y);Real二input(积分真值:);%输入积分真值fprintf(输入样本容量 101*101102*102:丫) V=zeros(1,2);V二inpu t(Vl:);%输入样本容量V(2)=input(V2:);for m=V(1):V(2)% 样本容量 10F110F2n=10 mx=r
10、andn(1,n);y=randn(1,n);for i=1:nt2(i)=I2(x(i),y(i);%示性及求和向量endy2二g2(x,y)*(2*pi).*exp(x2+y.八2)/2);Y2(j)=y2*t2/n; %单次实验样本均值endt=ones(1,10);EY=Y2*t/10;%十次均值D=abs(EY-Real); %绝对误差RD=D/Real;%绝对误差d=0;for i=1:10d二d+(Y2(i)-Real厂2;endd=d/(10-1);EY2(m-V(1)+1)=EY;%样本容量为10 m时的样本均值D2(m-V(1)+1)=D;%绝对误差RD2(m-V(1)+1
11、)=RD;%绝对误差MSE2(m-V(1)+1)=d;%方差end outf12=EY2;D2;RD2;MSE2; %存放样本数字特征%保存为fl2.m运行结果：%估计积分 口 ex2+y2dxdy，真值为 pi*(exp(l)-l)%5.3981x2 + y2 1m=f12输入二元被积函数如exp(x2+y.八2):exp(x2+y.八2)g2 =exp(x2+y.积分真值:pi*(exp(l)T)%5.3981输入样本容量 10X1*10X110X2*10X2:V1:1V2:4n =100100010000Real =5.3981EY2 =4.7702 5.1250 5.43175.4041D2 =0.6279 0.2732 0.03350.0060RD2 =0.1163MSE2 =0.05060.00623.89650.55640.0247m =4.77025.12505.43170.62790.27320.03350.11630.05060.00623.89650.55640.02470.00110.00175
