概率论与数理统计实验报告2500字

概率论与数理统计实验报告2500字

概率论与数理统计实验报告应物12班郭帅2110903026一、实验内容：用蒙特卡洛方法估计积分值?2??1用蒙特卡洛方法估计积分xsinxdx，edx和00?1?x2x2?y2?1??ex2?y2并将估计值与真dxdy的值，值进行比较。2用蒙特卡洛方法估计积分edx和0?x2x2?y2???的值，并对误差进行估计。二、要求：（1）针对要估计的积分选择适当的概率分布设计蒙特卡洛方法；（2）利用计算机产生所选分布的随机数以估计积分值；（3）进行重复试验，通过计算样本均值以评价估计的无偏性；通过计算均方误差（针对第1类题）或样本方差（针对第2类题）以评价估计结果的精度。三、目的：（1）能通过MATLAB或其他数学软件了解随机变量的概率密度、分布函数及其期望、方差、协方差等；（2）熟练使用MATLAB对样本进行基本统计，从而获取数据的基本信息；（3）能用MATLAB熟练进行样本的一元回归分析。蒙特卡洛方法：当所求解问题是某种随机事件出现的概率，或者是某个随机变量的期望值时，通过某种“实验”的方法，以这种事件出现的频率估计这一随机事件的概率，或者得到这个随机变量的某些数字特征，并将其作为问题的解。四、实验步骤：?2（1）xsinxdx0?方法：x在0至pi/2区间上随机取10000个数为均匀分布的简单随机样本,然后计算y的值一共计算二十次,即可用样本均值作为积分的估计值.Y=pi/2*x.*sin(x)y*f(x)即为被积函数??2?,x?[0,]f(x)???2?其他?0,clcclearx=rand(20,10000)*pi/2y=(pi/2)*x.*sin(x)a=sum(y,2)/10000u=sum(a,1)/20H=1E=abs(H-u)b=abs(H-u)^2D=sum(b,1)/19结果样本均值为u=0.9987E=0.0013D=8.971e-008=0.00000008971真值计算:clcclearsymsxf='x*sin(x)'int(f,x,0,pi/2)结果真值为1??(2)?0exdx2方法:x在负无穷到正无穷之间按标准正态分布取10000个样本，然后计算y值二十次，即可用样本均值估计积分值。clcclearx=randn(20,10000)y=exp((-x.^2)/2)*sqrt(2*pi)a=sum(y,2)/10000b=sum(a,1)/20u=b/2H=1/2*pi^(1/2)E=abs(H-u)c=abs(H-u)^2D=sum(c,1)/19结果样本均值为：0.8870误差=0.00076131样本方差=3.0505e-008=0.0000000.0505真值：symsxf='exp(-x^2)'int(f,x,0,inf)真值为：1/2*pi^(1/2)(3)x2?y2?12?1??ex2?y2dxdy原积分=?00??e?d?d?2方法:把?看做x,使x在0到1之间取10000个值坐在0至1区间上的均匀分布的随机简单样本,计算样本对应的y值,因为y与?无关,所以第二重积分可直接乘以2?得样本均值,一共计算二十次,即可用样本均值作为积分的估计值.y=exp(x^2)*xf(x)??clcclear1,x?[0,1]0,其他x=rand(20,10000)g=x.*exp(x.^2)a=sum(g,2)/10000b=sum(a,1)/20u=2*pi*bH=pi*(exp(1)-1)E=abs(H-u)c=abs(H-u)^2D=sum(c,1)/19结果样本均值为5.4140误差=0.0159样本方差=1.3294e-005=0.000013294真值计算:clcclearsymsxyf='exp(x^2)*x'int(int(f,x,0,1),y,0,2*pi)结果pi*(exp(1)-1)=5.3981二1(1)?e0x2dx方法：x在0至1区间上随机取10000个数为均匀分布的简单随机样本,然后计算y的值一共计算二十次,即可用样本均值作为积分的估计值.Y=exp(x^2)f(x)??clcclear1,x?[0,1]0,其他x=rand(20,10000)y=exp(x.^2)a=sum(y,2)/10000u=sum(a,1)/20b=(a-u).^2D=sum(b,1)/19U为样本均值,计算结果样本均值=1.4627D为样本方差,样本方差=1.0156e-005=0.000010156（2）x2?y2???方法:随机产生20组数据每组10000个,保证这些数据为均匀分布在单位园内的简单随机样本,然后计算z的值,即可用z在单位圆上均匀分布的期望值对积分制进行估计.Z=pi/(sqrt(1+x^4+y^4))?122?,x?y?1f(x)?????0,其他clcclearm=rand(20,10000)n=rand(20,10000)a=rand(20,10000)*2*pib=rand(20,10000)*2*pix=m.*cos(a)y=n.*sin(b)z=pi./(sqrt(1+x.^4+y.^4))b=sum(z,2)/10000u=sum(b,1)/20c=(b-u).^2D=sum(c,1)/19计算结果,U为样本均值,u=2.9621D为样本方差,D=2.6849e-006=0.0000026849

