Matlab生成随机数

2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
【例4.1-5】调用random函数生成10000×1的二项分布随机数向 量，然后作出频率直方图。其中二项分布的参数为 n=10，p=0.3
% 调用random函数生成10000行1列的随机数向量x，其元素服从二项分布 B(10,0.3) >> x = random('bino', 10, 0.3, 10000, 1); >> [fp, xp] = ecdf(x); % 计算经验累积概率分布函数值 >> ecdfhist(fp, xp, 50); % 绘制频率直方图 >> xlabel('二项分布（n = 10, p = 0.3）随机数'); % 为X轴加标签 >> ylabel('f(x)'); % 为Y轴加标签 1.5
19937 53
,1 253 ] 上
1) / 2
利用 Marsaglia's subtract with borrow 算法，在 MATLAB 5 至 7.3 中默认采用这个算法。 此算法能产生闭区间 [2 值
53
1492 理论上能产生 2 个不重复 ,1 253 ] 上的所有双精度值，
'seed'
利用乘同余算法（multiplicative congruential algorithm），是 MATLAB 4 中的默认算法。 此算法产生闭区间 [1/(2
31
1), 1 1/(231 1)] 上的双精度值，周期为 231 2
2013-7-26
©
谢中华, 天津科技大学数学系.
2013-7-26
©
谢中华, 天津科技大学数学系.
1. 创建RandStream类对象 RandStream类方法列表
方 法 说 明 RandStream 创建一个随机数流 创建多个独立的随机数流 获取 RandStream 类对象的属性 列出可用的随机数生成算法 获取默认随机数流 设置默认随机数流 重置一个流到它的初始内部状态 生成均匀分布伪随机数 生成标准正态分布伪随机数 生成离散均匀分布伪随机整数 生成一个随机排列
>> x = random('chi2', 10, 10000, 1); >> [fp, xp] = ecdf(x); % 计算经验累积概率分布函数值 >> ecdfhist(fp, xp, 50); % 绘制频率直方图 >> hold on >> t = linspace(0, max(x), 100); >> y = chi2pdf(t, 10); >> plot(t, y, 'r', 'linewidth', 3) >> xlabel('x ( \chi^2(10) )'); % 为X轴加标签 >> ylabel('f(x)'); % 为Y轴加标签 >> legend('频率直方图', '密度函数曲线') % 为图形加标注框
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
三、常见一元分布随机数
MATLAB统计工具箱中函数名以rnd三个字符结尾的 函数用来生成常见分布的随机数。例如： betarnd exprnd gamrnd lognrnd normrnd Beta分布 指数分布 Gamma分布 对数正态分布 正态分布
300 250 第一列 第二列 第三列
200
频数
150
100
50
0 40
50
60 70 80 90 正态分布随机数（ = 75, = 8）
100
110
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
【例4.1-4】调用normrnd函数生成1000×3的正态分布随机数矩 阵，其中第各列均值分别为0，15，40，标准差分别为1，2，3， 并作出各列的频数直方图
1
f(x)
0.5 0 0
1
2
3 4 5 6 7 8 二项分布（n = 10, p = 0.3）随机数
9
10
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
【例4.1-6】调用random函数生成10000×1的卡方分布随机数向 量，然后作出频率直方图，并与自由度为10的卡方分布的密度 函数曲线作比较。其中卡方分布的参数（自由度）为 10
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
四、任意一元分布随机数
1. 离散分布随机数
X p
x1 p1
x2 xn . p2 pn
定理 4-1 设一元随机变量 X 的分布函数为 F ( x) ， 令 Y F ( x) ，则 Y 服从 [0, 1] 上的均匀分布，记为 Y ~ U (0, 1) .
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
在MATLAB7.7以前的版本中，rand函数还可以这样调用： rand(method, s) s = rand(method) 其中method是字符串变量，它的可能取值如下表所列：
method 参数的取值 'twister' 说 明 利用 Mersenne Twister 算法（译为马特赛特旋转演算法，个人认为译为梅森旋转算法更 合适） 在 MATLAB7.4 及以后版本中默认采用这个算法。 ， 它是由 Makoto Matsumoto (松 本)和 Takuji Nishimura (西村)于 1997 年开发的。此算法产生闭区间 [2 的双精度值，周期为 (2 'state'
生成随机数
【例4.1-2】调用RandStream函数创建一个指定随机数生成算法 的RandStream类对象，然后利用对象的randn方法生成10×10 的标准正态分布随机数矩阵，并将矩阵按列拉长，画出频数直 方图
% 创建一个RandStream类对象s，其随机数生成器的算法为‘mlfg6331_64’， 初始种子为10，对象s的randn方法的算法为'Inversion' >> s = RandStream('mlfg6331_64', 'seed', 10, 'RandnAlg', 'Inversion'); % 调用对象s的randn方法生成10行10列的随机数矩阵x，其元素服从标准正 态分布 >> x = s.randn(10) >> y = x(:); % 将x按列拉长成一个列向量 >> hist(y) % 绘制频数直方图 >> xlabel('标准正态分布随机数'); % 为X轴加标签 >> ylabel(‘频数’); % 为Y轴加标签
©
谢中华, 天津科技大学数学系.
生成随机数
第一节 生成一元分布随机数
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
一、均匀分布随机数和标准正态分布随机数
1. rand函数 调用格式：
Y = rand Y = rand(n)
Y = rand(m,n)
Y = rand([m n]) Y = rand(m,n,p, …) Y = rand([m n p …]) Y = rand(size(A))
% 调用normrnd函数生成1000行3列的随机数矩阵x，其各列元素分别服从不同 的正态分布 >> x = normrnd(repmat([0 15 40], 1000, 1), repmat([1 2 3], 1000, 1), 1000, 3); >> hist(x, 50) % 绘制矩阵x每列的频数直方图 >> xlabel('正态分布随机数'); % 为X轴加标签 >> ylabel('频数'); % 为Y轴加标签 % 为图形加标注框 >> legend('\mu = 0, \sigma = 1','\mu = 15, \sigma = 2','\mu = 40, \sigma = 3')
生成随机数
create get list getDefaultStream setDefaultStream reset rand randn randi randperm
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
2. RandStream函数的调用方法 调用格式：
s = RandStream('gentype')
>> xvalue = [-2 -1 0 1 2]; % 定义向量xvalue >> xp = [0.05 0.2 0.5 0.2 0.05]; % 定义向量xp % 调用randsample函数生成100个服从指定离散分布的随机数 >> x = randsample(xvalue, 100, true, xp); >> reshape(x,[10 10]) % 调用randsrc函数生成10*10的服从指定离散分布的随机数矩阵 >> y = randsrc(10,10,[xvalue;xp])
生成随机数
2. randn函数 调用格式：
与rand函数类似
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
【例4.1-1】设置随机数生成器的算法为Mersenne Twister算法， 生成均匀分布随机数矩阵 % 设置随机数生成器的算法为Mersenne Twister算法，初始种 子为1 >> rand('twister',1); % 生成2行6列的随机数矩阵，其元素服从[0,1]上均匀分布 >> x1 = rand(2,6)
参数名 Seed RandnAlg 参数值 非负整数，默认值为 0 'Ziggurat'、'Polar' 或 'Inversion'，默认 值为'Ziggurat' 说 明
设定随机数种子，用来初始化随机数生成器 设定标准正态分布随机数函数 randn 所用算法
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
生成随机数
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
教 材
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
主要内容 生成一元分布随机数 生成多元分布随机数 蒙特卡洛方法
2013-7-26
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
二、RandStream类
MATLAB 7.7及以后的版本中，依然支持rand函数的 上述两种调用方式，但已经是过时的调用方式了，因 为MATLAB 7.7中对生成随机数作了重大调整，给出 了RandStream（随机数流）类，通过调用类的构造 函数并传递合适的参数可以创建类对象，然后调用类 对象的rand，randn，randi，randperm方法生成随 机数。
[…] = RandStream('gentype','param1',val1,'param2',val2,…) gentype有6个可能的取值：mcg16807、mlfg6331_64、 mrg32k3a、mt19937ar、shr3cong和swb2712，对应随机数生 成器的6个不同算法，默认值为mt19937ar
MATLAB统计工具箱中的randsample和randsrc函数
用来生成指定离散分布随机数
2013-7-26
©
谢中华, 天津科技大学数学系.
生成随机数
【例 4.1-7】 设离散总体 X 的分布列为
X p 2 1 0 1 2 0.05 0.2 0.5 0.2 0.05
分别调用 randsample 和 randsrc 函数生成 100 个服从该分 布的随机数。
poissrnd
randomΒιβλιοθήκη Baidu
2013-7-26
泊松分布
指定分布
©
谢中华, 天津科技大学数学系.
randsample 从有限总体中随机抽样
生成随机数
【例4.1-3】调用normrnd函数生成1000×3的正态分布随机数矩 阵，其中均值 为75，标准差为8 ，并作出各列的频数直方图
% 调用normrnd函数生成1000行3列的随机数矩阵x，其元素服从均值为75， 标准差为8的正态分布 >> x = normrnd(75, 8, 1000, 3); >> hist(x) % 绘制矩阵x每列的频数直方图 >> xlabel('正态分布随机数（\mu = 75, \sigma = 8）'); % 为X轴加标签 >> ylabel('频数'); % 为Y轴加标签 >> legend('第一列', '第二列', '第三列') % 为图形加标注框