生成特定分布随机数的方法

生成随机数是程序设计里常见的需求。一般的编程语言都会自带一个随机数生成函数,用于生成服从均匀分布的随机数。不过有时需要生成服从其它分布的随机数,例如高斯分布或指数分布等。有些编程语言已经有比较完善的实现,例如Python的NumPy。这篇文章介绍如何通过均匀分布随机数生成函数生成符合特定概率分布的随机数,主要介绍Inverse Ttransform和Acceptance-Rejection两种基础算法以及一些相关的衍生方法。下文我们均假设已经拥有一个可以生成0到1之间均匀分布的随机数生成函数,关于如何生成均匀分布等更底层的随机数生成理论,请参考其它资料,本文不做讨论。

基础算法

Inverse Transform Method

最简单的生成算法是Inverse Transform Method(下文简称ITM)。如果我们可以给出概率分布的累积分布函数(下文简称CDF)及其逆函数的解析表达式,则可以非常简单便捷的生成指定分布随机数。

ITM算法描述

  1. 生成一个服从均匀分布的随机数U~Uni(0,1)
  2. 设F(X)为指定分布的CDF,F?1(Y)是其逆函数。返回X=F?1(U)作为结果

ITM算法说明

这是一个非常简洁高效的算法,下面说明其原理及正确性。

我们通过图示可以更直观的明白算法的原理。下图是某概率分布的CDF:

我们从横轴上标注两点xa和xb,其CDF值分别为F(xa)和F(xb)。

由于U服从0到1之间的均匀分布,因此对于一次U的取样,U落入F(xa)和F(xb)之间的概率为:

P(U∈(F(xa),F(xb)))=F(xb)?F(xa)

而由于CDF都是单调非减函数,因此这个概率同时等于X落入xa和xb之间的概率,即:

P(U∈(F(xa),F(xb)))=P(F?1(U)∈(F?1(F(xa)),F?1(F(xb))))=P(X∈(xa,xb))=F(xb)?F(xa)

而根据CDF的定义,这刚好说明X服从以F(x)为CDF的分布,因此我们的生成算法是正确的。

ITM实现示例

下面以指数分布为例说明如何应用ITM。

首先我们需要求解CDF的逆函数。我们知道指数分布的CDF为

F(X)=1?e?λX

通过简单的代数运算,可以得到其逆函数为

F?1(Y)=?1λln(1?Y)

由于U服从从0到1的均匀分布蕴含着1?U服从同样的分布,因此在实际实现时可以用Y代替1?Y,得到:

F?1(Y)=?1λln(Y)

下面给出一个Python的实现示例程序。

  1. import random
  2. import math
  3. def exponential_rand(lam):
  4. if lam <=0:
  5. return-1
  6. U = random.uniform(0.0,1.0)
  7. return(-1.0/ lam)* math.log(U)

Acceptance-Rejection Method

一般来说ITM是一种很好的算法,简单且高效,如果可以使用的话,是第一选择。但是ITM有自身的局限性,就是要求必须能给出CDF逆函数的解析表达式,有些时候要做到这点比较困难,这限制了ITM的适用范围。

当无法给出CDF逆函数的解析表达式时,Acceptance-Rejection Method(下文简称ARM)是另外的选择。ARM的适用范围比ITM要大,只要给出概率密度函数(下文简称PDF)的解析表达式即可,而大多数常用分布的PDF是可以查到的。

ARM算法描述

  1. 设PDF为f(x)。首先生成一个均匀分布随机数X~Uni(xmin,xmax)
  2. 独立的生成另一个均匀分布随机数Y~Uni(ymin,ymax)
  3. 如果Y≤f(X),则返回X,否则回到第1步

ARM算法说明

通过一幅图可以清楚的看到ARM的工作原理。

ARM本质上是一种模拟方法,而非直接数学方法。它每次生成新的随机数后,通过另一个随机数来保证其被接受概率服从指定的PDF。

