# 概率分布采样与 NumPy 实践

[[MCS/day11/00 学习指南|day11 学习指南]] · [[MCS/00 MCS学习总览|十一天总览]] · [[MCS/day11/阅读核验与疑问|阅读核验与疑问]] · 上一篇：[[MCS/day11/08 伯努利样本、标准误与 Z 分数|伯努利样本、标准误与 Z 分数]] · 下一篇：[[MCS/day11/10 方差与取球概率例题|方差与取球概率例题]]

> 来源：[[MCS/day11 补充3|day11 补充3]]。用途：实践。原有推导、追问与解答保留；阅读时请结合本日核验说明。

Q 1 大数定律和中心极限定理区别?大数定律，样本均值趋近于期望

中心极限定理，样本均值的分布趋近于正态分布

Q 2 为什么赌场总是赢?

赌场本身有微弱优势，根据大数定律，长期重复必然盈利，赌客短期可赢，长期必输

Q 3 大数定律多少次有效?没有绝对数字，n 越大，误差越小

Q 4 概率为零的事件一定不发生吗?

不一定，连续概率空间中单点概率为零，但仍可能发生。

Q 5 如何用代码生成特定分布的数据?



numpyrandom 模块, np.random.normal, np.random.binomial
采用逆变换采样
 


## 用 Python 生成符合特定概率分布的数据 — 总结

### 一、核心思路

分两大类：**内置分布直接采样** 和 **自定义分布采样**。

### 二、常见方法一览

|需求|推荐方法|
|---|---|
|常见标准分布|numpy / scipy.stats 内置|
|自定义离散分布|rng.choice(p=...) 或 rv_discrete|
|自定义连续分布、CDF 可求|逆变换采样|
|自定义连续分布、CDF 难求|拒绝采样|
|多维相关分布|multivariate_normal / Copula|
|只有数据、无解析分布|KDE 或 bootstrap|
|复杂高维分布|MCMC（pymc、emcee）|

### 三、各方法要点

**1. NumPy / SciPy 内置分布（最常用）**

```python
rng = np.random.default_rng(42)
rng.uniform(0, 1, size=10000)
rng.normal(loc=0, scale=1, size=10000)
rng.binomial(n=10, p=0.5, size=10000)
rng.poisson(lam=3, size=10000)
rng.exponential(scale=2.0, size=10000)
from scipy import stats
stats.norm.rvs(loc=0, scale=1, size=10000, random_state=42)
stats.gamma.rvs(a=2.0, scale=1.0, size=10000)
stats.beta.rvs(a=2, b=5, size=10000)
```

建议优先用 scipy.stats，分布库最全，还能算 pdf、cdf、拟合参数。

**2. 离散自定义分布**

```python
values = ['A', 'B', 'C', 'D']
probs  = [0.5, 0.3, 0.15, 0.05]
samples = rng.choice(values, size=10000, p=probs)
# 或用 rv_discrete
custom = stats.rv_discrete(name='custom', values=([1,2,3,4], [0.1,0.2,0.3,0.4]))
samples = custom.rvs(size=10000, random_state=42)
```

**3. 连续自定义分布：逆变换采样**

原理：若 U ~ Uniform(0,1)，则 X = F⁻¹(U) 服从目标分布 F。

```python
x_grid = np.linspace(0, 4, 1000)
pdf_vals = custom_pdf(x_grid)
cdf_vals = np.cumsum(pdf_vals)
cdf_vals /= cdf_vals[-1]
inv_cdf = interp1d(cdf_vals, x_grid, bounds_error=False, fill_value=(0, 4))
u = rng.uniform(0, 1, size=10000)
samples = inv_cdf(u)
```

**4. 拒绝采样（CDF 难求时）**

```python
def rejection_sampling(n, x_min, x_max, M):
    samples = []
    while len(samples) < n:
        x = rng.uniform(x_min, x_max)
        y = rng.uniform(0, M)
        if y < target_pdf(x):
            samples.append(x)
    return np.array(samples)
```

**5. 多维 / 相关分布**

```python
mean = [0, 0]
cov  = [[1, 0.8], [0.8, 1]]
samples = rng.multivariate_normal(mean, cov, size=10000)
```

**6. 从真实数据学分布**

