从均匀分布到任意分布
生成器给出在 0 和 1 之间均匀分布的数。模型需要别的形状：胜率 30% 的硬币、两笔成交之间的等待时间、正态的日收益率。诀窍始终相同：对均匀分布的数做变换。最通用的是逆变换法：若 F 是你想要的分布的累积分布函数，那么对均匀的 u，F⁻¹(u) 恰好服从该分布。整个数组可以一次完成变换。
import math
import numpy as np
# Exponential waiting times with rate 2 (a mean of 0.5): F(x) = 1 - exp(-rate * x), so x = -ln(1 - u) / rate
u = np.array([0.1, 0.5, 0.9])
print((-np.log(1 - u) / 2).round(4).tolist()) # [0.0527, 0.3466, 1.1513] small u, short wait
rng = np.random.default_rng(42)
waits = -np.log(1 - rng.random(100000)) / 2
print(round(waits.mean(), 2)) # 0.5 the mean is 1 / rate离散结果与正态分布
对于离散的结果表，rng.choice(values, size=n, p=probs) 替你完成逆变换：按各自概率抽取每个值，比如某天以 20%、50%、30% 的概率为 -2%、0%、+3%。对正态分布，其分布函数的反函数没有简单公式，经典技巧是 Box-Muller 变换：由两个独立的均匀数 u1 和 u2，z = sqrt(-2 ln u1) * cos(2 pi u2) 就是一个标准正态抽样。实际中 rng.standard_normal(n) 更快，但自己写一次就会明白这里没有魔法。
rng = np.random.default_rng(3)
days = rng.choice(np.array([-0.02, 0.0, 0.03]), size=100000, p=[0.2, 0.5, 0.3])
print(round(days.mean(), 3), round(0.2 * -0.02 + 0.3 * 0.03, 3)) # 0.005 0.005 the sample mean matches the table
rng = np.random.default_rng(42)
u1, u2 = 1 - rng.random(100000), rng.random(100000) # 1 - random() is in (0, 1]: never log(0)
z = np.sqrt(-2 * np.log(u1)) * np.cos(2 * math.pi * u2)
print(abs(z.mean()) < 0.02, round(z.std(), 2)) # True 1.0 mean 0 and standard deviation 1从标准正态抽样到价格
标准正态抽样 z 是一个冲击；模型决定它如何推动价格。股票的常用模型几何布朗运动说，t 年后的价格是 S = S0 * exp((mu - 0.5 sigma²) t + sigma sqrt(t) z)。它是对数正态的：恒为正，并具有真实价格那样的右侧长尾。-0.5 sigma² 这一项校正了均值，使期望价格恰好以 mu 增长：E[S] = S0 * exp(mu t)。中位数低于均值，所以典型结果比平均结果差。
s0, mu, sigma, t = 100, 0.08, 0.2, 1
z = np.array([-1.0, 0.0, 1.0]) # a one-sigma loss, no shock, a one-sigma gain
prices = s0 * np.exp((mu - 0.5 * sigma ** 2) * t + sigma * math.sqrt(t) * z)
print(prices.round(2).tolist()) # [86.94, 106.18, 129.69] all positive, the upside is larger
draws = np.random.default_rng(0).standard_normal(100000)
simulated = s0 * np.exp((mu - 0.5 * sigma ** 2) * t + sigma * math.sqrt(t) * draws)
print(round(simulated.mean(), 1), round(s0 * math.exp(mu * t), 1)) # 108.3 108.3 the mean grows at mu
print(round(float(np.median(simulated)), 1)) # 106.1 the median is lower# NumPy ships ready-made samplers too; the transformations above show what they do inside
z = np.random.default_rng(0).standard_normal(100000)
print(abs(z.mean()) < 0.02, round(z.std(), 2)) # True 1.0 standard_normal: mean 0, standard deviation 1用已知事实检验采样器写完任何采样器后，抽取一个大样本，把它的均值、标准差或某个概率与精确值比较：速率为 2 的指数分布均值 0.5，正态分布为 0 和 1，价格为
S0 * exp(mu t)。符号错误或漏掉平方根会产生看似合理却错误的数，只有这样的检验才能发现。- 逆变换法：对均匀数应用分布函数的反函数；指数分布：
-ln(1 - u) / rate。 rng.choice(values, size=n, p=probs)从离散表抽样；Box-Mullersqrt(-2 ln u1) * cos(2 pi u2)生成正态抽样。t年后的价格：S0 * exp((mu - 0.5 sigma²) t + sigma sqrt(t) z)，对数正态；务必用精确均值检验大样本。
下一节：蒙特卡洛估计，用随机抽样计算没有公式的平均值和概率。
实验：评分内容
- 用
sample_discrete从结果表中抽样25 分 - 用
exponential_from_uniform反转分布函数25 分 - 用
box_muller生成钟形分布的数25 分 - 用
terminal_prices把冲击转换为未来价格25 分