显然ARM从效率上不如ITM,但是其适应性更广,在无法得到CDF的逆函数时,ARM是不错的选择。

ARM实现示例

下面使用ARM实现一个能产生标准正态分布的随机数生成函数。

首先我们要得到标准正态分布的PDF,其数学表示为:

f(x)=1√2πe?x22

为了方便,这里我会直接使用SciPy来计算其PDF。

程序如下。

  1. import random
  2. import scipy.stats as ss
  3. def standard_normal_rand():
  4. whileTrue:
  5. X = random.uniform(-3.0,3.0)
  6. Y = random.uniform(0.0,0.5)
  7. if Y < ss.norm.pdf(X):
  8. return X

注意:标准正态分布的x取值范围从理论上说是(?∞,∞),但是当离开均值点很远后,其概率密度可忽略不计。这里只取(?3.0,3.0),实际使用时可以根据具体需要扩大这个取值范围。

衍生算法

组合算法

当目标分布可以用其它分布经过四则运算表示时,可以使用组合算法生成对应随机数。

最常见的就是某分布可以表示成多个独立同分布(下文简称IID)随机变量之和。例如二项分布可以表示成多个0-1分布之和,Erlang分布可以由多个IID的指数分布得出。

以Erlang分布为例说明如何生成这类随机数。

设X1,X2,?,Xk为服从0到1均匀分布的IID随机数,则?1λlnX1,?1λlnX2,?,?1λlnXk为服从指数分布的IID随机数,因此

X=?1λlnX1?1λlnX2???1λlnXk=?1λlnk∏i=1Xi~Erl(k,λ)

所以生成Erlang分布随机数的算法如下:

  1. 生成X1,X2,?,Xk~Uni(0,1)
  2. 返回?1λln∏ki=1Xi

这类分布的随机数生成算法很直观,就是先生成相关的n个IID随机数,然后带入简单求和公式或其它四则公式得出最终随机数。其数学理论基础是卷积理论,稍微有些复杂,这里不再讨论,有兴趣的同学可以查阅相关资料。

生成具有相关性的随机数

现在考虑生成多维随机数,以最简单的二维随机数为例。

如果两个维度的随机数是相互独立的,那么只要分别生成两个列就可以了。但是如果要求两列具有一定的相关系数,则需要做一些特殊处理。

下列算法可以生成两列具有相关系数ρ的随机数。

  1. 生成IID随机变量X和Y
  2. 计算X′=ρX+√1?ρ2Y
  3. 返回(X,X′)

可以这样验证其正确性:

再分享一下我老师大神的人工智能教程吧。零基础!通俗易懂!风趣幽默!还带黄段子!希望你也加入到我们人工智能的队伍中来!https://blog.csdn.net/jiangjunshow

原文地址:https://www.cnblogs.com/skiwnchiwns/p/10345428.html

时间: 2024-10-17 09:39:44

生成特定分布随机数的方法的相关文章

C#生成不重复随机数的方法

