Python中的隨機取樣和概率分佈(二)

orion發表於2021-12-06
Python隨機概率分佈

在上一篇博文《Python中的隨機取樣和概率分佈(一)》(連結：https://www.cnblogs.com/orion-orion/p/15647408.html)中，我們介紹了Python中最簡單的隨機取樣函式。接下來我們更進一步，來看看如何從一個概率分佈中取樣，我們以幾個機器學習中最常用的概率分佈為例。

1. 二項（binomial）/伯努利（Bernoulli）分佈

1.1 概率質量函式(pmf)

\[P(X = x;\space n, \space p)=\left(\begin{array}{c}n \\ x\end{array}\right) p^{x}(1-p)^{n-x}\\ x=0,1,2,...n; \space 0\leqslant p \leqslant 1 \]

當\(n=1\)時，則取到下列極限情況，是為引數為\(p\)的二項分佈：

\[P(X = x;\space p)= p^{x}(1-p)^{1-x}\\ x=0,1; \space 0\leqslant p \leqslant 1 \]

二項分佈\(P(X = x;\space n, \space p)\)可以表示進行獨立重複試驗\(n\)次，每次有兩成功和失敗可能結果（分別對應概率\(p\)和\(1-p\)），共成功\(x\)次的概率。

1.2 函式原型

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

引數：
n: int or array_like of ints   對應分佈函式中的引數 n，>=0，浮點數會被截斷為整形。
p: float or array_like of floats   對應分佈函式引數\(p\), >=0並且<=1。
size: int or tuple of ints, optional   如果給定形狀為\((m, n, k)\)，那麼\(m\times n \times k\)個隨機樣本會從中抽取。預設為None，即返回一個一個標量隨機樣本。

返回：
out: ndarray or scalar  從帶引數的概率分佈中採的隨機樣本，每個樣本表示獨立重複實驗\(n\)次中成功的次數。

1.3 使用樣例

設進行獨立重複實驗10次，每次成功概率為0.5，取樣樣本表示總共的成功次數（相當於扔10次硬幣，正面朝上的次數）。總共採20個樣本。

import numpy as np
n, p = 10, .5  
s = np.random.binomial(n, p, 20)
print(s) # [4 5 6 5 4 2 4 6 7 2 4 4 2 4 4 7 6 3 5 6]

可以粗略的看到，樣本幾乎都在5周圍上下波動。
我們來看一個有趣的例子。一家公司鑽了9口井，每口井成功的概率為0.1，所有井都失敗了，發生這種情況的概率是多少？
我們總共取樣2000次，來看下產生0結果的概率。

s = sum(np.random.binomial(9, 0.1, 20000) == 0)/20000.
print(s) # 0.3823

可見，所有井失敗的概率為0.3823，這個概率還是蠻大的。

2. 多項（multinomial）分佈

2.1 概率質量函式(pmf)

\[P(\bm{X} = \bm{x};\space n, \space \bm{p})=\frac{n !}{x_{1} ! \cdots x_{k} !} p_{1}^{x_{1}} \cdots p_{k}^{x_{k}}\\ \bm{x}=(x_1, x_2, ..., x_k), x_i \in \{0, ..., n\}, \space \sum_{i}{x_i}=n; \\ \bm{p}=(p_1, p_2, ..., p_k), 0\leqslant p_i \leqslant 1, \space \sum_{i}{p_i}=1 \]

當\(k=2\)時，則取到下列極限情況，是為引數為\(n\), \(p\)的二項分佈：

\[P(X = x;\space n, \space p)=\frac{n !}{x !(n-x) !} p^{x}(1-p)^{n-x}\\ x=0,1,2,...n; \space 0\leqslant p \leqslant 1 \]

也就是說，多項分散式二項分佈的推廣：仍然是獨立重複實驗\(n\)次，但每次不只有成功和失敗兩種結果，而是\(k\)種可能的結果，每種結果的概率為\(p_i\)。多項分佈是一個隨機向量的分佈，\(\bm{x}=(x_1, x_2, ..., x_k)\)意為第\(i\)種結果出現\(x_i\)次，\(P(\bm{X} = \bm{x};\space n, \space p)\)也就表示第\(i\)種結果出現\(x_i\)次的概率。

2.2 函式原型

random.multinomial(n, pvals, size=None)

引數：
n: int   對應分佈函式中的引數 n。
pvals: sequence of floats   對應分佈函式引數\(\bm{p}\), 其長度等於可能的結果數\(k\)，並且有\(0 \leqslant p_i \leqslant 1\)。
size: int or tuple of ints, optional   為輸出形狀大小，因為採出的每個樣本是一個隨機向量，預設最後一維會自動加上\(k\)，如果給定形狀為\((m, n)\)，那麼\(m\times n\)個維度為\(k\)的隨機向量會從中抽取。預設為None，即返回一個一個\(k\)維的隨機向量。

返回：
out: ndarray   從帶引數的概率分佈中採的隨機向量，長度為可能的結果數\(k\)，如果沒有給定 size，則shape為 (k,)。

2.3 使用樣例

設進行獨立重複實驗20次，每次情況的概率為1/6，取樣出的隨機向量表示每種情況出現次數（相當於扔20次六面骰子，點數為0, 1, 2, ..., 5出現的次數）。總共採1個樣本。

s = np.random.multinomial(20, [1/6.]*6, size=1)
print(s) # [[4 2 2 3 5 4]]

當然，如果不指定size，它直接就會返回一個一維向量了

s = np.random.multinomial(20, [1/6.]*6)
print(s) # [4 1 4 3 5 3]

如果像進行多次取樣，改變 size即可：

s = np.random.multinomial(20, [1/6.]*6, size=(2, 2))
print(s)
# [[[4 3 4 2 6 1]
#   [5 2 1 6 3 3]]

#  [[5 4 1 1 6 3]
#   [2 5 2 5 4 2]]]

這個函式在論文<sup>[1]</sup>的實現程式碼<sup>[2]</sup>中用來設定每一個 client分得的樣本數：

for cluster_id in range(n_clusters): 
    weights = np.random.dirichlet(alpha=alpha * np.ones(n_clients))
    clients_counts[cluster_id] = np.random.multinomial(clusters_sizes[cluster_id], weights)
    # 一共扔clusters_sizes[cluster_id]次篩子，該函式返回骰子落在某個client上各多少次，也就對應著該client應該分得的樣本數

3.均勻（uniform）分佈

3.1 概率密度函式(pdf)

\[p(x; \space a, \space b)=\frac{1}{b-a} \]

均勻分佈可用於隨機地從連續區間\([a, b)\)內進行取樣。

3.2 函式原型

random.uniform(low=0.0, high=1.0, size=None)

引數：
low: float or array_like of floats, optional   對應分佈函式中的下界引數 a，預設為0。
high: float or array_like of floats   對應分佈函式中的下界引數 b，預設為1.0。
size: int or tuple of ints, optional   為輸出形狀大小，如果給定形狀為\((m, n, k)\)，那麼\(m\times n\times k\)的樣本會從中抽取。預設為None，即返回一個單一標量。

返回：
out: ndarray or scalar   從帶引數的均勻分佈中採的隨機樣本

3.3 使用樣例

s = np.random.uniform(-1,0,10)
print(s)
# [-0.9479594  -0.86158902 -0.63754099 -0.0883407  -0.92845644 -0.11148294
#  -0.19826197 -0.77396765 -0.26809953 -0.74734785]

4. 狄利克雷(Dirichlet)分佈

4.1 概率密度函式(pdf)

\[P(\bm{x}; \bm{\alpha}) \propto \prod_{i=1}^{k} x_{i}^{\alpha_{i}-1} \\ \bm{x}=(x_1,x_2,...,x_k),\quad x_i > 0 , \quad \sum_{i=1}^k x_i = 1\\ \bm{\alpha} = (\alpha_1,\alpha_2,..., \alpha_k). \quad \alpha_i > 0 \]

4.2 函式原型

random.dirichlet(alpha, size=None)

引數：
alpha: sequence of floats, length k   對應分佈函式中的引數向量 \(\alpha\)，長度為\(k\)。
size: int or tuple of ints, optional   為輸出形狀大小，因為採出的每個樣本是一個隨機向量，預設最後一維會自動加上\(k\)，如果給定形狀為\((m, n)\)，那麼\(m\times n\)個維度為\(k\)的隨機向量會從中抽取。預設為None，即返回一個一個\(k\)維的隨機向量。

返回：
out: ndarray   採出的樣本，大小為\((size, k)\)。

4.3 使用樣例

設\(\bm{\alpha}=(10, 5, 3)\)(意味著\(k=3\))，\(size=(2, 2)\)，則採出的樣本為\(2\times 2\)個維度為\(k=3\)的隨機向量。

s = np.random.dirichlet((10, 5, 3), size=(2, 2))
print(s)
# [[[0.82327647 0.09820451 0.07851902]
#   [0.50861077 0.4503409  0.04104833]]

#  [[0.31843167 0.22436547 0.45720285]
#   [0.40981943 0.40349597 0.1866846 ]]]

這個函式在論文[1]的實現程式碼[2]中用來生成符合狄利克雷分佈的權重向量

for cluster_id in range(n_clusters): 
    # 為每個client生成一個權重向量，文章中分佈引數alpha每一維都相同
    weights = np.random.dirichlet(alpha=alpha * np.ones(n_clients))
    clients_counts[cluster_id] = np.random.multinomial(clusters_sizes[cluster_id], weights)

參考文獻

相關文章