用蒙特卡洛方法估计积分方法及matlab编程实现

1、精品用蒙特卡洛方法估计积分方法及matlab编程实现专业班级：材料43学生姓名：王宏辉学号：2140201060指导教师：李耀武完成时间：2016年6月8日感谢下载载用蒙特卡洛方法估计积分方法及matlab编程实现实验内容：i用蒙特卡洛方法估计积分2xsin xdxo2 一22e-xdx 和ex ydxdy 的值,0x2 y2 1并将估计值与真值进行比较2用蒙特卡洛方法估计积分12exdx和0差进行估计。要求：(1)针对要估计的积分选择适当的概率分布设计蒙特卡洛方法；(2)利用计算机产生所选分布的随机数以估计积分值；(3)进行重复试验，通过计算样本均值以评价估计的无偏性；通过计算均方误差(针对

2、第1类题)或样本方差(针对第2类题)以评价估计结果的精度目的：(1)能通过MATLAB或其他数学软件了解随机变量的概率密度、分布函数及其期望、方差、协方差等；(2)熟练使用MATLAB对样本进行基本统计，从而获取数据的基本信息;(3)能用MATLAB熟练进行样本的一元回归分析。实验原理:蒙特卡洛方法估计积分值，总的思想是将积分改写为某个随机变量的数学期望，借助相应的随机数，利用样本均值估计数学期望，从而估计相应的积分值。具体操作如下:bb()b，一般地，积分Sg(x)dx改与成Sf(x)dxh(x)f(x)dx的形aaf(x)a式，(其中为f(x)一随机变量X的概率密度函数，且f(x)的支持域

3、x|f(x)0(a,b),h(x)吗);4Y=h(X),则积分S=E(Y);利f(x)用matlab软件，编程产生随机变量X的随机数，在由yh(x)I(x),I(x)1,x(；：)，得到随机变量Y的随机数，求出样本均值，以此估计积分值。积分Sg(x,y)dxdy的求法与上述方法类似，在此不赘述。A概率密度函数的选取:0(a,b),为使方法普一重积分，由于要求f(x)的支持域x|f(x)遍适用，考虑到标准正态分布概率密度函数1与f (x) 2je e 支持域为1x2R，故选用f(x)7re2。22支持域为R2oxy类似的，二重积分选用f(x,y)2e2精品估计评价：进行重复试验，通过计算样本均值

4、以评价估计的无偏性；通过计算均方误（针对第1类题，积得出）或样本方差（针对第2类题，积不出）以评价估计结果的精度。程序设计：依据问题分四类：第一类一重积分；第一类二重积分；第二类一重积分，第二类二重积分，相应程序设计成四类。为了使程序具有一般性以及方便以后使用：一重积分，程序保存为一个.m文本，被积函数，积分区间均采用键盘输入；二重积分，程序主体保存为一个.m文本，被积函数键盘输入，示性函数用function语句构造，求不同区域二重积分，只需改变function函数内容。编程完整解决用蒙特卡洛方法估计一重、二重积分值问题。程序代码及运行结果：第一类一重积分程序代码：%构造示性函数functio

5、nI=I1（x,a,b）ifx=a&x二bI=1;elseend%保存为I1.m%第一类一重积分，程序主体：%保存为f11.mfunctionoutf11=f11()g1=input(输入一元被积函数如x.*sin(x):,s)%输入被积函数g1=inline(g1);a=input(输入积分下界a:);%输入积分上下限b=input(输入积分上界b:);Real=input(积分真值:);输入积分真值fprintf(输入样本容量1071-102:)V=zeros(1,2);V(1)=input(V1:);%输入样本容量V(2)=input(V2:);form!=丫(1)出样本容量10Am1-

6、10Am2n=10Amforj=1:10x=randn(1,n);fori=1:nt1(i)=I1(x(i),a,b);%示性及求和向量endy1=g1(x)*(pi*2)八0.5).*exp(x.A2/2);Y1(j)=y1*t1/n;%单次实验样本均值endt=ones(1,10);EY=Y1*t/10;%十次均值D=abs(EY-Real);%绝对误差RD=D/Real;%绝对误差d=0;fo门=1:10d=d+(Y1(i)-Real)A2;endd=d/(10-1);EY1(m-V(1)+1)=EY;%样本容量为10Am时的样本均D1(m-V(1)+1)=D;%绝对误差RD1(m-V(

7、1)+1)=RD;%绝对误差MSE1(m-V(1)+1)=d;%方差endReal,EY1,D1,RD1,MSE1outf11=EY1;D1;RD1;MSE1;%存放样本数字特征%保存为f11.m运行结果：2%估计积分xsinxdx,积分真值为10m=f11输入一元被积函数如x.*sin(x):x.*sin(x)g1=x.*sin(x)输入积分下界a:0输入积分上界b:pi/2积分真值：1输入样本容量10AV1-10AV2:V1:1V2:5n=10感谢下载载100100010000100000Real=EY1=1.26351.00881.00661.0109D1=0.26350.00880.0

8、0660.0109RD1=0.26350.00880.00660.0109MSE1=0.64390.02050.00280.00061.00180.00180.00180.0001m=1.01091.00181.26351.00881.00660.26350.00880.00660.01090.00180.26350.00880.00660.01090.00180.64390.02050.00280.00060.00012%估计积分e-xdx真值为0.88620M=f11输入一元被积函数如x.*sin(x):exp(-xJ2)g1=exp(-x.A2)输入积分下界a:0输入积分上界b:+inf

9、积分真值钟0.5/2%0.8862输入样本容量10AV1-10AV2:V1:1V2:4n=10n=100n=1000n=10000精品MSE1 =感谢下载载Real=0.8862EY1=0.93330.90770.88730.8871D1=0.04700.02150.00100.0009RD1=0.05310.02430.00120.0010精品0.19270.01120.00160.0000M=0.93330.90770.88730.88710.04700.02150.00100.00090.05310.02430.00120.00100.19270.01120.00160.0000感谢下载

