Python中幂律分布的拟合

Question

我使用不同的python来适应数据集中的密度函数。此数据集由从1秒开始的正时间值组成。

我测试了来自scipy.statistics和powerlaw库的不同密度函数，以及使用scipy.optimize的函数curve_fit()测试自己的函数。

到目前为止，当拟合下列修正的幂律函数时，我得到了最好的结果：

def funct(x, alpha, x0):
    return((x+x0)**(-alpha))

我的代码如下：

bins = range(1,int(s_distrib.max())+2,1)
y_data, x_data = np.histogram(s_distrib, bins=bins, density=True)
x_data = x_data[:-1]

param_bounds=([0,-np.inf],[np.inf,np.inf])
fit = opt.curve_fit(funct,
                    x_data,
                    y_data,
                    bounds=param_bounds) # you can pass guess for the parameters/errors
alpha,x0 = fit[0]
print(fit[0])

C = 1/integrate.quad(lambda t: funct(t,alpha,x0),1,np.inf)[0]

# Calculate fitted PDF and error with fit in distribution
pdf = [C*funct(x,alpha,x0) for x in x_data]
sse = np.sum(np.power(y_data - pdf, 2.0))
print(sse)

fig, ax = plt.subplots(figsize=(6,4))
ax.loglog(x_data, y_data, basex=10, basey=10,linestyle='None',  marker='.')
ax.loglog(x_data, pdf, basex=10, basey=10,linestyle='None',  marker='.')

拟合返回的值为x0为8.48，alpha为1.40。在日志图中，数据和fit图如下所示：

​

​

• 我的第一个问题是技术。为什么在opt.curve_fit函数中将(x+x0)更改为(x-x0)时，在funct中得到以下警告和错误？由于我对x0的界限是(-inf，+inf)，所以我期望拟合值返回-8.48。

/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:3: RuntimeWarning:除以倒数中遇到的零，这与ipykernel包是分开的，所以我们可以避免进行导入，直到ValueError:残留物在初始点上没有限制。

• ，我的其他问题是理论上的。(x+x0)^(-alpha)是一个标准发行版吗？x0值代表什么，如何物理地解释这个8.48s值？根据我的理解，这意味着我的分布对应于移动幂律分布？当数据与幂律相匹配时，x0是否与经典需要的xmin值相对应？
• 关于这个xmin值，我知道只考虑大于这个阈值的数据是有意义的，因为拟合过程可以描述分布的尾部。但是，我想知道用一个分布来描述整个数据的标准方法是什么，这个分布在xmin之后是幂律，在xmin之前是其他的。

这是很多问题，因为我非常不熟悉的主题，任何评论和回答，即使是部分，将非常感谢！

Answer 1

(x+x0)^(-alpha)是一个标准发行版吗？

要回答第二个问题，是的，它是标准发行版，称为Zipf分布。它是用Python/NumPy 也是实现的。

x0值代表什么？

这是移位参数。在标准参数之上的任何分布(如Zipf中的功率参数)都可能有移位和缩放参数，这基本上意味着您的X值是在不同起源点的不同单元中测量的。

关于这个xmin值，我知道只考虑大于这个阈值的数据是有意义的，因为拟合过程可以描述分布的尾部。

这就是Zipf定律的定义，从0到无限。改变它意味着你的起源会不一样