直觉和示例代码
这是一系列关于< Power Laws和Fat Tails的第三篇文章。在之前的文章中，我们探讨了如何从实证数据中检测到幂律分布。虽然这种技术很有用，但是fat tails并不仅仅是将数据拟合到幂律分布中。在这篇文章中，我将解析4种我们可以量化fat tails的方法，并分享分析真实数据的示例Python代码。
注：如果您对Power Law分布或Fat Tail等术语不熟悉，请查看这篇文章作为入门。
在本系列的第一篇文章中，我们介绍了“fat tails”的概念，这描述了罕见事件对分布的总体统计数据的影响程度。我们通过Pareto分布看到了fat tails的一个极端例子，例如，80%的销售由20%的客户产生（而50%的销售仅由1%的客户产生）。
Pareto、Power Laws和Fat Tails
统计学没有教给你的东西
towardsdatascience.com
尽管Pareto（以及更一般的幂律）分布为我们提供了一个明显的fat tails示例，但这是一个更一般的概念，范围从具有薄尾（如高斯分布）到具有非常fat-tailed（如 Pareto 80-20）的谱。
这种fat-tailedness的观点为我们提供了比简单将其标记为Power Law的灵活和精确的分类数据的方式。然而，这引发了一个问题：我们如何定义fat-tailedness？
量化Fat Tails的4种方法
虽然没有“真实”的fat-tailedness度量，但在实践中我们可以使用一些启发式（即经验法则）来量化fat-tailed数据。在这里，我们介绍4种这样的启发式方法。我们首先从概念上介绍每种技术，然后深入示例Python代码。
启发式1：幂律尾指数
最fat的fat tails出现在幂律分布中，其中幂律的尾部指数（即α）越小，尾部越fat，如下图所示。
观察到较小的尾指数意味着更fat的tails自然地激发了我们使用α来量化fat tails。在实践中，这归结为将幂律分布拟合到给定的数据集并提取估计的α值。
虽然这是一种直接的方法，但它有一个明显的限制。也就是说，当处理与幂律差异较大的数据时，该方法会失效。
启发式二：峰度（即非高斯性）
与长尾相反的是“瘦尾”（即“罕见事件如此罕见，以至于可以忽略不计”）。一个瘦尾示例是所钟爱的高斯分布，其中一个6 sigma离均值的事件发生的概率约为千亿分之一。
这启发我们可通过衡量数据的“非高斯性”来量化长尾。我们可以通过所谓的非高斯性度量来实现这一点。虽然我们可以设计很多这样的度量，但最常用的是峰度，其定义如下所示。
![根据引用详见[1]和[2]定义的峰度图。图像由作者提供。](https://ai.miximages.com/miro.medium.com/v2/resize:fit:640/format:webp/1*DcgMTRkGVtoE4VdpKhh98A.png "用Python量化Fat Tails的四种方法 四海 第4张-四海吧")
峰度受到远离中心（即尾部）的值的影响。因此，峰度越大，尾部越厚。
这个度量在所有矩尺度是有限的情况下效果良好[3]。然而，它有一个主要限制，即峰度不能对某些分布进行定义，例如α=<4的帕累托分布，这使得它在许多长尾数据上无用。
启发式三：对数正态分布的σ
在本系列的之前文章中，我们讨论了由以下概率密度函数定义的对数正态分布。
![对数正态分布的概率密度函数[4]。图像由作者提供。](https://ai.miximages.com/miro.medium.com/v2/resize:fit:640/format:webp/1*Xb79fUujEQHf1AWi4yU5_A.png "用Python量化Fat Tails的四种方法 四海 第5张-四海吧")
正如我们之前所见，对于较小的σ，该分布可能类似于高斯分布，但对于较大的σ则类似于帕累托分布。这自然提供了另一种量化长尾的方式，其中σ越大，尾部越厚。
我们可以像启发式一一样得到这个度量。即，我们将一个对数正态分布适应到我们的数据中并提取拟合的σ值。尽管这是一个简单的过程，但它（就像启发式一一样）在对数正态拟合不能很好地解释底层数据时会失效。
启发式四：塔勒布的κ
先前的启发式（H）从特定的分布开始（即H1 – 幂律和H3 – 对数正态）。然而，在实践中，我们的数据很少精确遵循任何特定的分布。
此外，使用这些度量进行比较时，评估遵循定性不同分布的两个变量可能存在问题。例如，使用幂律的尾指数来比较类似帕累托和类似高斯的数据可能没有太大意义，因为幂律无法很好地适应类似高斯的数据。
这促使我们使用非分布特定的长尾度量。塔勒布在引用[3]中提出了一种这样的度量。提出的指标（κ）适用于具有有限均值的单峰数据，并且取值介于0和1之间，其中0表示数据是极度瘦尾，1表示数据是极度长尾。其定义如下所示。
![Taleb's κ度量定义[3]。图片作者：笔者。](https://ai.miximages.com/miro.medium.com/v2/resize:fit:640/format:webp/1*2A-bJKuEMpt0kIVX6hILig.png "用Python量化Fat Tails的四种方法 四海 第6张-四海吧")
这个度量比较两个样本（比如，Sₙ₀和Sₙ），其中Sₙ是从特定分布中抽取的n个样本的总和。例如，如果我们评估一个高斯分布并选择n=100，我们会从高斯分布中抽取100个样本，并将它们全部相加以创建S₁₀₀。
上述表达式中的M(n)表示平均绝对偏差（mean absolute deviation），根据下面的方程定义。这个度量反映了平均值周围的离散度，通常比标准差更具鲁棒性[3][5]。
![从κ方程中定义的平均绝对差[3]。图片作者：笔者。](https://ai.miximages.com/miro.medium.com/v2/resize:fit:640/format:webp/1*X47ZBi7QKxpX_xjgnQ_d1w.png "用Python量化Fat Tails的四种方法 四海 第7张-四海吧")
为了简化，我们可以选择n₀=1，得到下面的表达式。
这里的关键术语是M(n)/M(1)，其中M(n)量化了n个样本的平均值周围的离散度（某个分布的样本）。
对于尾部纤细的分布，M(30)将与M(1)相对较接近，因为数据通常接近均值。因此，M(30)/M(1) ~ 1。
然而，对于尾部较厚的数据，M(30)将比M(1)大得多。因此，M(30)/M(1) >> 1。下面的示例中，左图显示了高斯分布和右图显示了Pareto分布的离散度如何随n的增加而扩展。
因此，对于尾部较厚的数据，κ方程中的分母将大于分子，使得右侧的第二项较小，最终导致κ的值较大。
如果这对你来说有点太数学了，下面是要点：大κ表示尾部较厚，小κ表示尾部较纤细。
示例代码：量化（现实世界的）社交媒体数据的尾部厚度
掌握这些概念后，让我们看看在实践中使用这些启发式方法的样子。在这里，我们将使用上述每种方法来分析本系列的前一篇文章中的相同数据。
这些数据来自我的社交媒体账户，包括每月在VoAGI上获得的关注者数量、每个YouTube视频的收益以及LinkedIn上的每日印象。数据和代码可以在GitHub仓库上免费获取。
YouTube-Blog/power laws/3-测量尾部的重尾现象在主文中 · ShawhinT/YouTube-Blog
补充YouTube视频和VoAGI上的博客帖子的代码。- YouTube-Blog/power laws/3-测量尾部的重尾现象在主文中…
github.com
我们首先导入一些有用的库。
import matplotlib.pyplot as pltimport pandas as pdimport numpy as npimport powerlawfrom scipy.stats import kurtosis接下来，我们将加载每个数据集并将它们存储在一个字典中。
filename_list = ['VoAGI-followers', 'YT-earnings', 'LI-impressions']df_dict = {}for filename in filename_list: df = pd.read_csv('data/'+filename+'.csv') df = df.set_index(df.columns[0]) # 设置索引 df_dict[filename] = df此时，查看数据总是一个好主意。我们可以通过绘制直方图和打印每个数据集的前5条记录来实现。
for filename in filename_list: df = df_dict[filename] # 绘制直方图（以下函数在GitHub上的笔记本中定义） plot_histograms(df.iloc[:,0][df.iloc[:,0]>0], filename, filename.split('-')[1]) plt.savefig("images/"+filename+"_histograms.png") # 打印前5条记录 print("百分比最高的5条记录") print((df.iloc[:,0]/df.iloc[:,0].sum()).sort_values(ascending=False)[:5]) print("")
根据上述直方图，每个数据集都在某种程度上显示了重尾现象。让我们看看百分比最高的5条记录，以进一步了解情况。
从这个视角来看，VoAGI关注者似乎是最重尾的，有60%的关注者来自于仅有的2个月。YouTube收入也明显呈现重尾现象，约60%的收入来自于仅有的4个视频。LinkedIn印象似乎是最不重尾的。
虽然我们可以通过查看数据来定性地了解尾部重尾现象，但通过我们的4个启发式方法，让我们将这个问题定量化。
启发式方法1：幂律尾指数
要为每个数据集获取α值，我们可以使用powerlaw库，就像我们在上一篇文章中所做的那样。在下面的代码块中完成此操作，我们使用循环对每个数据集执行拟合并打印参数估计结果。
for filename in filename_list: df = df_dict[filename] # 执行Power Law拟合 results = powerlaw.Fit(df.iloc[:,0]) # 打印结果 print("") print(filename) print("-"*len(filename)) print("Power Law拟合") print("a = " + str(results.power_law.alpha-1)) print("xmin = " + str(results.power_law.xmin)) print("")
上面的结果与我们的定性评估一致，VoAGI的关注者是最疲软的，其次是YouTube的收入和LinkedIn的印象（记住，较小的α意味着尾部更粗）。
启发式2：峰度
计算峰度的一种简单方法是使用现成的实现。在这里，我使用Scipy库，并以与之前类似的方式打印结果。
for filename in filename_list: df = df_dict[filename] # 打印结果 print(filename) print("-"*len(filename)) print("峰度 = " + str(kurtosis(df.iloc[:,0], fisher=True))) print("")
峰度给出的结果与启发式1不同。根据此度量，尾部的粗细排序如下：LinkedIn > VoAGI > YouTube。
然而，这些结果应该带有一些怀疑态度。正如我们在上面的Power Law拟合中看到的，所有3个数据集都符合α<4的幂律分布，这意味着峰度是无穷大的。因此，虽然计算返回了一个值，但对这些数字持怀疑态度可能是明智的。
启发式3：对数正态分布的σ
我们可以再次使用powerlaw库，与启发式1类似地获得σ估计。以下是它的样子。
for filename in filename_list: df = df_dict[filename] # 执行Power Law拟合 results = powerlaw.Fit(df.iloc[:,0]) # 打印结果 print("") print(filename) print("-"*len(filename)) print("对数正态分布拟合") print("mu = " + str(results.lognormal.mu)) print("sigma = " + str(results.lognormal.sigma)) print("")
从上面的σ值可以看出，所有拟合都暗示数据是粗尾的，其中VoAGI的关注者和LinkedIn的印象具有类似的σ估计值。另一方面，YouTube的收入具有显着较大的σ值，暗示有（更）粗尾。
然而，令人怀疑的是，拟合估计了一个负的μ，这可能意味着对数正态分布可能无法很好地解释数据。
启发式4：塔勒布的κ
由于我找不到现成的用于计算κ的Python实现（我也没很用心地找），所以这个计算需要一些额外的步骤。即我们需要定义3个辅助函数，如下所示。
def mean_abs_deviation(S): """ 计算输入样本S的平均绝对偏差 """ M = np.mean(np.abs(S - np.mean(S))) return Mdef generate_n_sample(X,n): """ 从数组X生成n个长度为len(X)的随机样本的函数 """ # 初始化样本 S_n=0 for i in range(n): # 从X中随机选择长度为len(X)的观测，将其添加到样本中 S_n = S_n + X[np.random.randint(len(X), size=int(np.round(len(X))))] return S_ndef kappa(X,n): """ 根据此处描述的Taleb的kappa度量从n0=1计算得到：https://arxiv.org/abs/1802.05495 注意：K_1n = kappa(1,n) = 2 - ((log(n)-log(1))/log(M_n/M_1))，其中M_n表示n个随机样本的平均绝对偏差 """ S_1 = X S_n = generate_n_sample(X,n) M_1 = mean_abs_deviation(S_1) M_n = mean_abs_deviation(S_n) K_1n = 2 - (np.log(n)/np.log(M_n/M_1)) return K_1n第一个函数mean_abs_deviation()计算了前面定义的平均绝对偏差。
接下来，我们需要一种方法来生成并求和n个来自我们的经验数据的样本。在这里，我采用了一种简单的方法，随机采样一个输入数组（X）n次并将样本相加。
最后，我将mean_abs_deviation(S)和generate_n_sample(X,n)结合起来，实现了之前定义的κ计算，并为每个数据集计算了κ值。
n = 100 # 在kappa计算中包含的样本数量for filename in filename_list: df = df_dict[filename] # 打印结果 print(filename) print("-"*len(filename)) print("kappa_1n = " + str(kappa(df.iloc[:,0].to_numpy(), n))) print("")
上述结果给了我们另一个故事。然而，考虑到此计算的隐式随机性（请回顾generate_n_sample()的定义）以及我们处理的数据具有fat tails特征，点估计（即仅运行一次计算）是不可靠的。
因此，我运行了相同的计算1000次，并打印了每个数据集的平均κ(1,100)。
num_runs = 1_000kappa_dict = {}for filename in filename_list: df = df_dict[filename] kappa_list = [] for i in range(num_runs): kappa_list.append(kappa(df.iloc[:,0].to_numpy(), n)) kappa_dict[filename] = np.array(kappa_list) print(filename) print("-"*len(filename)) print("mean kappa_1n = " + str(np.mean(kappa_dict[filename]))) print("")
这些更稳定的结果表明，VoAGI关注者是最具fat-tailed特征的，其次是LinkedIn Impressions和YouTube收入。
注意：可以将这些值与参考文献[3]中的第III表进行比较，以更好地了解每个κ值。即，这些值与α在2到3之间的Pareto分布相当。
尽管每个启发式方法都讲述了略有不同的故事，但所有迹象都指向了VoAGI关注者的fat-tailed特征是这3个数据集中最为明显的。
结论
将数据标记为长尾型（或非长尾型）可能会很诱人，但长尾性存在于一个光谱中。在这里，我们分解了四个用于量化数据长尾性的启发式方法。
尽管每种方法都有其局限性，但它们为从业者提供了比较实证数据长尾性的定量方式。
资源
支持: 给我买杯咖啡 ☕️
数据企业家
数据领域的企业家社区。👉 加入Discord！
VoAGI.com
[1] Scipy的峰度
[2] Scipy的Moment函数
[3] arXiv:1802.05495 [stat.ME]
[4] https://en.wikipedia.org/wiki/Log-normal_distribution
[5] Pham-Gia, T., & Hung, T. (2001). The mean and median absolute deviations. Mathematical and Computer Modelling, 34(7–8), 921–936. https://doi.org/10.1016/S0895-7177(01)00109-1