```python
# KDE
kde = stats.gaussian_kde(data)
samples = kde.resample(10000, seed=42).flatten()
# bootstrap
samples = np.random.choice(data, size=10000, replace=True)
```

**7. 验证生成结果**

```python
plt.hist(samples, bins=50, density=True, alpha=0.6)
stat, p = stats.kstest(samples, 'norm', args=(0, 1))
print(f'KS stat={stat:.4f}, p={p:.4f}')  # p > 0.05 不能拒绝"符合该分布"
```

### 四、总结口诀

1. 能用内置就用内置：NumPy、SciPy 覆盖 90% 常见需求。
    
2. 离散自定义：rng.choice(p=...) 最方便。
    
3. 连续自定义：优先逆变换采样；CDF 难求就用拒绝采样。
    
4. 高维 / 复杂：上 MCMC 或 Copula。
    
5. 生成后务必验证：画直方图 + KS 检验。
    

> 一句话：先看有没有现成分布，没有就"逆变换"或"拒绝采样"，再复杂就上 MCMC。

---

## np.random.binomial — 总结

### 一、作用

生成**二项分布**随机数。

场景：做 n 次独立试验，每次只有成功/失败两种结果，成功概率为 p，问 n 次里成功了多少次。

记作：**B(n, p)**

经典例子：

- 抛 10 次硬币，正面出现几次？（n=10, p=0.5）
    
- 射击 20 发，命中几发？（n=20, p=命中率）
    
- 100 个产品里，有几个次品？（n=100, p=次品率）
    

### 二、函数签名

```text
np.random.binomial(n, p, size=None)
```

|参数|含义|
|---|---|
|n|每次试验的试验次数|
|p|每次试验成功的概率（0 到 1）|
|size|生成多少个结果（输出形状）|

返回值：每次"n 次试验中成功的次数"，取值是 0 到 n 之间的整数。

### 三、代码示例

```python
import numpy as np
rng = np.random.default_rng(42)
result = rng.binomial(n=10, p=0.5)
print(result)   # 例如 6
results = rng.binomial(n=10, p=0.5, size=10000)
print(results[:10])        # [6 5 5 7 4 6 5 6 4 5 ...]
print(results.mean())      # 约等于 n*p = 5
print(results.std())       # 约等于 sqrt(n*p*(1-p)) ≈ 1.58
```

### 四、关键理解：n 和 size 的区别

```python
rng.binomial(n=10, p=0.5, size=5)
```

- n=10：每做一轮实验，抛 10 次硬币
    
- size=5：这样的"一整轮实验"做 5 轮
    
- 输出：5 个数，每个数代表那一轮 10 次里成功了几次
    

注意：n=1 时，二项分布退化为伯努利分布。

```python
rng.binomial(n=1, p=0.3, size=10)  # 输出 10 个 0 或 1
```

### 五、和别的分布的关系

|分布|关系|
|---|---|
|伯努利分布|二项分布在 n=1 时的特例|
|多项分布|二项分布的推广（结果从 2 种变成 k 种）|
|泊松分布|n 很大、p 很小时，二项分布近似泊松分布|
|正态分布|n 很大时，二项分布近似正态分布（中心极限定理）|

### 六、直观图示

```python
results = rng.binomial(n=10, p=0.5, size=10000)
plt.hist(results, bins=range(0, 12), align='left', density=True, edgecolor='black')
plt.xlabel('成功次数')
plt.ylabel('概率')
plt.title('二项分布 B(10, 0.5)')
plt.show()
```

会看到中间高、两边低的钟形（近似正态），峰值在 n*p = 5 附近。

### 七、理论公式

```text
P(X = k) = C(n, k) * p^k * (1-p)^(n-k)
```

- 期望：E[X] = np
    
- 方差：Var(X) = np(1-p)
    

### 八、总结

1. np.random.binomial(n, p, size) 生成二项分布随机数。
    
2. n = 每次实验的试验次数，p = 每次成功概率，size = 生成多少个结果。
    
3. 返回值是成功次数（0 到 n 的整数）。
    
4. n=1 时退化为伯努利分布；n 很大时近似正态分布。
    
5. 典型用途：模拟抛硬币、命中次数、次品数等"多次二选一试验的成功计数"。
    

> 一句话：它是"做 n 次成败试验、数成功几次"的模拟器。