|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
随机数序列在数值分析和概率统计中占有非常重要的地位,因为使用蒙特卡罗模拟方法的前提就是要求很多足够多的,真正的随机数。matlab是基于某种算法,通过rand函数来产生随机数的。从随机数的定义看,rand函数产生的序列不是随机数,是伪随机数。但我们在使用蒙特卡罗模拟方法算法时,不可能成千上万次的去投掷硬币来产生随机数,所以要考察matlab产生的伪随机数能不能当随机数使用。
<P dir=ltr style="MARGIN-RIGHT: 0px"> 考察的方法是:</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 1:利用rand函数产生200个伪随机数,分别统计出有多少个伪随机数小数点后的5位数字中奇数出现的频数为0,有多少个伪随机数小数点后的5位数字中奇数出现的频数为1,有多少个伪随机数小数点后的5位数字中奇数出现的频数为2,一直统计到有多少个伪随机数小数点后的5位数字中奇数出现的频数为5。 这样就获得一个包含6个元素的行向量s,其元素依次为小数点后5位包含奇数个数为0或1或2……或5的伪随机数的个数。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 2:给出两组源于实际观测的数据。一组是记录某个医院相继出生的1000个婴儿的随机序列,记5个连续出生的婴儿为一组,这样共有200组,统计分别含有0,1,2……5个男婴的组数,获得向量a;另一组是从一个装有500黑球和500个白球的口袋里,每次有放回的抽取一个球,供抽取1000次,记连续抽样5次为一组,这样共有200组,统计分别含有0,1,2……5个白球的组数,获得向量b。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 3:对向量a,b,s进行自由度为5卡方检验,分别获得以向量a,b,s为代表的三组数据的χ<SUP>2</SUP>值。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 4:将上述步骤重复1000次,每次向量a,b都是不变的,但每次的s向量都不同。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 5:计算1000组s对应的1000个χ<SUP>2</SUP>值的最大值,最小值,平均值,对1000组a,b同样如此。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 6:如果假设显著水平为0.05,那么自由度为5的χ<SUP>2</SUP>分布临界值是11.1,所以还要计算1000组s对应的1000个χ<SUP>2</SUP>值中大于11.1的数值所占的分率。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 总结果如下:</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> a b s(基于matlab)</P>
<P dir=ltr style="MARGIN-RIGHT: 0px">平均值:2.2240 5.0400 5.0038</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"><BR>最大值:2.2240 5.0400 19.3760</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"><BR>最小值: 2.2240 5.0400 0.2560</P>
<P dir=ltr style="MARGIN-RIGHT: 0px">1000个χ<SUP>2</SUP>值中大于11.1的数值所占的分率</P>
<P dir=ltr style="MARGIN-RIGHT: 0px">x = 0.0520</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 从平均值的计算结果看,matlab产生的伪随机数的随机程度和从口袋摸球相当,所以随机性满足要求。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 从最大值结果看,基于matlab的伪随机序列产生的χ<SUP>2</SUP>值最大达到19.3760,大于显著水平0.05,自由度为5的χ<SUP>2</SUP>分布的临界值,似乎有些序列不够随机。但考虑到χ<SUP>2</SUP>分布中,总有0.05的概率,使得χ<SUP>2</SUP>值大于11.1,所以验算了基于matlab的伪随机序列产生的χ<SUP>2</SUP>值中大于11.1的数值占的分率,这个分率是0.0520,非常接近0.05。</P>
<P dir=ltr style="MARGIN-RIGHT: 0px"> 所以,<B>matlab产生的伪随机序列可以作为真正的随机序列使用。</B></P>
<P dir=ltr style="MARGIN-RIGHT: 0px"><B> </B>matlab程序如下:</P>
<BLOCKQUOTE dir=ltr style="MARGIN-RIGHT: 0px">
<P dir=ltr style="MARGIN-RIGHT: 0px">clear<BR>rr=[];<BR>for l=1:1000<BR>p=[];<BR>rand('seed',prod(clock))<BR>r=fix(rand(200)*100000);<BR>for i=1:length(r)<BR> m=r(i); s=[];<BR> for j=1:5<BR> s=[s rem(m,10)];<BR> m=round((m-rem(m,10))/10);<BR> end<BR> p=[p;s];<BR>end<BR>s=zeros(1,6);<BR>for i=1:size(p,1)<BR> k=length(find(rem(p(i,:),2*ones(1,size(p,2)))));<BR> s(k+1)=s(k+1)+1;<BR>end<BR>a=[5 27 64 65 30 9];<BR>b=[4 34 65 70 22 5];<BR>c=[];<BR>for i=0:5<BR> p=combine_m(5,i)*0.5^i*0.5^(5-i); <BR> c=[c 200*p];<BR>end<BR>ka2=sum([([a;b;s]-[c;c;c]).^2./[c;c;c]]');<BR>rr=[rr;ka2];<BR>end<BR>ave=mean(rr)<BR>mx=max(rr)<BR>mi=min(rr)<BR>x=length(find(rr(:,3)>11.1))/length(rr(:,3))</P></BLOCKQUOTE> |
|