Datawhale 第十九期 Numpy下 之 Task02：随机抽样

一、随机抽样模块

Numpy.random模块对Python内置的random进行了补充，增加了一些用于高效生成多种概率分布的样本值函数，比如正态分布、泊松分布等。

numpy.random.seed(seed=None)

seed() 用于指定随机数生成时所用算法开始的整数值，如果使用相同的seed() 值，则每次生成的随机数都相同，如果不设置这个值，则系统根据时间来自己选择这个值，此时每次生成的随机数因时间差异而不同。

当对数据进行预处理时，往往需要加入新的操作或改变处理策略，假如伴随着随机操作，需要指定唯一的随机种子，防止由于随机的差异影响后续的结果。

二、离散型随机变量

2.1 二项分布

二项分布可以用于单次实验只有两种结果，各结果对应的概率相等的多次实验的概率问题。比如抛硬币10次正面6次的概率等类似的问题。

二项分布概率函数的代码表示：

binom.pmf(k) = choose(n, k) pk (1-p)(n-k)

二项分布概率函数的数学表示：

numpy.random.binomial(n, p, size=None)

表示对一个二项分布进行采样，size表示采样的次数，n表示做了n重伯努利试验，p表示成功的概率，函数的返回值表示n中成功的次数

【例1】野外正在进行9（n=9）口石油勘探井的发掘工作，每一口井能够开发出油的概率是0.1（p=0.1）。请问，最终所有的勘探井都勘探失败的概率？

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats

np.random.seed(20201124)
n = 9
p = 0.1
size = 500000
x = np.random.binomial(n, p, size)
print('np.sum(x==0)/size=', np.sum(x == 0) / size)
# np.sum(x==0)/size= 0.38806

plt.hist(x)
plt.show()

s = stats.binom.pmf(range(10), n, p)
print('np.around(s,3)=',np.around(s, 3))
#np.around(s,3)= [0.387 0.387 0.172 0.045 0.007 0.001 0.    0.    0.    0.   ]

【例2】模拟投硬币，投2次，请问两次都为正面的概率?

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

np.random.seed(20201124)
n = 2
p = 0.5
size = 50000
x = np.random.binomial(n, p, size)
print("np.sum(x == 0) / size=", np.sum(x == 0) / size)  
# np.sum(x == 0) / size= 0.249814
print("np.sum(x == 1) / size=", np.sum(x == 1) / size)  
# np.sum(x == 1) / size= 0.500848
print("np.sum(x == 2) / size)=", np.sum(x == 2) / size)  
# np.sum(x == 2) / size)= 0.249338

plt.hist(x, density=True)
plt.show()

s = stats.binom.pmf(range(n + 1), n, p)
print("np.around(s, 3)=", np.around(s, 3))
# np.around(s, 3)= [0.25 0.5  0.25]

2.2 泊松分布

numpy.random.poisson(lam=1.0, size=None) Draw samples from a Poisson distribution.

表示对一个泊松分布进行采样，size表示采样的次数，lam表示一个单位内发生事件的平均值，函数的返回值表示一个单位内事件发生的次数。

【例3】假定某航空公司预定票处平均每小时接到42次订票电话，那么10分钟内恰好接到6次电话的概率是多少？

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

np.random.seed(20201124)
lam = 42 / 6
size = 500000
x = np.random.poisson(lam, size)
print("np.sum(x == 6) / size=", np.sum(x == 6) / size)  
# np.sum(x == 6) / size= 0.148998

plt.hist(x)
plt.show()

x = stats.poisson.pmf(6, lam)
print("x=", x)  # x= 0.1490027796743377

2.3 超几何分布

numpy.random.hypergeometric(ngood, nbad, nsample, size=None) Draw samples from a Hypergeometric distribution.

表示对一个超几何分布进行采样，size表示采样的次数，ngood表示总体中具有成功标志的元素个数，nbad表示总体中不具有成功标志的元素个数，ngood+nbad表示总体样本容量，nsample表示抽取元素的次数（小于或等于总体样本容量），函数的返回值表示抽取nsample个元素中具有成功标识的元素个数。

【例4】一共20只动物里有7只是狗，抽取12只有3只狗的概率（无放回抽样）。

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

np.random.seed(20201124)
size = 500000
x = np.random.hypergeometric(ngood=7, nbad=13, nsample=12, size=size)
print("np.sum(x == 3) / size=", np.sum(x == 3) / size)  
# np.sum(x == 3) / size= 0.19872

