The following code is adapted from an example invented by @KieranDelamotte.
We try to fit a Gumbel distribution on a sample that is actually obtained from a mixture of a Gumbel and a Fréchet distribution.
GumbelFactory does not manage to find a sensible set of parameters to try to fit the sample.
WeibullMaxFactory, on the other hand, manages to find a Weibull distribution that kind of fits, but it has crazy parameters… because it is actually very close to a Gumbel. In fact, when put into a GEV, the GEV \xi parameter is 0: it is a Gumbel distribution.
The bottom line is that WeibullMaxFactory manages to find a better Gumbel than GumbelFactory.
The difference is that WeibullMaxFactory maximizes (by default) the likelihood whereas GumbelFactory fits the first and second order moments to the empirical moments of the sample. Likelihood maximization is more robust to (this particular) distribution misspecification than the method of moments.
# Try to fit a WeibullMax on the sample.
# WeibullMaxFactory maximizes the likelihood.
weibullmax = ot.WeibullMaxFactory().build(sample)
print("WeibullMaxFactory result: ", weibullmax)
# Express the estimated WeibullMax as a GEV
gev = ot.GeneralizedExtremeValue(weibullmax)
print("GEV actual distribution: ", gev.getActualDistribution())
# Recompute the GEV actual distribution
# This WeibullMax turns out to be a Gumbel
# because the GEV xi parameter is zero.
print("GEV xi parameter = ", gev.getXi())
gev.setParameter(gev.getParameter())
print("More sensible GEV actual distribution: ", gev.getActualDistribution())
# Try to directly estimate a Gumbel distribution.
# GumbelFactory uses the method of moments.
gumbel = ot.GumbelFactory().build(sample)
print("GumbelFactory result: ", gumbel)
graph = ot.HistogramFactory().build(sample).drawPDF()
graph.add(gev.drawPDF())
graph.add(gumbel.drawPDF())
graph.setColors(["green", "red", "blue"])
graph.setLegends(["Sample", "MLE", "Moments"])
graph.setTitle("Estimated Gumbel distribution")
view = viewer.View(graph)
axes = view.getAxes()
_ = axes[0].set_xlim(-20.0, 20.0)
Output:
WeibullMaxFactory result:
WeibullMax(beta = 4.25531e+06, alpha = 1.45411e+06, gamma = 4.25531e+06)
GEV actual distribution:
WeibullMax(beta = 4.25531e+06, alpha = 1.45411e+06, gamma = 4.25531e+06)
GEV xi parameter = -6.877037385642863e-07
More sensible GEV actual distribution:
Gumbel(beta = 2.92639, gamma = 2.5121)
GumbelFactory result: Gumbel(beta = 15.1936, gamma = -3.88787)
