用蒙特卡洛方法估计积分方法及matlab编程实现
matlab用蒙特卡洛方法进行概率和分位计算
matlab用蒙特卡洛方法进行概率和分位计算【主题】matlab用蒙特卡洛方法进行概率和分位计算【序号1】引言在概率和统计领域，计算概率和分位数一直是一个重要的课题。
传统的方法可能在计算复杂的分布时显得力不从心，而蒙特卡洛方法却能够以随机模拟的方式来解决这些问题。
本文将介绍如何使用MATLAB来进行概率和分位计算，重点讨论如何利用蒙特卡洛方法来进行模拟，以及如何在MATLAB环境中实现这一过程。
【序号2】MATLAB中的蒙特卡洛方法MATLAB作为一个强大的数值计算工具，提供了丰富的函数和工具箱，可以方便地进行概率和统计计算。
在MATLAB中，蒙特卡洛方法可以通过随机数生成函数和循环结构来实现。
我们需要生成符合指定分布的随机数样本，然后利用这些样本进行模拟计算，最终得到所需的概率和分位数结果。
【序号3】随机数生成在MATLAB中，可以利用内置的随机数生成函数来生成符合某个特定分布的随机数样本。
可以使用randn函数来生成符合正态分布的随机数样本，使用rand函数来生成在[0,1]区间均匀分布的随机数样本。
除了内置函数，MATLAB还提供了更多灵活的工具箱，可以生成更加复杂的分布样本，如指数分布、泊松分布等。
【序号4】模拟计算一旦得到了符合特定分布的随机数样本，就可以利用这些样本进行模拟计算。
以正态分布为例，我们可以利用蒙特卡洛方法来估计在某个区间内的概率，或者计算某个分位数的取值。
通过多次模拟，取平均值可以得到一个较为准确的估计结果。
在MATLAB中，可以利用循环结构和向量化的方式来高效地实现这一过程，并得到稳健可靠的结果。
【序号5】具体案例下面通过一个具体案例来展示如何在MATLAB中使用蒙特卡洛方法进行概率和分位计算。
假设我们需要计算标准正态分布的概率P(-1<Z<1)和95%的分位数。
我们可以利用randn函数生成一组标准正态分布的随机数样本，然后利用循环结构来进行模拟计算。
我们得到了P(-1<Z<1)约等于0.6827和95%的分位数约等于1.645，这些结果可以帮助我们更好地理解正态分布的性质。
今天刚好在研究蒙特卡洛方法计算定积分
今天刚好在研究蒙特卡洛方法计算定积分，下面是用随机点方法做的，希望对你有所帮助：%%%%%%%Matlab程序%%%%%方法一%%%%%%%%%%%%%clc;clear;n=10^4; %% 计算精度count=0; %% count为频数a=0; %% 积分下限b=1; %% 积分上限k=unifrnd(a,b,[2 n]); %% 积分区间内产生2n个随机数（均匀分布）A=k(1,:);B=k(2,:);for (i=1:n)x(i)=A(i);y(i)=B(i);end%%开始计算函数for i=1:nif y(i)<=fun(x(i)) %%ma函数，可定义不同的function计算不同的积分count=count+1;endendj=count/n%%%%%%%%计算结果为%%%%%%%%%j=0.743600000%%%%%%%%Matlab程序%%%%%方法二%%%%%%%%%%%%MC随机点法clear;clc;format long;a=0 ;b=1 ;n=10^2 ;k1=unifrnd(a,b,[1 n]);for i=1:nc1(i)=fun(k1(i));endc=min(c1-c1(n/2));d=max(c1+c1(n/2));k2=unifrnd(c,d,[1 n]);count=0;for i=1:nx(i)=k1(i);y(i)=k2(i);if 0<y(i)&0<fun(x(i))&y(i)<fun(x(i))count=count+1;else if fun(x(i))<y(i)&fun(x(i))<0&y(i)<0count=count+1;endendendu=(count/n)*((b-a)*(d-c))%%%%%%%%%%%%计算结果为%%%%%%%%%%%%u=0.75868947评价：1.方法二适合几乎所有的不确定积分面积上下限的定积分，精确性更好。
基于matlab环境下蒙特卡罗法的实现
基于Matlab 环境下蒙特卡罗法的实现建筑与土木工程2011级 201121022 温秋平针对应用蒙特卡罗对连续型分布采取直接抽样法解决结构可靠度所遇到的困难，提出利用MATLAB 其强大数值计算功能来解决此类问题。
利用MATLAB 进行蒙特卡罗抽样模拟，在一定程度上减少了对连续型分布采用直接抽样时的困难，大大提高了计算效率。
1.蒙特卡罗法蒙特卡洛方法是以数理统计原理为基础的，又称随机模拟方法，是随着电子计算机的发展而逐步发展起不来的一种独特的数值方法。
用蒙特卡洛方法来研究事件的随机性是结构可靠度分析的一个重要方面。
蒙特卡洛方法的优点是，它回避了结构可靠度分析中的数学困难，不需要考虑结构极限状态曲面的复杂性，只需要得到结构的响应即可；缺点是计算虽大，因此目前还不作为一种常规的结构可靠度分析的方法来使用，只适用于一些情况复杂的结构，由于其具有相对较高的精度，常用于结构可靠度各种近似方法计算精度的检验和计算结果的校核。
直接抽样方法是蒙特卡洛分析最基本的一种方法，对于基本随机变量12(,,,)n X X X X =，其概率密度函数为()f x ，对应结构某一状态的功能函数为()Z g x =。
将随机样本值序列X 代入功能函数()Z g x =，若Z<0，则模拟的结构失效一次。
若总的模拟数为N ，功能函数Z<0的次数为f n ，则结构失效概率f P 的估计值ˆfP 为： ˆf fn P N= (1.1) 由伯努利大数定理：lim ()1f f N nP P Nε→∞-<= (1.2) 可得ˆfP 以概率收敛于f P 。
失效概率的同样可以表达为：[()]()f P I g x f x dx +∞-∞=⎰(1.3)其中[()]I g x 为()g x 的示性函数，即：1 ()0[()]0 ()0g x I g x g x <⎧=⎨≥⎩ (1.4)则结构失效概率f P 的估计值ˆf P 为：11ˆ[()]Nffii n P I g x NN===∑ (1.5)对于结构可靠度问题，其对应的结构失效概率的数量级通常为371010--。
matlab蒙特卡洛法求定积分
文章标题：探索matlab中的蒙特卡洛法求定积分在数学和计算科学中，求解定积分是一个常见的问题。
传统的数值积分方法中，蒙特卡洛法是一种非常有趣和强大的方法，能够对一些特殊的不易求解的定积分问题提供解决方案。
而在matlab这一强大的数学计算软件中，蒙特卡洛法同样有着广泛的应用。
1. 什么是蒙特卡洛法？蒙特卡洛法是一种基于随机采样的数值积分方法，其核心思想是利用随机抽样的方法逼近定积分的值。
具体来说，对于给定的函数$f(x)$以及区间$[a, b]$，蒙特卡洛法通过对函数在该区间上进行随机采样，并利用采样点的平均值来逼近定积分的值。
2. 在matlab中应用蒙特卡洛法在matlab中，可以利用蒙特卡洛法求解定积分问题。
通过生成服从均匀分布的随机数，并代入原函数，然后求解采样点的平均值，可以得到定积分的近似值。
matlab内置了丰富的数学计算和随机数生成函数，能够方便地实现蒙特卡洛法的计算。
3. 实例分析：使用matlab进行蒙特卡洛法求解定积分假设我们要求解函数$f(x)=x^2$在区间$[0, 1]$上的定积分，即$$\int_{0}^{1} x^2 \, dx$$我们可以在matlab中编写如下代码：```matlabN = 1000000; % 设定采样点的个数X = rand(1, N); % 生成均匀分布的随机数Y = X.^2; % 代入原函数integral_value = mean(Y); % 求解采样点的平均值```通过上述代码，我们得到了定积分的近似值integral_value。
在这个例子中，我们利用蒙特卡洛法求得了定积分的近似值。
4. 总结与展望通过本文的介绍，我们对matlab中蒙特卡洛法求解定积分的方法有了初步的了解。
蒙特卡洛法作为一种基于随机采样的数值积分方法，在matlab中有着广泛的应用。
在实际应用中，我们可以根据定积分的具体问题来灵活选择采样点的个数，并结合matlab强大的数学计算能力，在求解定积分问题中取得更加准确的结果。
蒙特卡洛方法 matlab代码
蒙特卡洛方法 matlab代码蒙特卡洛方法是一种随机模拟的数值计算方法，广泛应用于金融、物理、计算机科学、统计学等领域。
在这里，我们将介绍如何用matlab实现蒙特卡洛方法。
本文主要内容包括：蒙特卡洛方法的基本原理、常见应用、matlab代码实现、实例应用等。
一、蒙特卡洛方法基本原理蒙特卡洛方法是一种基于统计学的数值计算方法，其基本原理是使用随机数模拟复杂系统的行为，从而获得数值上的解决方案。
它的核心思想是，通过大量的重复实验来模拟随机过程，最终得到与实际结果相似的概率性解决方案。
蒙特卡洛方法有许多不同的应用，例如解决随机过程和数学物理问题、评估金融和投资风险等。
二、蒙特卡洛方法常见应用蒙特卡洛方法在许多领域都有广泛应用。
以下是一些蒙特卡洛方法的应用：1. 金融和投资风险评估：通过模拟资产价格的随机行为，可以估计资产组合的波动性和风险。
2. 物理建模：用来计算复杂系统的行为，例如氢气中原子之间的相互作用。
3. 工程设计：用于模拟复杂系统的行为，例如机械系统的振动行为。
4. 统计学：用于估算总体参数的不确定性和置信区间。
蒙特卡洛方法的实现过程包括以下几个步骤:1. 选择模型：选择适合模型，其模型应足够灵活，可以处理不同类型的数据和数据格式。
2. 生成随机数：应生成适量的随机数。
这些随机数可以是具有不同分布的数值，例如正态分布、均匀分布等。
3. 执行计算：为获得数值上的解决方案，应对随机数进行计算。
1. 步骤1：选择模型例如，我们要计算正态分布的均值。
这是一种非常基本的蒙特卡洛问题。
2. 步骤2：生成随机数我们可以用matlab生成伪随机数，例如normal随机数：n = 10000;data = normrnd(0, 1, n, 1);在这里，我们将生成10000个均值为0，标准差为1的正态分布随机数。
3. 步骤3：执行计算为了计算正态分布的均值，我们可以使用matlab的mean函数：ans = mean(data)输出将是这些随机数的平均值。
Matlab笔记——蒙特卡罗方法018
18. 蒙特卡罗方法（一）概述一、原理蒙特卡罗（Monte Carlo ）方法，是一种基于“随机数”的计算机随机模拟方法，通过大量随机试验，利用概率统计理论解决问题的一种数值方法。
其理论依据是：大数定律、中心极限定理（估计误差）。
常用来解决如下问题：1. 求某个事件的概率，基于“频率的极限是概率”；2. 可以描述为“随机变量的函数的数学期望”的问题，用随机变量若干个具体观察值的函数的算术平均值代替。
一般形式为求积分：[]()()()d ba E f X f x x x ϕ=⎰ X 为自变量（随机变量），定义域为[a,b], f (x )为被积函数，()x ϕ为概率密度函数。
蒙特卡罗法步骤为：(1) 依据概率分布()x ϕ不断生成N 个随机数x , 依次记为x 1, …, x N , 并计算f (x i );(2) 用f (x i )的算术平均值1()Nii f x N =∑近似估计上述积分类比：“积分”同“求和”，“d x ”同“1/N ”，“()x ϕ”同“服从()x ϕ分布的随机数”；(3) 停止条件：至足够大的N 停止；或者误差小于某值停止。
3. 利用随机数模拟各种分布的随机现象，进而解决实际问题。
二、优缺点优点：能够比较逼真地描述具有随机性质的事物的特点及物理实验过程；受几何条件限制小；收敛速度与问题的维数无关；误差容易确定。
缺点：收敛速度慢；误差具有概率性；进行模拟的前提是各输入变量是相互独立的。
三、应用随机模拟实验，随机最优化问题，含有大量不确定因素的复杂决策系统进行风险模拟分析（金融产品定价、期权）。
（二）用蒙特卡罗法求事件概率一、著名的“三门问题”源自博弈论的一个数学游戏：参赛者面前有三扇关闭的门，其中一扇门的后面藏有一辆汽车，而两扇门的后面各藏有一只山羊。
参赛者从三扇门中随机选取一扇，若选中后面有车的那扇门就可以赢得该汽车。
当参赛者选定了一扇门，但尚未开启它的时候，节目主持人会从剩下的两扇门中打开一扇藏有山羊的门，然后问参赛者要不要更换自己的选择，选取另一扇仍然关着的门。
概率实验报告_蒙特卡洛积分
本科实验报告实验名称：《概率与统计》随机模拟实验随机模拟实验实验一设随机变量X 的分布律为-i P{X=i}=2,i=1,2,3......试产生该分部的随机数1000个，并作出频率直方图。
一、实验原理采用直接抽样法：定理：设U 是服从[0,1]上的均匀分布的随机变量，则随机变量-1()Y F U =与X 有相同的分布函数-1()Y F U =（为F(x)的逆函数），即-1()Y F U =的分部函数为()F x .二、题目分析易得题中X 的分布函数为1()1- ,1,0,1,2,3, (2i)F x i x i i =≤≤+=若用ceil 表示对小数向正无穷方向取整，则F(x)的反函数为产生服从[0,1]上的均匀分布的随机变量a ，则m=F -1(a)则为题中需要产生的随 机数。
三、MATLAB 实现f=[]; i=1;while i<=1000a=unifrnd(0,1); %产生随机数a ，服从【0,1】上的均匀分布 m=log(1-a)/log(1/2);b=ceil(m); %对m 向正无穷取整 f=[f,b]; i=i+1; enddisplay(f);[n,xout]=hist(f); bar(xout,n/1000,1)产生的随机数（取1000个中的20个）如下：-1ln(1-)()1ln()2a F a ceil ⎡⎤⎢⎥=⎢⎥⎢⎥⎣⎦频率分布直方图实验二设随机变量X 的密度函数为24,0,()0,0x xe x f x x -⎧>=⎨≤⎩试产生该分布的随机数1000个，并作出频率直方图 一、实验原理取舍抽样方法，当分布函数的逆函数难以求出时，可采用此方法。
取舍抽样算法的流程为：（1） 选取一个参考分布，其选取原则，一是该分布的随机样本容易产生；二是存在常数C ，使得()()f x Cg x ≤。
（2） 产生参考分布()g x 的随机样本0x ； （3） 独立产生[0,1]上的均匀分布随机数0u ；（4） 若000()()u Cg x f x ≤，则保留x 0,作为所需的随机样本；否则舍弃。
用蒙特卡洛方法估计积分值
用蒙特卡洛方法估计积分值一.实验目的：1.初步了解蒙特卡洛算法及其用途2.利用蒙特卡洛算法计算积分值并与其真实之进行比较二.实验原理：做Monte Carlo 时，求解积分的一般形式是：X 为自变量，它应该是随机的，定义域为(x0, x1)，f(x)为被积函数，ψ(x)是x 的概率密度。
Monte Carlo 步骤：1．依据概率分布ψ(x)不断生成随机数x, 并计算f(x)由于随机数性质，每次生成的x 的值都是不确定的，为区分起见，我们可以给生成的x 赋予下标。
如x i 表示生成的第i 个x 。
生成了多少个x ，就可以计算出多少个f(x)的值2．误差分析Monte Carlo 方法得到的结果是随机变量，因此，在给出点估计后，还需要给出此估计值的波动程度及区间估计。
严格的误差分析首先要从证明收敛性出发，再计算理论方差，最后用样本方差来替代理论方差。
三．实验内容：第一题 估计积分，并将估计值与真值进行比较(1).∫322dx x将被积函数展开成[2,3]上的均匀分布，1)(=x f x 2)(x x h =真值为6.3333程序内容运算结果（2）xdx x sin 20∫π将被积函数展开为[0,2π]上的均匀分布，π2)(=x f x ，x x x h sin 2)(π=真值为1程序内容运算结果（3）∫+∞−02dxe x将被积函数展开为（0, ∞+）上参数为2的卡方分布，f(x)=,ℎ()=2程序内容运算结果：第二题估计积分，并对误差进行估计（1）.∫将被积函数展开为(0,1)均匀分布，f(x)=1,ℎ()=程序内容运算结果(2)∫√将被积函数展开为(0,10)上的均匀分布，f(x)=0.1,ℎ()=√ 程序内容四.实验总结：通过实验了解了概率论在值估计及生活中的运用。
通过对上面五个问题的求解知道，由于随机数的任意性，虽然计算机每次的运行结果都是不一样的，但是结果往往与理论值偏差不大。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
用蒙特卡洛方法估计积分方法及matlab编程实现专业班级：材料43学生姓名：王宏辉学号：2140201060指导教师：李耀武完成时间：2016年6月8日用蒙特卡洛方法估计积分 方法及matlab 编程实现实验内容：1用蒙特卡洛方法估计积分 20sin x xdx π⎰，2-0x e dx +∞⎰和22221xy x y e dxdy ++≤⎰⎰的值，并将估计值与真值进行比较。
2用蒙特卡洛方法估计积分 210x e dx ⎰和22x y +≤⎰⎰的值，并对误差进行估计。
要求：（1）针对要估计的积分选择适当的概率分布设计蒙特卡洛方法；（2）利用计算机产生所选分布的随机数以估计积分值； （3）进行重复试验，通过计算样本均值以评价估计的无偏性；通过计算均方误差（针对第1类题）或样本方差（针对第2类题）以评价估计结果的精度。
目的：（1）能通过 MATLAB 或其他数学软件了解随机变量的概率密度、分布函数及其期望、方差、协方差等；（2） 熟练使用 MATLAB 对样本进行基本统计，从而获取数据的基本信息；（3） 能用 MATLAB 熟练进行样本的一元回归分析。
实验原理：蒙特卡洛方法估计积分值，总的思想是将积分改写为某个随机变量的数学期望，借助相应的随机数，利用样本均值估计数学期望，从而估计相应的积分值。
具体操作如下：一般地，积分⎰=bdx x g a )(S 改写成⎰⎰==bb dx f h dx f g aa )(x )(x )(x f(x))(x S 的形式，（其中为)f(x 一随机变量X 的概率密度函数，且)f(x 的支持域）（｝｛b f ,a 0)(x |x ⊇>)，f(x ))(x )(x g h =）;令Y=h(X)，则积分S=E （Y ）；利用matlab 软件，编程产生随机变量X 的随机数，在由⎩⎨⎧∉∈==)b (a,,)b (a,,01I(x) ,)(x )(x y x x I h ，得到随机变量Y 的随机数，求出样本均值，以此估计积分值。
积分⎰⎰=Adxdy g S )y (x,的求法与上述方法类似，在此不赘述。
概率密度函数的选取：一重积分，由于要求)f(x 的支持域）（｝｛b f ,a 0)(x |x ⊇>，为使方法普遍适用，考虑到标准正态分布概率密度函数22e 21)(x xf -=π支持域为R ，故选用22e 21)(x x f -=π。
类似的，二重积分选用22221)y (x,y x e f +-=π，支持域为2R 。
估计评价：进行重复试验，通过计算样本均值以评价估计的无偏性；通过计算均方误（针对第1类题，积得出）或样本方差（针对第2类题，积不出）以评价估计结果的精度。
程序设计：依据问题分四类：第一类一重积分；第一类二重积分；第二类一重积分，第二类二重积分，相应程序设计成四类。
为了使程序具有一般性以及方便以后使用：一重积分，程序保存为一个.m文本，被积函数，积分区间均采用键盘输入；二重积分，程序主体保存为一个.m文本，被积函数键盘输入，示性函数用function 语句构造，求不同区域二重积分，只需改变function 函数内容。
编程完整解决用蒙特卡洛方法估计一重、二重积分值问题。
程序代码及运行结果：第一类一重积分程序代码：%%%构造示性函数function I=I1(x,a,b)if x>=a&&x<=bI=1;elseI=0;%保存为I1.m%%%%%%%%%%%%%%%%%%第一类一重积分,程序主体：%保存为f11.mfunction outf11=f11()g1=input('输入一元被积函数如x.*sin(x):','s')%输入被积函数g1=inline(g1);a=input('输入积分下界a:');%输入积分上下限b=input('输入积分上界b:');Real=input('积分真值:');%输入积分真值fprintf('输入样本容量10^V1--10^V2:\r')V=zeros(1,2);V(1)=input('V1:');%输入样本容量V(2)=input('V2:');for m=V(1):V(2)%样本容量10^m1--10^m2n=10^mfor j=1:10x=randn(1,n);for i=1:nt1(i)=I1(x(i),a,b);%示性及求和向量y1=g1(x)*((pi*2)^0.5).*exp(x.^2/2);Y1(j)=y1*t1'/n; %单次实验样本均值endt=ones(1,10);EY=Y1*t'/10; %十次均值D=abs(EY-Real); %绝对误差RD=D/Real; %绝对误差d=0;for i=1:10d=d+(Y1(i)-Real)^2;endd=d/(10-1);EY1(m-V(1)+1)=EY; %样本容量为10^m时的样本均值D1(m-V(1)+1)=D; %绝对误差RD1(m-V(1)+1)=RD; %绝对误差MSE1(m-V(1)+1)=d; %方差endReal,EY1,D1,RD1,MSE1outf11=[EY1;D1;RD1;MSE1]; %存放样本数字特征%保存为f11.m运行结果：%估计积分2sinx xdxπ⎰，积分真值为1 m=f11输入一元被积函数如x.*sin(x):x.*sin(x) g1 =x.*sin(x)输入积分下界a:0输入积分上界b:pi/2积分真值:1输入样本容量10^V1--10^V2:V1:1V2:5n =10n =100n =1000 n =10000 n =100000 Real =1EY1 =1.2635 1.0088 1.0066 1.0109 1.0018 D1 =0.2635 0.0088 0.0066 0.0109 0.0018 RD1 =0.2635 0.0088 0.0066 0.0109 0.0018 MSE1 =0.6439 0.0205 0.0028 0.0006 0.0001m=1.2635 1.0088 1.0066 1.0109 1.00180.2635 0.0088 0.0066 0.0109 0.00180.2635 0.0088 0.0066 0.0109 0.00180.6439 0.0205 0.0028 0.0006 0.0001%估计积分 2-0x edx +∞⎰真值为0.8862M=f11输入一元被积函数如x.*sin(x):exp(-x.^2) g1 =exp(-x.^2)输入积分下界a:0 输入积分上界b:+inf 积分真值:pi^0.5/2%0.8862 输入样本容量 10^V1--10^V2:V1:1V2:4n =10n =100n =1000 n =10000Real =0.8862EY1 =0.9333 0.9077 0.8873 0.8871 D1 =0.0470 0.0215 0.0010 0.0009 RD1 =0.0531 0.0243 0.0012 0.0010 MSE1 =0.1927 0.0112 0.0016 0.0000 M =0.9333 0.9077 0.8873 0.88710.0470 0.0215 0.0010 0.00090.0531 0.0243 0.0012 0.00100.1927 0.0112 0.0016 0.0000第一类二重积分程序代码：%%%构造示性函数，求不同区域上积分只需更改示性函数function I=I2(x,y)if x^2+y^2<=1I=1;elseI=0;end%保存为I2.m%第一类二重积分程序主体%保存为f12.mfunction outf12=f12()g2=input('输入二元被积函数如exp(x.^2+y.^2):','s')%输入被积函数g2=inline(g2,'x','y');Real=input('积分真值:');%输入积分真值fprintf('输入样本容量10^V1*10^V1--10^V2*10^V2:\r') V=zeros(1,2);V(1)=input('V1:');%输入样本容量V(2)=input('V2:');for m=V(1):V(2)%样本容量10^m1--10^m2n=10^mfor j=1:10x=randn(1,n);y=randn(1,n);for i=1:nt2(i)=I2(x(i),y(i));%示性及求和向量endy2=g2(x,y)*(2*pi).*exp((x.^2+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)=RD; %绝对误差MSE2(m-V(1)+1)=d; %方差 endReal,EY2,D2,RD2,MSE2outf12=[EY2;D2;RD2;MSE2]; %存放样本数字特征 %保存为f12.m 运行结果： %估计积分 22221xy x y e dxdy ++≤⎰⎰，真值为pi*(exp(1)-1)%5.3981m=f12输入二元被积函数如exp(x.^2+y.^2):exp(x.^2+y.^2) g2 =exp(x.^2+y.^2)积分真值:pi*(exp(1)-1)%5.3981输入样本容量 10^V1*10^V1--10^V2*10^V2: V1:1 V2:4 n =10n =100n =1000 n =10000 Real =5.3981EY2 =4.77025.1250 5.4317 5.4041 D2 =0.6279 0.2732 0.0335 0.0060RD2 =0.1163 0.0506 0.0062 0.0011 MSE2 =3.8965 0.5564 0.0247 0.0017m =4.77025.1250 5.4317 5.40410.6279 0.2732 0.0335 0.00600.1163 0.0506 0.0062 0.00113.8965 0.5564 0.0247 0.0017第二类一重积分程序代码：%%%构造示性函数function I=I1(x,a,b)if x>=a&&x<=bI=1;elseI=0;end%保存为I1.m%第二类一重积分程序主体%程序保存为f21.mfunction outf21=f21()g1=input('输入一元被积函数如exp(x.^2):','s')%输入被积函数g1=inline(g1);a=input('输入积分下界a:');%输入积分上下限b=input('输入积分上界b:');fprintf('输入样本容量10^V1--10^V2:\r')V=zeros(1,2);V(1)=input('V1:');%输入样本容量V(2)=input('V2:');for m=V(1):V(2)%样本容量10^m1--10^m2n=10^mfor j=1:10x=randn(1,n);for i=1:nt1(i)=I1(x(i),a,b);%示性及求和向量endy1=g1(x)*((pi*2)^0.5).*exp(x.^2/2);Y1(j)=y1*t1'/n; %单次实验样本均值endt=ones(1,10);EY=Y1*t'/10; %十次均值d=0;for i=1:10d=d+(Y1(i)-EY)^2;endd=d/(10-1);EY1(m-V(1)+1)=EY; %样本容量为10^m时的样本均值MSE1(m-V(1)+1)=d; %方差endEY1,MSE1outf21=[EY1;MSE1]; %存放样本数字特征%%%%程序保存为f21.m 运行结果：%估计积分21xe dxm=f21输入一元被积函数如exp(x.^2):exp(x.^2) g1 =exp(x.^2)输入积分下界a:0输入积分上界b:1输入样本容量10^V1--10^V2:V1:1V2:4n =10n =100n =1000n =10000EY1 =2.0782 1.6583 1.5029 1.4590 MSE1 =0.4315 0.0889 0.0057 0.0008m =2.0782 1.6583 1.5029 1.45900.4315 0.0889 0.0057 0.0008%用matlab 指令求积分21xe dxf=inline('exp(x.^2)')f =Inline function:f(x) = exp(x.^2)>> S=quadl(f,0,1)S =1.4627第二类二重积分程序代码：%%%构造示性函数，求不同区域上积分只需更改示性函数function I=I2(x,y)if x^2+y^2<=1I=1;elseI=0;end%保存为I2.m%第二类二重积分函数主体%，程序保存为f22.mfunction outf22=f22()g2=input('输入二元被积函数如1./(1+x.^4+y.^4).^0.5:','s')%输入被积函数g2=inline(g2,'x','y');fprintf('输入样本容量10^V1*10^V1--10^V2*10^V2:\r')V=zeros(1,2);V(1)=input('V1:');%输入样本容量V(2)=input('V2:');for m=V(1):V(2)%样本容量10^m1--10^m2n=10^mfor j=1:10x=randn(1,n);y=randn(1,n);for i=1:nt2(i)=I2(x(i),y(i));%示性及求和向量endy2=g2(x,y)*(2*pi).*exp((x.^2+y.^2)/2);Y2(j)=y2*t2'/n; %单次实验样本均值endt=ones(1,10);EY=Y2*t'/10; %十次均值d=0;for i=1:10d=d+(Y2(i)-EY)^2;endd=d/(10-1);EY2(m-V(1)+1)=EY; %样本容量为10^m时的样本均值MSE2(m-V(1)+1)=d; %方差endEY2,MSE2outf22=[EY2;MSE2]; %存放样本数字特征%第二类二重积分，程序保存为f22.m运行结果：%估计积分⎰⎰22+≤x ym=f22输入二元被积函数如1./(1+x.^4+y.^4).^0.5:1./(1+x.^4+y.^4).^0.5g2 =1./(1+x.^4+y.^4).^0.5输入样本容量10^V1*10^V1--10^V2*10^V2:V1:1V2:4n =10n =100n =1000n =10000EY2 =3.0759 2.9699 2.8566 2.8269 MSE2 =1.3267 0.0900 0.0060 0.0014m =3.0759 2.9699 2.8566 2.82691.3267 0.0900 0.0060 0.0014实验结果整理：第一类一重积分：估计积分2sinx xdxπ⎰积分真值：1 积分估计值：1.0018样本容量：10 100 1000 10000 100000样本均值：1.2635 1.0088 1.0066 1.0109 1.0018 绝对误差：0.2635 0.0088 0.0066 0.0109 0.0018 相对误差：0.2635 0.0088 0.0066 0.0109 0.0018 均方误差：0.6439 0.0205 0.0028 0.0006 0.0001估计积分2-0x e dx +∞⎰积分真值：0.8862 积分估计值：0.8871 样本容量：101001000 10000样本均值：0.9333 0.9077 0.8873 0.8871 绝对误差：0.0470 0.0215 0.0010 0.0009 相对误差：0.0531 0.0243 0.0012 0.0010 均方误差：0.1927 0.0112 0.0016 0.0000第一类二重积分： 估计积分22221xy x y e dxdy ++≤⎰⎰积分真值：5.3981 积分估计值： 5.4041 样本容量：101001000 10000样本均值：4.7702 5.1250 5.4317 5.4041 绝对误差： 0.6279 0.2732 0.0335 0.0060 相对误差：0.1163 0.0506 0.0062 0.0011均方误差：3.8965 0.5564 0.0247 0.0017第二类一重积分： 估计积分 21x e dx ⎰积分估计值：1.4590 样本容量：101001000 10000样本均值：2.0782 1.6583 1.5029 1.4590 样本方差：0.4315 0.0889 0.0057 0.0008 用matlab 指令求得积分结果1.4627 第二类二重积分：估计积分22x y +≤⎰⎰积分估计值：2.8269 样本容量：101001000 10000样本均值：3.0759 2.9699 2.8566 2.8269 样本方差：1.3267 0.0900 0.0060 0.0014实验结果分析：从第一类积分看，以估计积分2sinx xdxπ⎰为例：积分真值：1 积分估计值：1.0018样本容量：10 100 1000 10000 100000样本均值：1.2635 1.0088 1.0066 1.0109 1.0018 绝对误差：0.2635 0.0088 0.0066 0.0109 0.0018 相对误差：0.2635 0.0088 0.0066 0.0109 0.0018 均方误差：0.6439 0.0205 0.0028 0.0006 0.0001 随着样本容量的增大，样本均值有接近积分真值的趋势，绝对误差、相对误差、均方误差呈减小趋势；随着样本容量的增大，样本均值有接近积分真值的趋势，说明估计具有无偏性；绝对误差、相对误差、均方误差呈减小趋势，说明增大样本容量能提高估计精度；验证了蒙特卡洛方法估计积分值的可行性，为后续估计第二类积分提供了参考。