Python

Question

我是 Python 的新手，我在互联网上四处寻找，但找不到任何可以帮助我解决问题的逻辑。

我在图表中有降水值，现在我需要根据图表中的这些值拟合 GEV 分布。每个值等于从 1974 年到 2017 年的一年中的最大值（因此总共有 43 个值）。

这些是值：

max_precip = [9.4, 38.0, 12.5, 35.3, 17.6, 12.9, 12.4, 19.6, 15.0, 13.2, 12.3, 16.9, 16.9, 29.4, 13.6, 11.1, 8.0, 16.6, 12.0, 13.1, 9.1, 9.7, 21.0, 11.2, 14.4, 18.8, 14.0, 19.9, 12.4, 10.8, 21.6, 15.4, 17.4, 14.8, 22.7, 11.5, 10.5, 11.8, 12.4, 16.6, 11.7, 12.9, 17.8]

我发现我需要使用gev.fit，所以我想使用以下：

t = np.linspace(1,43,43)
fit = gev.fit(max_precip,loc=3)
pdf = gev.pdf(t, *fit)
plt.plot(t,pdf)
plt.plot(t, max_precip, "o")

但这只会在图表中打印出 max_precip 的点，而不是 GEV 分布。

有人可以帮助我吗？抱歉，如果这个问题已经有人问过，我找不到类似的问题。

我使用了这些导入：

import csv
import matplotlib.pyplot as plt
import numpy as np

from dateutil.rrule import rrule, YEARLY
import datetime
from matplotlib.dates import DateFormatter
from scipy.stats import genextreme as gev
from scipy.stats import genpareto as gpd
from scipy.optimize import minimize

Answer 1

我试过拟合你的数据

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import genextreme as gev

def main(rvs):
    shape, loc, scale = gev.fit(rvs)
    return shape, loc, scale

if __name__ == '__main__':
    rvs = [9.4, 38.0, 12.5, 35.3, 17.6, 12.9, 12.4, 19.6, 15.0, 13.2, 12.3, 16.9, 16.9, 29.4, 13.6, 11.1, 8.0, 16.6, 12.0, 13.1, 9.1, 9.7, 21.0, 11.2, 14.4, 18.8, 14.0, 19.9, 12.4, 10.8, 21.6, 15.4, 17.4, 14.8, 22.7, 11.5, 10.5, 11.8, 12.4, 16.6, 11.7, 12.9, 17.8]

    shape, loc, scale = main(rvs)

    print(shape)
    print(loc)
    print(scale)

    l = loc + scale / shape

    xx = np.linspace(l+0.00001, l+0.00001+35, num=71)
    yy = gev.pdf(xx, shape, loc, scale)

    hist, bins = np.histogram(rvs, bins=12, range=(-0.5, 23.5), density=True)
    plt.bar(bins[:-1], hist, width = 2, align='edge')

    plt.plot(xx, yy, 'ro')
    plt.show()

但我得到的是

-0.21989526255575445
12.749780017954315
3.449061347316184

对于 shape、loc 和 scale。如果看GEV distribution as defined in scipy，shape为负数时，有效区间为[loc + scale/shape...+infinity]。我计算了后一个值，它等于

-2.935417290135696

应该可以...

Python3, 蟒蛇, scipy 1.1, Windows 10 64 位

更新

好的，我已经更新了代码并添加了绘图，看起来有点合理。这是你在找什么吗？基本上，技巧是对其进行直方图绘制并绘制与 PDF

重叠的密度箱

Answer 2

出于好奇，我尝试了 OpenTURNS

中可用的 GeneralizedExtremeValue Factory (GEV)

import openturns as ot
sample = ot.Sample([[p] for p in max_precip])

gev = ot.GeneralizedExtremeValueFactory().buildAsGeneralizedExtremeValue(sample)
print (gev)
>>> GeneralizedExtremeValue(mu=12.7497, sigma=3.44903, xi=0.219894)

我可以确认它给出了相同的结果。

Python - 根据这些值拟合 GEV 分布

Python - Fitting a GEV distribution from these values

plot

distribution