Generalized Hyperbolic Distribution#

The Generalized Hyperbolic Distribution is defined as the normal variance-mean mixture with Generalized Inverse Gaussian distribution as the mixing distribution. The “hyperbolic” characterization refers to the fact that the shape of the log-probability distribution can be described as a hyperbola. Hyperbolic distributions are sometime referred to as semi-fat tailed because their probability density decrease slower than “sub-hyperbolic” distributions (e.g. normal distribution, whose log-probability decreases quadratically), but faster than other “extreme value” distributions (e.g. pareto distribution, whose log-probability decreases logarithmically).

Functions#

Different parameterizations exist in the literature; SciPy implements the “4th parametrization” in Prause (1999).

f(x,p,a,b)=(a2b2)p/22πap0.5Kp(a2b2)ebx×Kp1/2(a1+x2)(1+x2)1/2p

for:

  • x,p(;)

  • |b|<a if p0

  • |b|a if p<0

  • Kp(.) denotes the modified Bessel function of the second kind and order p (scipy.special.kn)

The probability density above is defined in the “standardized” form. To shift and/or scale the distribution use the loc and scale parameters. Specifically, f(x,p,a,b,loc,scale) is identically equivalent to 1scalef(y,p,a,b) with y=1scale(xloc).

This parameterization derives from the original (λ,α,β,δ,μ) parameterization in Barndorff (1978) by setting:

  • λ=p

  • α=aδ=α^δ

  • β=bδ=β^δ

  • δ=scale

  • μ=location

Random variates for the scipy.stats.genhyperbolic can be efficiently sampled from the above-mentioned normal variance-mean mixture where scipy.stats.geninvgauss is parametrized as GIG(p=p,b=α^2β^2,loc=location,scale=1α^2β^2) so that: GH(p,α^,β^)=β^GIG+GIGN(0,1)

The “generalized” characterization suggests the fact that this distribution is a superclass of several other probability distribution, for instance:

  • f(p=ν/2,a=0,b=0,loc=0,scale=ν) has a Student’s t-distribution (scipy.stats.t) with ν degrees of freedom.

  • f(p=1,a=α^,b=β^,loc=μ,scale=δ) has a Hyperbolic Distribution.

  • f(p=1/2,a=α^,b=β^,loc=μ,scale=δ) has a Normal Inverse Gaussian Distribution (scipy.stats.norminvgauss).

  • f(p=1,a=δ,b=0,loc=μ,scale=δ) has a Laplace Distribution (scipy.stats.laplace) for δ0

Examples#

It is useful to understand how the parameters affect the shape of the distribution. While it is fairly straightforward to interpret the meaning of b as the skewness, understanding the difference between a and p is not as obvious, as both affect the kurtosis of the distribution. a can be interpreted as the speed of the decay of the probability density (where a>1 the asymptotic decay is faster than loge and vice versa) or - equivalently - as the slope of the log-probability hyperbola asymptote (where a>1 decay is faster than |1| and vice versa). p can be seen as the width of the shoulders of the probability density distribution (where p<1 results in narrow shoulders and vice versa) or - equivalently - as the shape of the log-probability hyperbola, which is convex for p<1 and concave otherwise.

import numpy as np
from matplotlib import pyplot as plt
from scipy import stats

p, a, b, loc, scale = 1, 1, 0, 0, 1
x = np.linspace(-10, 10, 100)

# plot GH for different values of p
plt.figure(0)
plt.title("Generalized Hyperbolic | -10 < p < 10")
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        label = 'GH(p=1, a=1, b=0, loc=0, scale=1)')
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'red', alpha = 0.5, label='GH(p>1, a=1, b=0, loc=0, scale=1)')
[plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'red', alpha = 0.1) for p in np.linspace(1, 10, 10)]
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'blue', alpha = 0.5, label='GH(p<1, a=1, b=0, loc=0, scale=1)')
[plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'blue', alpha = 0.1) for p in np.linspace(-10, 1, 10)]
plt.plot(x, stats.norm.pdf(x, loc, scale), label = 'N(loc=0, scale=1)')
plt.plot(x, stats.laplace.pdf(x, loc, scale), label = 'Laplace(loc=0, scale=1)')
plt.plot(x, stats.pareto.pdf(x+1, 1, loc, scale), label = 'Pareto(a=1, loc=0, scale=1)')
plt.ylim(1e-15, 1e2)
plt.yscale('log')
plt.legend(bbox_to_anchor=(1.1, 1))
plt.subplots_adjust(right=0.5)

# plot GH for different values of a
plt.figure(1)
plt.title("Generalized Hyperbolic | 0 < a < 10")
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        label = 'GH(p=1, a=1, b=0, loc=0, scale=1)')
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'blue', alpha = 0.5, label='GH(p=1, a>1, b=0, loc=0, scale=1)')
[plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'blue', alpha = 0.1) for a in np.linspace(1, 10, 10)]
plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'red', alpha = 0.5, label='GH(p=1, 0<a<1, b=0, loc=0, scale=1)')
[plt.plot(x, stats.genhyperbolic.pdf(x, p, a, b, loc, scale),
        color = 'red', alpha = 0.1) for a in np.linspace(0, 1, 10)]
plt.plot(x, stats.norm.pdf(x, loc, scale),  label = 'N(loc=0, scale=1)')
plt.plot(x, stats.laplace.pdf(x, loc, scale), label = 'Laplace(loc=0, scale=1)')
plt.plot(x, stats.pareto.pdf(x+1, 1, loc, scale), label = 'Pareto(a=1, loc=0, scale=1)')
plt.ylim(1e-15, 1e2)
plt.yscale('log')
plt.legend(bbox_to_anchor=(1.1, 1))
plt.subplots_adjust(right=0.5)

plt.show()

References#

Implementation: scipy.stats.genhyperbolic