第二篇：概率论与数理统计实验报告2400字概率论与数理统计实验报告——使用Mathematica工具软件实验人：郭鹏院系专业班级：物理学院应用物理系0801学号：U200810167试验时间：2009.12.23——2010-01-18实验题目：实验主题：使用mathematic工具研究非线性回归问题下面的数据来自对某种遗传特征的研究结果,一共有2723对数据,把它们分成846620.9117719.379819.36类后归纳为下表.频率分类变量x遗传性指标y579138.081021229.7607325.42324423.15120521.79研究者通过散点图认为y和x符合指数关系:y?aebx?c,其中a,b,c是参数.求参数a,b,c的最小二乘估计.实验前分析：因为y和x的关系不是能用Fit命令拟合的线性关系,也不能转换为线性回归模型.因此考虑用(1)多元微积分的方法求a,b,c的最小二乘估计;(2)非线性拟合命令NonlinearFit求a,b,c的最小二乘估计.实验过程：方法(1)微积分方法：输入Off[Genera1::spe11]Off[Genera1::spe111]Clear[x,y,a,b,c]dataset={{579,1,38.08},{1021,2,29.70},{607,3,25.42},{324,4,23.15},{120,5,21.79},{46,6,20.91},{17,7,19.37},{9,8,19.36}};y[x_]:=aExp[bx]+c下面一组命令先定义了曲线y?aebx?c与2723个数据点的垂直方向的距离平方和,记为g(a,b,c).再求g(a,b,c)对a,b,c的偏导数?g?g?g,,,?a?b?c分别记为ga,gb,gc.用FindRoot命令解三个偏导数等于零组成的方程组(求解a,b,c).其结果就是所要求的a,b,c的最小二乘估计.输入Clear[a,b,c,f,fa,fb,fc]g[a_,b_,c_]:=Sum[dataset[[i,1]]*(dataset[[i,3]]-a*Exp[dataset[[i,2]]*b]-c)^2,{i,1,Length[dataset]}]ga[a_,b_,c_]=D[g[a,b,c],a];gb[a_,b_,c_]=D[g[a,b,c],b];gc[a_,b_,c_]=D[g[a,b,c],c];Clear[a,b,c]oursolution=FindRoot[{ga[a,b,c]==0,gb[a,b,c]==0,gc[a,b,c]==0},{a,40.},{b,-1.},{c,20.}](*40是a的初值,-1是b的初值,20是c的初值*)则输出{a->33.2221,b->-0.626855,c->20.2913}再输入yhat[x_]=y[x]/.oursolution则输出20.2913+33.2221e?0.626855x这就是y和x的最佳拟合关系.输入以下命令可以得到拟合函数和数据点的图形:p1=Plot[yhat[x],{x,0,12},PlotRange->{15,55},DisplayFunction->Identity];pts=Table[{dataset[[i,2]],dataset[[i,3]]},{i,1,Length[dataset]}];p2=ListPlot[pts,PlotStyle->PointSize[.01],DisplayFunction->Identity];Show[p1,p2,DisplayFunction->$DisplayFunction];则输出图2.4.如图一方法(2)直接用非线性拟合命令NonlinearFit方法输入data2=Flatten[Table[Table[{dataset[[j,2]],dataset[[j,3]]},{i,dataset[[j,1]]}],{j,1,Length[dataset]}],1];(*把数据集恢复成2723个数对的形式*)<<Statistics`w=NonlinearFit[data2,a*Exp[b*x]+c,{x},{{a,40},{b,-1},{c,20}}]则输出?0.626855x20.2913?33.2221e这个结果与(1)的结果完全相同.这里同样要注意的是参数a,b,c必须选择合适的初值.?i)2.输入如果要评价回归效果,则只要求出2723个数据的残差平方和?(yi?yyest=Table[yhat[dataset[[i,2]]],{i,1,Length[dataset]}];yact=Table[dataset[[i,3]],{i,1,Length[dataset]}];w