幂律分布拟合神器——Python库powerlaw

导语

具有长尾特征的分布往往一目了然，但实际拟合过程却可能遇到各种各样的问题。本文将为读者介绍2014年由新加坡科技设计大学和麻省理工研究者联合发布的python库：powerlaw，专门适用于幂律等长尾特征分布的拟合，解决拟合烦恼。

胡一冰 | 作者

邓一雪 | 编辑

期刊来源：PLOS ONE

论文标题：

powerlaw: A Python Package for Analysis of Heavy- Tailed Distributions

论文网址：

https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0085777

1. 幂律分布简要回顾

幂律分布作为长尾分布的一种，满足可以用来印证生活中很多有趣的现象，比如最为大家熟知的“二八定律”“富者越富”等等。早在上世纪早期就先后有克莱伯定律、Zipf定律等幂律现象的发现，而千禧年之际BA无标度网络的正式提出更是掀起一波“幂律热潮”，大家越来越接受“幂律分布是复杂系统中的普适现象”这一观点。

图1. 人们常说“20%的人掌握着80%的财富”就是典型的幂律分布，其具有尖峰、长尾特征。

2. powerlaw库拟合效果

图2. 三个不同数据集用powerlaw拟合后的可视化图像

图2显示了可视化、拟合和评估长尾分布的基本要素，每组子图将在后续部分中进行进一步的详细描述。

图1中三列子图包含三个示例数据集和拟合结果，分别表示幂律良好拟合、中等拟合和较差拟合。第一列的数据，也是最好的拟合数据集可能是所有幂律分布中最著名和最可靠的：英语单词使用频率。具体使用的数据是Herman Melville的小说《Moby Dick》中词语的使用频率。第二列的数据是线虫每个神经元的连接数。最后一个拟合不佳的数据是1984年至2002年间受停电影响的美国人数。

图2中的子图A表示了三个样本数据集的概率密度函数，子图B显示了只有少数分布的尾部可能遵循幂律，子图C表示了如何将幂律拟合与其他长尾分布进行比较。

2. powerlaw库基本操作介绍

可视化

下面以1984年至2002年间受停电影响的美国人数数据为例进行操作说明。

>import powerlaw
>fit = powerlaw.Fit(data)  #对导入的数据data进行幂律拟合
>fit.power_law.alpha  #拟合得到幂律分布幂指数alpha2.273  
>fit.distribution_compare(‘power_law’, ‘exponential’)  #与指数分布对比得到对数似然比和p值
(12.755, 0.152)

powerlaw库能容易绘制概率密度函数（PDF）、累积分布函数（CDF）和互补累积分布函数（CCDF）。对应的计算工作通过pdf、cdf和ccdf完成的，而绘图则是结合matplotlib通过plot_pdf、plot_cdf和plot_ccdf命令完成。

>powerlaw.plot_pdf(data, color = ‘b’)  
>powerlaw.plot_pdf(data, linear_bins = True, color = ‘r’)

该命令中的linear_bins参数控制对数坐标轴与否。线性刻度不利于幂律尾部“大值”的观察，采用对数刻度增加在尾部观察一系列数据的可能性。图2子图A显示了选择对数刻度如何极大地提高数据分布的可视化。在这里默认为对数刻度，但也可以用linear_bins = True进行线性划分。

关于PDF、CDF和CCDF的拟合示例如下，具体选择需结合拟合需求。如果要在同一画布中呈现多条拟合曲线，需要传递带有ax的matplotlib对象。例如下方同时画出神经元连接数据的CCDF（红色）和PDF（蓝色）曲线。

