文章目录
- 介绍
- 代码实例
介绍
数学建模中很多时候，我们有某个随机变量 X X X 的若干样本 X 1 , X 2 , ⋯ , X n X_1,X_2,\cdots,X_n X1,X2,⋯,Xn，想要还原随机变量 X X X 的概率密度函数 f ( x ) f(x) f(x)。诚然，高斯核密度估计可以做到这一点，这里我们介绍一个新的 Python 库 Fitter。
Fitter 的原理是，用指定的某些分布去匹配连续性随机变量 X X X 的观测值，对于每一种分布都拟合出响应的参数、计算出均方误差。一般选择均方误差最小的分布作为拟合结果。核心方法Fitter的函数原型为：
Fitter(data, xmin=None, xmax=None, bins=100, distributions=None, verbose=True, timeout=10)
其中：
data是已知的 X X X 的观测值集合，或者说是样本点，即上面的 X 1 , ⋯ , X n X_1,\cdots,X_n X1,⋯,Xn。xmin表示忽略观测值data中小于xmin的数据。xmax表示忽略观测值data中大于xmax的数据。bins是累计直方图的组数（直方图里有几个矩形），值越大拟合曲线越平滑。verbose表示是否显示 Fitter 拟合过程中的一些提示信息，设置verbose=False将会使你的控制台更干净。timeout表示对于一种分布尝试的最长时间。因为各种各样的原因，给出的数据data可能不适合用某种分布去拟合，当 Fitter 在这种分布上拟合的时间超过timeout秒时，将会自动放弃拟合这种分布。distributions是一个列表，用来告诉 Fitter你想尝试哪些概率密度函数（比如正态分布， t t t 分布， β \beta β 分布等）。默认情况下会尝试所有的 scipy 分布（ 106 106 106 种）。常用分布distributions=['norm', 't', 'laplace', 'cauchy', 'chi2', 'expon', 'exponpow', 'gamma', ' lognorm', 'uniform']。（比如norm是正态分布，t是 t t t 分布，laplace是拉普拉斯分布，cauchy是柯西分布……）
参考文献：Python fitter包：拟合数据样本的分布_fitter的80多个分布有哪些-CSDN博客
代码实例
from scipy import stats
import numpy as np
import matplotlib.pyplot as plt
from fitter import Fitter# 直接生成 2000 个 beta 分布的样本点
data = np.array(np.random.beta(a=2, b=5, size=2000))# Euclid 字体是 latex 公式使用的字体
plt.rcParams['font.family'] = 'Euclid'# 尝试使用 distributions 里的分布去拟合样本点
f = Fitter(data, distributions=['norm', 't', 'laplace','beta','lognorm'])f.fit()print(f.summary()) # 返回排序好的分布拟合质量（拟合效果从好到坏）,并绘制数据分布和Nbest分布
"""sumsquare_error aic ... ks_statistic ks_pvalue
beta 19.103771 164.735384 ... 0.016081 6.791340e-01
lognorm 22.106143 143.794008 ... 0.026787 1.113541e-01
t 30.674133 248.908752 ... 0.060242 9.434693e-07
norm 30.674533 246.915057 ... 0.060244 9.426322e-07
laplace 52.044523 184.287028 ... 0.068199 1.562740e-08
可以看到用 beta 分布去拟合，均方误差最小。
"""print(f.df_errors) # 返回这些分布的拟合质量（均方根误差的和）
"""
sumsquare_error aic ... ks_statistic ks_pvalue
norm 30.674533 246.915057 ... 0.060244 9.426322e-07
laplace 52.044523 184.287028 ... 0.068199 1.562740e-08
beta 19.103771 164.735384 ... 0.016081 6.791340e-01
lognorm 22.106143 143.794008 ... 0.026787 1.113541e-01
t 30.674133 248.908752 ... 0.060242 9.434693e-07
"""print(f.fitted_param) # 返回拟合分布的参数
"""
{'norm': (0.2796987317013659, 0.15616306964918625),
'laplace': (0.26198195036801175, 0.12567921087793762),
'beta': (1.9506594116888034, 4.9951825006926605, 0.0024740416891275907, 0.9862799904927366),
'lognorm': (0.31353430150353023, -0.21970151525943354, 0.47572327635046574),
't': (850256.2207528697, 0.2796983453401272, 0.15616418590110978)}
"""print(f.fitted_pdf) # 使用最适合数据分布的分布参数生成的概率密度
"""
结果是返回从分布到 array 的字典。array 是一系列概率密度函数的取值，用来作概率密度曲线。
"""print(f.get_best(method='sumsquare_error')) # 返回最佳拟合分布及其参数
"""
{'beta': {'a': 1.9506594116888034, 'b': 4.9951825006926605, 'loc': 0.0024740416891275907, 'scale': 0.9862799904927366}
}
"""f.hist() # 绘制组数=bins的标准化直方图
# 下面这个方法绘制分布的概率密度函数，按照方差从小到大取前 Nbest 个。
# 因为前面的 f.summary() 也会画图，所以就不启用这个方法。
# 顺带一提，前面的 f.summary() 也可以设置 Nbest 参数。
# f.plot_pdf(names=None, Nbest=3, lw=2)
plt.show()
画出来的图如下所示， β \beta β 分布是均方误差最小的分布，可以认为是最合适的分布。（因为数据就是按照 β \beta β 分布生成的，所以这个结果很合理。）同时，对数正态分布也比较合适， t t t 分布次之。
注意：
Fitter方法的distributions参数决定了需要拟合哪些分布，默认拟合里面所有的 106 106 106 种分布。默认情况下会跑比较长的时间，而且往往给出最佳的分布是一个你从来没见过的分布。不建议使用默认情况。f.summary()默认参数Nbest=5，给出 5 5 5 个最好拟合的结果。如果想要所有分布的结果，请设置Nbest=len(f.distributions)，即将f.summary()改为f.summary(Nbest=len(f.distributions))。
