随机数产生的问题
这2天发现这个问题讨论较多,所以就搜索了一些资料。发现自己之前的理解有些很有问题,同时欢迎大家继续讨论。先澄清一下几个容易弄错的地方(也不一定全对)(1)用计算机产生的是“伪随机数”。用投色子计数的方法产生真正的随机数 , 但电脑若也这样做 , 将会占用大量内存 ; 用噪声发生器或放射性物质也可产生真正的随机数 , 但不可重复 . 而用数学方法产生最适合计算机 , 这就是周期有限 , 易重复的 ” 伪随机数 ”
(2)随机数的产生需要有一个随机的种子,因为用计算机产生的随机数是通过递推的方法得来的,必须有一个初始值。
(3)用同一台电脑,且在初始值和递推方法相同的情况下,可以产生相同的随机序列(由于以前每次使用randn或者rand得到都是不同值,所以曾经误以为相同的seed无法产生相同的序列)
一 matlab里产生随机数的方法
matlab里和随机数有关的函数:
(1) rand:产生均值为0.5、幅度在0~1之间的伪随机数
(2) randn:产生均值为0、方差为1的高斯白噪声
(3) randperm(n):产生1到n的均匀分布随机序列
(4) normrnd(a,b,c,d):产生均值为a、方差为b大小为cXd的随机矩阵
还有很多的扩展函数,不再一一列出。不过他们都调用的是rand或者randn函数,由此可见在matlab里rand和randn是产生随机数的关键所在。看来只有看他们的源文件了
function = randn(varargin)
%%%help 文档的内容略去%%%
if nargout == 0
builtin('randn', varargin{:});
else
= builtin('randn', varargin{:});
end
从这里也看不出到底是怎么产生的,就只看到builtin。而builtin函数的源文件是这样的:
%BUILTINExecute built-in function from overloaded method.
% BUILTIN is used in methods that overload built-in functions to execute
% the original built-in function. If F is a string containing the name
% of a built-in function then BUILTIN(F,x1,...,xn) evaluates that
% function at the given arguments.
%
% BUILTIN(...) is the same as FEVAL(...) except that it will call the
% original built-in version of the function even if an overloaded one
% exists (for this to work, you must never overload BUILTIN).
%
% = BUILTIN(F,x1,...,xn) returns multiple output arguments.
%
% See also FEVAL.
% Copyright 1984-2002 The MathWorks, Inc.
% $Revision: 1.9 $$Date: 2002/04/15 04:16:04 $
% Built-in function.
后来发现matlab基本函数的源文件都是这么一个结构
function = functionname(varargin)
%%%help 文档的内容%%%
if nargout == 0
builtin('functionname', varargin{:});
else
= builtin('functionname', varargin{:});
end
其中一些是.m文件一些是.bi文件,数学里最基本的运算就是加减乘除(对于计算机来说就只有加法器和乘法器)。那么最简单的sin函数如何用四则运算求值?反正我搞不清楚,sin.m内容也是一个builtin函数。扯了这么多,只是为了说明在使用matlab基本函数的时候,很多情况下看源文件并不能知道其中具体用的是什么数值方法。
有了rand和randn就可以产生轻松产生均匀分布和正态分布的随机数了
(1)产生在区间服从均匀分布随机序列的方法
(b-a)*rand(m,n)+a
>> 3*rand(2)+2
ans =
2.8166 2.0458
2.5964 4.2404
(2)产生服从正态分布的随机数
>> randn('state',2)
>> a=normrnd(0,1,1,6)
a =
1.7491 0.1326 0.3252 -0.7938 0.3149 -0.5273
>> randn('state',2)
>> b=randn(1,6)
b =
1.7491 0.1326 0.3252 -0.7938 0.3149 -0.5273
>> randn('state',2)
>> c=randn(2,3)
c =
1.7491 0.3252 0.3149
0.1326 -0.7938 -0.5273
>> d=randn(2,3)
d=
0.9323 -2.0457 1.7411
1.1647 -0.6444 0.4868
>> mean(a)
ans =
0.2001
--------------------------
>> randn(1,2)
ans =
1.0488 1.4886
>> randn(1,2)
ans =
1.2705 -1.8561
---------------------------
上边几个典型的例子可以看出:
(1)如果不设置种子,那么种子会“随机”变化。每次使用randn就会得到不同的结果(c和d)
(2)种子相同时可以得到相同的结果,如果是矩阵那么只是将产生的随机数按列重构(a、b、c)
(3)randn无法准确保证均值为0,小样本的时候尤为明显。去均值后可以严格保证均值为0,但是个人觉得意义不大。
(4)在不同的计算里得到的结果也可能有差别,特别是不同的操作系统。大家可以试一下这个语句
randn('state',2);randn(1,6)看看结果,我电脑每次都一样的
二 编程产生随机数的方法
在使用rand和randn的时候,用的具体方法我一直没有找到。那么如果要自己编程该如何实现呢?由于我对这块也不熟悉,通过搜到的资料来看,相关的方法和研究文献应该相当的多。所以如果是做这方面的朋友可以多找找看,里边深层的东西还很多。这里给出2个链接,权做抛砖引玉。里边有不同随机序列的不同生成方法,我就不啰唆了,感兴趣的朋友自己去看。搞研究嘛,遇到问题,要多动手,多搜索。
http://www.blog.edu.cn/user2/33128/archives/2006/1242828.shtml
http://blog.csdn.net/EmilMatthew/archive/2006/04/21/672276.aspx
最后欢迎各位版友发表自己的看法和问题,不讨论也就不容易发现问题和错误
[ 本帖最后由 eight 于 2007-9-5 15:02 编辑 ] 关于这个语句: randn('state',2);randn(1,6)
在我的电脑里面,也是一样的结果。
只要是输入“randn('state',j); randn(m,n)"
出来的结果一定是相同的。
j,m,n变成什么都一样。
但如果随后再生成一个randn(m,n),而不重新设定seed,就是与第一个序列不同的了,不过同一个seed每次生成的第二个序列也都是一样的。
这大约说明randn的伪随机过程还是按照某种固定模式迭代的吧。
[ 本帖最后由 不化顽石 于 2007-9-5 16:21 编辑 ] 另外,关于花兄说的”去均值后可以严格保证均值为0,但是个人觉得意义不大。“我想补充一点:
我现在正在做的一个工作是与Monte Carlo Simulation相关的,它要求生成标准白噪声样本,是因为在随后的过程中,要用到“样本标准差为1”来进一步计算。
所以,如果在这种情况下严格地对已有序列标准化还是必须的,而不能完全相信randn函数,虽然它在方差上造成的偏差要远小于均值,然而不知道前面这一点小偏差在后面会造成什么结果。 原帖由 不化顽石 于 2007-9-5 16:20 发表 http://www.chinavib.com/forum/images/common/back.gif
关于这个语句: randn('state',2);randn(1,6)
在我的电脑里面,也是一样的结果。
只要是输入“randn('state',j); randn(m,n)"
出来的结果一定是相同的。
j,m,n变成什么都一样。
但如果随后再生成一个randn ...
请你再试试,种子不同时。得到的结果肯定是不同的
>>randn('state',2);randn(1,6)
ans =
1.7491 0.1326 0.3252 -0.7938 0.3149 -0.5273
>> randn('state',3);randn(1,6)
ans =
0.9280 0.1733 -0.6916 -0.7230 -0.5744 -0.3077
随机数的产生就是根据种子(用来决定初始值),然后按照一定的规则外推。具体的方法有很多,可以看我附的链接。而randn用的是什么方法,我也弄不清楚,因为没法看到其真正的源文件
回复 #3 不化顽石 的帖子
谢谢你的意见,这方面的资料应该是比较多的。你可以看看我在“高斯信号和非高斯信号”帖子里提到的那篇文章[ 本帖最后由 花如月 于 2007-9-5 16:41 编辑 ] 我可能没写清楚。
种子不同的时候,确实不一样,我的意思是,只要确定了种子,生成的就是一样的了。 另外,关于花兄你的那篇贴子,我有一点疑问,就是下面两句:
subplot(121),plot(y); %此处可以看到高斯白噪声的频谱还是白噪声
subplot(122),plot(real(Y)); %高斯函数的傅立叶变换还是高斯函数
这里是不是有点问题?按说频谱不是这么看的。
回复 #7 不化顽石 的帖子
那个是有些不太准确,我已改过。用abs好些,不过结果基本是一样的。那个结论有严格的数学证明,可以推导出来的另外,我对频谱其实了解也不多少,只是看到别人都是这么用的。你觉得应该怎么看频谱呢?
[ 本帖最后由 花如月 于 2007-9-5 17:11 编辑 ] 现在用得比较多的,其实是功率谱。
你所说的abs(Y),应该是能量谱。
我是做大气的,接触小波分析要比纯粹的FFT多一些,在小波中,单看实部图或者虚部图,从原理上讲,都是不能反映出序列的振荡性质的。
当然这部分并不是我的方向,所以了解的也许不是很透彻,期待再交流。
页:
[1]