plt.hist(x, bins=8)
plt.show()

"""
M 为总体容量
n 为总体中具有成功标志的元素的个数
N,k 表示抽取N个元素有k个是成功元素
"""
x = range(8)
s = stats.hypergeom.pmf(k=x, M=20, n=7, N=12)
print("np.round(s, 3))=", np.round(s, 3))
# np.round(s, 3))= [0.    0.004 0.048 0.199 0.358 0.286 0.095 0.01 ]

三、连续型随机变量

3.1 均匀分布

numpy.random.uniform(low=0.0, high=1.0, size=None) Draw samples from a uniform distribution.

样品均匀分布在半开间隔[低，高）（包括低，但不包括高）。换句话说，给定间隔内的任何值都同样可能由均匀绘制。

【例5】在low到high范围内，创建大小为size的均匀分布的随机数。

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats

np.random.seed(20201124)
a = 0
b = 100
size = 500000
x = np.random.uniform(a, b, size=size)
print(np.all(x >= 0))  # True
print(np.all(x < 100))  # True
y = (np.sum(x < 50) - np.sum(x < 10)) / size
print("y=", y)  # y= 0.400444

plt.hist(x, bins=20)
plt.show()

a = stats.uniform.cdf(10, 0, 100)
b = stats.uniform.cdf(50, 0, 100)
print("b-a=",b - a)  # b-a=0.4

作为uniform()的特列，可以得到[0,1)之间的均匀分布的随机数。

numpy.random.rand(d0, d1, ..., dn) Random values in a given shape.

创建给定形状的数组，并用来自 +0、1 上统一分布的随机样本填充它。

【例6】根据指定大小产生[0,1)之间均匀分布的随机数。

import numpy as np

np.random.seed(20201124)
print("np.random.rand()=",np.random.rand())
# np.random.rand()= 0.16595762053434393

print("np.random.rand(5)=", np.random.rand(5))
# np.random.rand(5)= [0.18713139 0.16141027 0.08701148 0.89495402 0.44014387]

print("np.random.rand(4,3)=", np.random.rand(4, 3))
#np.random.rand(4,3)= [[0.45409661 0.15475296 0.70117664]
#[0.19671135 0.62914882 0.60799848]
#[0.69312633 0.07419549 0.46756484]
#[0.94894489 0.51506329 0.03977621]]

np.random.seed(20201124)
print("np.random.uniform()=", np.random.uniform())  
#np.random.uniform()= 0.16595762053434393
print("np.random.uniform(size=5)=", np.random.uniform(size=5))
#np.random.uniform(size=5)= [0.18713139 0.16141027 0.08701148 0.89495402 0.44014387]
print("np.random.uniform(size=(4,3))=", np.random.uniform(size=(4, 3)))
#np.random.uniform(size=(4,3))= [[0.45409661 0.15475296 0.70117664]
#[0.19671135 0.62914882 0.60799848]
#[0.69312633 0.07419549 0.46756484]
#[0.94894489 0.51506329 0.03977621]]

作为uniform的另一特例，可以得到[low,high)之间均匀分布的随机整数。

numpy.random.randint(low, high=None, size=None, dtype='l') Return random integers from low (inclusive) to high (exclusive).

【例7】若high不为None时，取[low,high)之间随机整数，否则取值[0,low)之间随机整数。

import numpy as np

np.random.seed(20201124)
x = np.random.randint(2, size=10)
print(x)
# [0 0 0 1 0 1 0 0 0 0]

x = np.random.randint(1, size=10)
print(x)
# [0 0 0 0 0 0 0 0 0 0]

x = np.random.randint(5, size=(2, 4))
print(x)
# [[3 3 0 1]
#  [1 1 0 1]]

x = np.random.randint(1, 10, [3, 4])
print(x)
# [[2 1 7 7]
#  [7 2 4 6]
#  [8 7 2 8]]

3.2 正态分布

numpy.random.randn(d0, d1, ..., dn) Return a sample (or samples) from the “standard normal” distribution.

【例】根据指定大小产生满足标准正态分布的数组（均值为0，标准差为1）。

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats

np.random.seed(20201124)
size = 50000
x = np.random.randn(size)
y1 = (np.sum(x < 1) - np.sum(x < -1)) / size
y2 = (np.sum(x < 2) - np.sum(x < -2)) / size
y3 = (np.sum(x < 3) - np.sum(x < -3)) / size
print("y1=", y1)  # y1= 0.68092
print("y2=", y2)  # y2= 0.9557
print("y3=", y3)  # y3= 0.99776

plt.hist(x, bins=20)
plt.show()

y1 = stats.norm.cdf(1) - stats.norm.cdf(-1)
y2 = stats.norm.cdf(2) - stats.norm.cdf(-2)
y3 = stats.norm.cdf(3) - stats.norm.cdf(-3)
print("y1=", y1)  # y1= 0.6826894921370859
print("y2=", y2)  # y2= 0.9544997361036416
print("y3=", y3)  # y3= 0.9973002039367398

另外，可以指定分布以及所需参数来进行随机，例如高斯分布中的mu和sigma。

numpy.random.normal(loc=0.0, scale=1.0, size=None) Draw random samples from a normal (Gaussian) distribution.

normal()为创建均值为 loc（mu），标准差为 scale（sigma），大小为 size 的数组。

sigma * np.random.randn(...) + mu

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(20201124)
x = 0.5 * np.random.randn(2, 4) + 5
print("x=",x)
# x= [[4.79670884 4.78295086 4.92900556 5.46844893]
#[4.40392523 4.92074704 4.53036121 5.31151958]]

#np.random.seed(20200614)
mu = 5
sigma = 0.5
x = np.random.normal(mu, sigma, (2, 4))
print("x=",x)
# x= [[5.66929259 5.80036633 4.76410781 5.1069904 ]
#[5.32335602 4.97663833 4.71316951 5.00938807]]

size = 50000
x = np.random.normal(mu, sigma, size)

print("np.mean(x)=", np.mean(x))  
# np.mean(x)= 4.9992922958417765
print("np.std(x, ddof=1)=", np.std(x, ddof=1))
# np.std(x, ddof=1)= 0.49928427081544485

plt.hist(x, bins=20)
plt.show()

3.3 指数分布

numpy.random.exponential(scale=1.0, size=None) Draw samples from an exponential distribution.

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats

np.random.seed(20200614)
lam = 7
size = 50000
x = np.random.exponential(1 / lam, size)
y1 = (np.sum(x < 1 / 7)) / size
y2 = (np.sum(x < 2 / 7)) / size
y3 = (np.sum(x < 3 / 7)) / size
print("y1=", y1)  # y1= 0.63218
print("y2=", y2)  # y2= 0.86518
print("y3=", y3)  # y3= 0.95056

plt.hist(x, bins=20)
plt.show()

y1 = stats.expon.cdf(1 / 7, scale=1 / lam)
y2 = stats.expon.cdf(2 / 7, scale=1 / lam)
y3 = stats.expon.cdf(3 / 7, scale=1 / lam)
print("y1=", y1)  # y1= 0.6321205588285577
print("y2=", y2)  # y2= 0.8646647167633873
print("y3=", y3)  # y3= 0.950212931632136

3.4 其它随机函数

numpy.random.choice(a, size=None, replace=True, p=None) Generates a random sample from a given 1-D array.

从序列中获取元素，若a为整数，元素取值从np.range(a)中随机获取；若a为数组，取值从a数组元素中随机获取。该函数还可以控制生成数组中的元素是否重复replace，以及选取元素的概率p。

import numpy as np

np.random.seed(20201124)
x = np.random.choice(10, 3)
print("x=", x)  # x= [7 6 2]

x = np.random.choice(10, 3, p=[0.05, 0, 0.05, 0.9, 0, 0, 0, 0, 0, 0])
print("x=",x)  # x=[3 3 3]

x = np.random.choice(10, 3, replace=False, p=[0.05, 0, 0.05, 0.9, 0, 0, 0, 0, 0, 0])
print("x=",x)  # x=[3 0 2]

aa_milne_arr = ['pooh', 'rabbit', 'piglet', 'Christopher']
x = np.random.choice(aa_milne_arr, 5, p=[0.5, 0.1, 0.1, 0.3])
print("x=", x)  # x= ['piglet' 'pooh' 'Christopher' 'pooh' 'piglet']

np.random.seed(20200614)
x = np.random.randint(0, 10, 3)
print("x=",x)  # x=[2 0 1]

参考文献
https://www.jianshu.com/p/63434ad5ea64