python - 对于相关分布采样，是否有 scipy _norm_pdf 的快速替代方案？

Question

我已经为蒙特卡罗模拟拟合了一系列 SciPy 连续分布，并且希望从这些分布中获取大量样本。但是，我希望能够获取相关样本，例如 i第一个样本从每个分布中获取例如第 90 个百分位数。

在这样做时，我发现了 SciPy 性能的一个怪癖:

# very fast way to many uncorrelated samples of length n
for shape, loc, scale, in distro_props:
    sp.stats.norm.rvs(*shape, loc=loc, scale=scale, size=n)

# verrrrryyyyy slow way to take correlated samples of length n
correlate = np.random.uniform(size=n)
for shape, loc, scale, in distro_props:
    sp.stats.norm.ppf(correlate, *shape, loc=loc, scale=scale)

大多数关于此的结果都声称这些 SciPy 发行版的缓慢，如果来自类型检查等包装。然而，当我分析代码时，大部分时间都花在了底层数学函数 [_continuous_distns.py:179(_norm_pdf)] 上。 1 .此外，它与 n 一起缩放，这意味着它在内部循环遍历每个元素。

SciPy docs on rv_continuous几乎似乎建议子类应该为了性能而覆盖它，但是我将monkeypatch到SciPy以加速他们的ppf似乎很奇怪。我只是根据 ppf 公式计算法线，但我也使用对数正态和偏斜正态，这更难实现。

那么，在 Python 中计算正态、对数正态和偏态正态分布的快速 ppf 的最佳方法是什么？或者更广泛地说，从几个这样的分布中获取相关样本？

Answer 1

如果您只需要普通的 ppf ，这么慢确实让人费解，不过可以用scipy.special.erfinv反而:

x = np.random.uniform(0,1,100)
np.allclose(special.erfinv(2*x-1)*np.sqrt(2),stats.norm().ppf(x))
# True
timeit(lambda:stats.norm().ppf(x),number=1000)
# 0.7717257660115138
timeit(lambda:special.erfinv(2*x-1)*np.sqrt(2),number=1000)
# 0.015020604943856597

编辑:
lognormal和 triangle也是直截了当的:

c = np.random.uniform()

np.allclose(np.exp(c*special.erfinv(2*x-1)*np.sqrt(2)),stats.lognorm(c).ppf(x))
# True

np.allclose(((1-np.sqrt(1-(x-c)/((x>c)-c)))*((x>c)-c))+c,stats.triang(c).ppf(x))
# True

skew normal 我不够熟悉，不幸的是。