在使用Random类生成随机数时,我们可能会碰到生成随机数重复的问题. 比如我们要生成6位数字验证码,虽然也是使用Random,但是可能出现111111,999999这样的情况. 这是因为在实例化Random类时,如果随机种子不填写,默认是以时间线作为种子进行伪随机运算,当计算运行速度过快时,导致所有的随机种子都是一个值. 解决的方法也很简单,我们使用Guid的哈希码作为种子值,就不会重复了,代码如下: 1 public class RandomHelper 2 { 3 /// <summary

一个php生成16位随机数的代码(两种方法)

一个php生成16位随机数的代码,php生成随机数的二种方法. 方法1<?php$a = mt_rand(10000000,99999999);$b = mt_rand(10000000,99999999);echo $a.$b; 方法2:<?php$a = range(0,9);for($i=0;$i<16;$i++){$b[] = array_rand($a);} // www.yuju100.comvar_dump(join("",$b));//结果string

数组应用实例(生成并打印随机数和统计随机数的分布)

一.生成并打印随机数 代码如下: <span style="font-size:18px;">#include <stdio.h> #include <stdlib.h> #define N 20 int a[N]; void gen_random(int upper_bound) //生成随机数在0-upper_bound之间 { int i; for(i = 0;i<N;i++) { a[i] = rand()%upper_bound; }

js生成随机数的方法实例总结 [收藏]

js生成随机数的方法实例总结 js生成随机数主要用到了内置的Math对象的random()方法.用法如:Math.random().它返回的是一个 0 ~ 1 之间的随机数.有了这么一个方法,那生成任意随机数就好理解了.比如实际中我们可能会有如下的需要: (1)生成一个 0 - 100 之间的随机整数,那么则可以: parseInt(100*Math.random()); 注意:因为Math.random()的返回值是包括0和1的,所以这里是有生成0和100的可能性的. (2)生成一个从 m -

C#生成互不相同随机数的实现方法

本文实例讲述了C#生成互不相同随机数的实现方法,在进行C#应用程序设计时非常具有实用价值.本文详细讲述了其功能的实现过程.分享给大家供大家参考之用.具体方法如下: 一般来说,用C#生成足够随机的互不相同的随机数 Dotnet.Frameword中提供了一个专门产生随机数的类System.Random,计算机并不能产生完全随机的数字,它生成的数字被称为伪随机数,它是以相同的概率从一组有限的数字中选取的,所选的数字并不具有完全的随机性,但就实用而言,其随机程度已经足够了. 在使用随机数时,要先初始化

已知一个函数rand5()能够生成1-5的随机数,请给出一个函数,该函数能够生成1-7的随机数。

这是朋友去笔试的一道题,有点考智商,当时我还很自信的说 random5+random5/2  不就可以了 他说不行,然后我就在网上搜了一下 有一道类似的题目 题目: 已知一个函数rand7()能够生成1-7的随机数,请给出一个函数,该函数能够生成1-10的随机数. 思路: 假如已知一个函数能够生成1-49的随机数,那么如何以此生成1-10的随机数呢? 解法: 该解法基于一种叫做拒绝采样的方法.主要思想是只要产生一个目标范围内的随机数,则直接返回.如果产生的随机数不在目标范围内,则丢弃该值,重新取

C# Random 生成不重复随机数

Random 类 命名空间:System 表示伪随机数生成器,一种能够产生满足某些随机性统计要求的数字序列的设备. 伪随机数是以相同的概率从一组有限的数字中选取的.所选数字并不具有完全的随机性,因为它们是用一种确定的数学算法选择的,但是从实用的角度而言,其随机程度已足够了. 伪随机数的生成是从种子值开始.如果反复使用同一个种子,就会生成相同的数字系列.产生不同序列的一种方法是使种子值与时间相关,从而对于 Random 的每个新实例,都会产生不同的系列.默认情况下,Random 类的无参数构造函数

Linq技术三:Linq to Object 和生成数据表的扩展方法

这篇来谈论一下Linq第三个方面的应用:Linq to Object,只要是继承了IEnumerable或IQueryable接口的Object都能使用Linq特性进行操作.在操作过程当中可能很多人都觉得不好调试不能实时地观察结果数据集,想把IQuery的Linq查询语句转换成数据表DataTable,要怎么实现转换呢?来看一下. 先来说一场景解释一下为什么需要用Linq来解决一些问题,能解决一些什么样的问题,相对于SQL,DataTable等一些传统操作方式有哪些优势? 场景:目前主要数据源有

java 产生随机数的方法

有三种方法: Math.random():这个方法返回一个[0.0, 1.0)的一个随机double型数.它实际是调用Random类的nextDouble()方法.只不过Math类使用的是一个静态随机数生成器(即new Random()),是线程安全的一个方法,所以多个线程共用一个随机数生成器.如果很多线程都在频繁的使用随机数生成器,那么还是为每个线程分配一个随机数生成器比较好.此外在J2ME中好像不支持这个方法. Random类:java提供的强大的伪随机数生成器.有两个构造方法和主要的六类获