从python中的对数正态分布生成随机数

Question

我需要从Python中的对数正态分布生成伪随机数.问题是我从对数正态分布的模式和标准偏差开始.我没有对数正态分布的均值或中值,也没有基础正态分布的任何参数.

numpy.random.lognormal取基础正态分布的均值和标准差.我试图从我的参数计算这些,但最后用四次函数.它有一个解决方案,但我希望有一个更简单的方法来做到这一点.

scipy.stats.lognorm采取我不理解的参数.我不是母语为英语的人,文档没有意义.

你能帮我吗？

Answer 1

您具有对数正态分布的模式和标准偏差.要使用rvs()SciPy的的方法lognorm,你必须在参数形状参数方面的分布s,这是标准偏差sigma的基础正态分布的,而且scale,这是exp(mu),这里mu是底层分布的均值.

您指出进行此重新参数化需要求解四次多项式.为此,我们可以使用numpy.poly1d该类.该类的实例具有roots属性.

一个小代数表明这exp(sigma**2)是多项式的唯一正实根

x**4 - x**3 - (stddev/mode)**2 = 0

其中stddev和mode是给定的对数正态分布的标准偏差和模式,对于该解决方案,scale(即exp(mu))是

scale = mode*x

这是一个将模式和标准偏差转换为形状和比例的函数:

def lognorm_params(mode, stddev):
    """
    Given the mode and std. dev. of the log-normal distribution, this function
    returns the shape and scale parameters for scipy's parameterization of the
    distribution.
    """
    p = np.poly1d([1, -1, 0, 0, -(stddev/mode)**2])
    r = p.roots
    sol = r[(r.imag == 0) & (r.real > 0)].real
    shape = np.sqrt(np.log(sol))
    scale = mode * sol
    return shape, scale

例如,

In [155]: mode = 123

In [156]: stddev = 99

In [157]: sigma, scale = lognorm_params(mode, stddev)

使用计算的参数生成样本:

In [158]: from scipy.stats import lognorm

In [159]: sample = lognorm.rvs(sigma, 0, scale, size=1000000)

这是样本的标准偏差:

In [160]: np.std(sample)
Out[160]: 99.12048952171304

这里有一些matplotlib代码用于绘制样本的直方图,并在绘制样本的分布模式下绘制垂直线:

In [176]: tmp = plt.hist(sample, normed=True, bins=1000, alpha=0.6, color='c', ec='c')

In [177]: plt.xlim(0, 600)
Out[177]: (0, 600)

In [178]: plt.axvline(mode)
Out[178]: <matplotlib.lines.Line2D at 0x12c5a12e8>

直方图:

如果要使用numpy.random.lognormal()而不是生成样本scipy.stats.lognorm.rvs(),可以执行以下操作:

In [200]: sigma, scale = lognorm_params(mode, stddev)

In [201]: mu = np.log(scale)

In [202]: sample = np.random.lognormal(mu, sigma, size=1000000)

In [203]: np.std(sample)
Out[203]: 99.078297384090902

我没有研究过鲁棒poly1d的roots算法是多少,所以一定要测试各种可能的输入值.或者,您可以使用scipy中的求解器来求解上述多项式x.您可以使用以下方法绑定解

max(sqrt(stddev/mode), 1) <= x <= sqrt(stddev/mode) + 1