尺度参数在scipy分布中产生了位移。

Question

我正在使用 SciPy的levy_stable 分布做一些事情。文档说：

levy_stable.pdf(x, alpha, beta, loc, scale) 与 y = (x-loc) / scale 及 levy_stable.pdf(y, alpha, beta) / scale 等价。

那么，我理解 levy_stable.pdf(x, alpha=1, beta=1, loc=0, scale=10) 应该是 levy_stable.pdf(x, alpha=1, beta=1, loc=0, scale=1) 的拉伸版本，但不会有位移。然而，我得到了这个：

请注意，在scale=1时，最大值为负值，但在scale>=5时，最大值为正值，这意味着它已经发生了偏移。为什么会这样？是一个错误还是我漏掉了什么？在其他更简单的分布中，比如norm，loc和scale参数按照其名称所指示的方式进行平移和拉伸分布。实际上，对于the moyal distribution，它应该是对另一个分布的近似，它的行为也符合预期：

生成这些图表的代码如下：

from scipy.stats import levy_stable, moyal
import numpy
import plotly.express as px
import pandas

data = []
for loc in [0]:
    for scale in [1,5,10]:
        distributions = {
            'levy_stable': levy_stable(alpha=1,beta=1,loc=loc,scale=scale), 
            'moyal': moyal(loc=loc,scale=scale),
        }
        for name,dist in distributions.items():
            x = numpy.linspace(loc-5*scale,loc+30*scale,999)
            df = pandas.DataFrame(
                {
                    'x': x,
                    'x/scale': x/scale,
                    'y': dist.pdf(x),
                }
            )
            df['loc'] = loc
            df['scale'] = scale
            df['distribution'] = name
            data.append(df)
data = pandas.concat(data)

fig = px.line(
    data,
    x = 'x/scale',
    y = 'y',
    facet_row = 'scale',
    color = 'distribution',
)
fig.update_yaxes(matches=None)
fig.write_html(f'deleteme/levy_stable_bug.html')

- user171780

Answer 1

这看起来像是一个错误。我在这里提交了一个错误报告：https://github.com/scipy/scipy/issues/19140

Levy稳定分布的PDF计算很复杂。没有闭式解，SciPy实现了几种不同的计算方法。SciPy中的levy_stable实现有一个额外的属性叫做parameterization，可以设置来控制内部参数化。默认值是levy_stable.parameterization = 'S1'，而这个版本有一个错误。作为一种解决方法，添加以下代码：

levy_stable.parameterization = 'S0'

将代码修改为选择替代参数化方式。目前看来，S0参数化方式似乎没有这个bug。

（到目前为止，我只在SciPy 1.11.2版本上进行了测试。）

- Warren Weckesser

---

使用'S0'参数化对我来说似乎没有改变，它完全一样。 - user171780

奇怪，这对我有效（但我没有运行你的脚本）。我的初步尝试没有成功，但那是因为我拼错了 parameterization。如果你正在使用一个IDE或交互式shell，请尝试完全退出程序，并用一个普通的 python 命令来运行你的代码。你正在使用哪个版本的SciPy？如果你不知道，你可以用 import scipy; print(scipy.__version__) 来检查。 - Warren Weckesser

Answer 2

levy_stable的值是正确的，但是文档遗漏了一个关键事实，即如果alpha == 1，则所声称的等价性不适用于S1参数化。

在S1参数化中，如果alpha == 1，则levy_stable.pdf(x, alpha, beta, loc, scale)等同于levy_stable.pdf(y, alpha, beta) / scale，其中y = (x - loc - 2 * beta * scale * np.log(scale) / np.pi) / scale。请参见Stable Distributions: Models for Heavy Tailed Data，第8-9页。

代码演示

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import levy_stable


alpha = 1
beta = 1
loc = 1
scale = 3


x = np.linspace(-5, 10, 1000)
y1 = levy_stable.pdf(x, alpha, beta, loc=loc, scale=scale)
y2 = levy_stable.pdf((x-loc)/scale, alpha, beta)/scale

plt.plot(x, y1, '--', linewidth=1,
         label='pdf(x, alpha, beta, loc=loc, scale=scale)')
plt.plot(x, y2, alpha=0.4, linewidth=3,
         label='pdf((x-loc)/scale, alpha, beta)/scale')

# Calculate scaled y.
y = (x - loc - 2 * beta * scale * np.log(scale) / np.pi) / scale
y3 = levy_stable.pdf(y, alpha, beta)/scale
plt.plot(x, y3, alpha=0.4, linewidth=3,
         label='pdf(y, alpha, beta)/scale')

plt.legend(shadow=True, framealpha=1)
plt.grid(alpha=0.25)
plt.savefig('levy_stable_plot_resolved.png')

- TMBailey