10、载第一类二重积分程序代码：%构造示性函数，求不同区域上积分只需更改示性函数functionI=I2(x,y)ifxA2+yA2=a&xS=quadl(f,0,1)S=1.4627感谢下载载第二类二重积分程序代码：%构造示性函数，求不同区域上积分只需更改示性函数functionI=I2(x,y)ifxA2+yA2=1I=1;elseI=0;end%保存为I2.m%第二类二重积分函数主体%,程序保存为f22.mfunctionoutf22=f22()g2=input(输入二元被积函数如1./(1+乂八4+丫八4)八0.5:,&)%输入被积函数g2=inline(g2,x,y);fprintf(输入

11、样本容量10A丫1*10八丫1-10八丫2*10八丫2：)V=zeros(1,2);V(1)=input(V1:);%输入样本容量V(2)=input(V2:);form=V(1):V(2)%样本容量10Am1-10Am2n=10Amforj=1:10x=randn(1,n);y=randn(1,n);fori=1:nt2(i)=I2(x(i),y(i);%示性及求和向量endy2=g2(x,y)*(2*pi).*exp(x.A2+y.A2)/2);Y2(j)=y2*t2/n;%单次实验样本均值endt=ones(1,10);EY=Y2*t/10;%十次均值d=0;fo门=1:10d=d+(Y

12、2(i)-EY)A2;endd=d/(10-1);时的样本均值EY2(m-V(1)+1)=EY;%样本容量为10AmMSE2(m-V(1)+1)=d;%方差endEY2,MSE2outf22=EY2;MSE2;%存放样本数字特征%第二类二重积分，程序保存为f22.m运行结果：%估计积分,1=dxdyx2y2i,1x4y4m=f22输入二元被积函数如1./(1+x.A4+y.A4).A0.5:1./(1+x.A4+y.A4).A0.5g2=1./(1+x.A4+yA4)A0.5输入样本容量10A丫1*10八丫1-10八丫2*10八丫2:V1:1V2:4n=10n=精品100感谢下载载100010

13、000EY2=3.07592.96992.85662.8269MSE2=1.32670.09000.00600.0014精品3.07592.96991.32670.09002.85662.82690.00600.0014感谢下载载实验结果整理：第一类一重积分：2估计积分xsinxdx0积分真值：1积分估计值：1.0018样本容量：101001000100001000001.26351.00881.00661.01090.26350.00880.00660.01090.26350.00880.00660.01090.64390.02050.00280.00062e-xdx0.8862积分估计值：

14、0.8871101001000100000.93330.90770.88730.88710.04700.02150.00100.00090.05310.02430.00120.00100.19270.01120.00160.0000样本均值:绝对误差:相对误差:均方误差:估计积分0积分真值:样本容量:样本均值:绝对误差:相对误差:均方误差:1.00180.00180.00180.0001第一类二重积分:估计积分x222exydxdyy21积分真值:5.3981积分估计值:5.4041样本容量:10100100010000样本均值:4.77025.12505.43175.4041绝对误差:0.6

15、2790.27320.03350.0060相对误差:0.11630.05060.00620.0011均方误差：3.89650.55640.02470.0017100001.45900.0008100002.82690.0014第二类一重积分：1-2估计积分exdx0积分估计值：1.4590样本容量：101001000样本均值：2.07821.65831.5029样本方差：0.43150.08890.0057用matlab指令求得积分结果1.4627第二类二重积分:估计积分1-dxdyx2y21,1x4y4积分估计值：2.8269样本容量：101001000样本均值：3.07592.96992.

16、8566样本方差：1.32670.09000.0060实验结果分析:2从第一类积分看，以估计积分xsinxdx为例:0积分真值：1样本容量：10样本均值：1.2635绝对误差：0.2635相对误差：0.2635均方误差：0.6439积分估计值：1.001810010001.00881.00660.00880.00660.00880.00660.02050.0028100001000001.01091.00180.01090.00180.01090.00180.00060.0001随着样本容量的增大，样本均值有接近积分真值的趋势，绝对误差、相对误差、均方误差呈减小趋势；随着样本容量的增大，样本均

17、值有接近积分真值的趋势，说明估计具有无偏性；绝对误差、相对误差、均方误差呈减小趋势，说明增大样本容量能提高估计精度；验证了蒙特卡洛方法估计积分值的可行性，为后续估计第二类积分提供了1一、一一一2-从第二类积分看，以估计积分exdx为例:0积分估计值：1.4590样本容量:样本均值:样本方差:102.07820.43151001.65830.088910001.50290.0057100001.45900.0008用matlab指令求得积分结果1.4627由于积分真值未知，无法直接比较估计值与积分值值；但随样本容量增大，样本方差减小，间接反映了估计精度的提高。蒙特卡洛方法估计值1.4590相比用matlab指令求得的积分结果1.4627,绝对偏差0.0038,相对偏差0.0025。蒙特卡洛方法估计值与用matlab指令求得的积分结果相互验证。总结与讨论：蒙特卡洛方法是基于随机数的一种统计方法。蒙特卡洛方法估计积分值，总的思想是将积分改写为某个随机变量的数学期望，借助相应的随机数，利用样本均值估计数学期望，从而估计相应的积分值。为使方法具有一般性，概率密度函数一重积分选择了1x21xLif(x)-y=e2,二重积分选用f(x,y)2e2。程序设计方面，本着使程序具有一般性以及方便以后使用的原则，依据问题分四类：第一类一重积分；第一类二重积分；第二类一重积分，