>fig2 = fit.plot_pdf(color = ‘b’, linewidth = 2)>fit.power_law.plot_pdf(color = ‘b’, linestyle = ‘–’, ax = fig2)
>fit.plot_ccdf(color = ‘r’, linewidth = 2, ax = fig2)
>fit.power_law.plot_ccdf(color = ‘r’, linestyle = ‘–’, ax = fig2

       图3. 神经元连接数据的CCDF（红色）和PDF（蓝色）曲线。

在绘图之外，PDF、CDF和CCDF信息也可用。拟合对象返回拟合数据和排序数据（CDF）和bin的边缘概率（PDF）。分布对象默认使用全部拟合数据，也可以指定具体的数据范围。

>x, y = fit.cdf()  #幂律分布CDF拟合
>bin_edges, probability = fit.pdf()
>y = fit.lognormal.cdf(data = [300, 350])  #对数正态分布CDF拟合
>y = fit.lognormal.pdf()

拟合范围

拟合幂律的第一步是确定拟合数据的哪一部分。长尾分布的有趣特征是尾部及其属性，所以如果数据的初始小值不遵循幂律分布，用户有可能会选择忽略它们。powerlaw库会为用户提供最优的幂律拟合最小值，当然用户也可以选择指定最小值或给最小值框定范围，不同的最小值选择往往会带来不一样的拟合结果，如下方代码所示。

>fit = powerlaw.Fit(data)
>fit.xmin  #得到最优的幂律拟合最小值230.000
>fit.power_law.alpha
2.273
>fit.power_law.D  #选择数据和拟合之间的Kolmogorov-Smirnov距离D0.061
>fit = powerlaw.Fit(data, xmin = 1.0)  #指定x_min进行幂律拟合
>fit.xmin
1.0
>fit.fixed_xmin
True
>fit.power_law.alpha
1.220
>fit.power_law.D
0.376
>fit = powerlaw.Fit(data, xmin = (250.0, 300.0))  #给定x_min范围
>fit.given_xmin
(250.000, 300.000)
>fit.xmin
272.0

另外，一些领域也希望分布有一个精确的上界x_max，因为有时考虑实际情况一些数据可能超出了理论极限（例如在天体物理学中，速度的分布在光速下可能有一个上限）。powerlaw库同样能满足这部分用户的需求。

>fit = powerlaw.Fit(data, xmax = 10000.0)
>fit.xmax
10000.0
>fit.fixed_xmax
True

       图4. 不同的x_min和x_max选择可能导致不同的拟合结果。

这里的x_min和x_max默认使用KS距离D最小的值。然而，这个距离D对分布尾部的差异明显不敏感，而尾部正是幂律的大部分有趣行为发生的地方。可能需要使用其他指标，如Kuiper或AndersonDarling，它们在测量分布之间的距离时给予尾部额外的权重，用户可以根据实际需要进行调整。

>fit = powerlaw.Fit(data, xmin_distance = ‘D’)
>fit = powerlaw.Fit(data, xmin_distance = ‘V’)
>fit = powerlaw.Fit(data, xmin_distance = ‘Asquare’)

离散与连续数据

为了适合幂律和其他分布的连续形式，数据默认是连续的。然而，许多数据是离散的。离散分布的概率分布不能精确地与连续分布相提并论。离散（整数）分布通过适当的规范化，可以在初始化时指定。不过离散形式的概率分布通常比连续形式更难计算，因此某些操作运行可能会比较慢。

>fit = powerlaw.Fit(data, xmin = 230.0)
>fit.discrete
False
>fit = powerlaw.Fit(data, xmin = 230.0, discrete = True)
>fit.discrete
True

与其他分布比较

从创建的Fit对象中，用户可以很容易地访问评估长尾分布所需的所有统计分析。Fit对象也能适用于其他可能的分布形式。每个分布都有适合该分布的最佳参数，可通过参数名称或更通用的“parameter1”访问。

>fit.power_law.alpha
2.273
>fit.power_law.parameter1
2.273
>fit.power_law.parameter1_name
>fit.lognormal.mu
0.154
>fit.lognormal.parameter1_name
‘mu’
>fit.lognormal.parameter2_name
‘sigma’
>fit.lognormal.parameter3_name = = NoneTrue

另外，得到参数之后还需要评价哪个分布优度更好。每个分布的拟合优度可以考虑单独或通过相互比较得到（这里使用KS 检验和对数似然比来确定）。比如下方代码直接对比幂律分布和指数分布的拟合效果。返回得到的R是两个候选分布之间的对数似然比，如果数据更有可能在第一个分布中，这个数字将为正，如果数据更有可能在第二个分布中，这个数字将为负。p是显著性指标。

>R, p = fit.distribution_compare(‘power_law’, ‘exponential’,
normalized_ratio = True)
>print R, p
16.384  0.024   #幂律分布比指数分布更优

一般认为，指数分布是长尾分布中绝对最小备选。因此，如果当下结果并不比指数分布更适合，就根本没有理由考虑其他的长尾分布，更别说幂律分布了。fit对象的fit.supported_distribution中支持的分布列表包括：幂律、截断幂律、指数、对数正态、威布尔和伽马分布。distribution_compare可以使用这些分布中的任何一个。

参考资料：

1. 该论文中的拟合数据和代码均可见：https://github.com/jeffalstott/powerlaw

2. 关于幂律分布和复杂网络的相关联系可以参考：集智之前关于幂律分布的一些长文链接。

4. 关于幂律分布拟合及检验原理可以参考集智学园陈清华老师授课视频：

https://campus.swarma.org/course/1747 。

课程推荐

课程地址：https://campus.swarma.org/course/647

推荐阅读

点击“阅读原文”，进入集智学园听课