pdqr module

All distributions are parameterized by the mean, standard deviation (sigma), skew and shape parameters. For the Generalized Hyperbolic distribution, there is an additional parameter lamda, representing a real parameter in the Generalized Inverse Gaussian (GIG) distribution of the Normal mean variance mixture representation of the distribution. For each distribution (*), there is a function for the cumulative distribution (p*), density (d*), quantile (q*) and random sampling (r*). Additional functions for calculating the skewness and kurtosis given the skew and shape parameters are also available.

Normal distribution

The Normal Distribution is a spherical distribution described completely by it first two moments, the mean and variance. Formally, the random variable \(x\) is said to be normally distributed with mean \(\mu\) and variance \(\sigma^2\) (both of which may be time varying), with density given by,

\[ \begin{equation} f\left( x \right) = \frac{{{e^{\frac{{ - 0.5{{\left( {x - \mu } \right)}^2}}} {{{\sigma ^2}}}}}}} {{\sigma \sqrt {2\pi } }}. \end{equation} \]

where \(\mu \in \mathbb{R}\) and \(\sigma \in \mathbb R_{> 0}\).

pdqr.dnorm(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, log: bool = False) → NDArray[float64]

(Normal) Probability Density Function (dnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• log – whether to return the log density

Return type:

a numpy array

pdqr.pnorm(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, lower_tail: bool = True) → NDArray[float64]

(Normal) Cumulative Probability Function (pnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qnorm(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, lower_tail: bool = True) → NDArray[float64]

(Normal) Quantile Function (qnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rnorm(n: int = 1, mu: float = 0, sigma: float = 1, seed: int | None = None) → NDArray[float64]

(Normal) Random Number Function (rnorm)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Student’s T distribution

The Student’s t distribution is described completely by a shape parameter \(\nu\), but for standardization we proceed by using its 3 parameter representation as follows:

\[ \begin{equation} f\left( x \right) = \frac{{\Gamma \left( {\frac{{\nu + 1}} {2}} \right)}} {{\sqrt {\beta \nu \pi } \Gamma \left( {\frac{\nu } {2}} \right)}}{\left( {1 + \frac{{{{\left( {x - \alpha } \right)}^2}}} {{\beta \nu }}} \right)^{ - \left( {\frac{{\nu + 1}} {2}} \right)}} \end{equation} \]

where \(\alpha\), \(\beta\), and \(\nu\) are the location, scale and shape parameters respectively, and \(\Gamma\) is the Gamma function. Similar to the GED distribution described later, this is a unimodal and symmetric distribution where the location parameter \(\alpha\) is the mean (and mode) of the distribution while the variance is:

\[ \begin{equation} Var\left( x \right) = \frac{{\beta \nu }}{{\left( {\nu - 2} \right)}}. \end{equation} \]

For the purposes of standardization we require that:

\[\begin{split} \begin{gathered} Var(x) = \frac{{\beta \nu }} {{\left( {\nu - 2} \right)}} = 1 \hfill \\ \therefore \beta = \frac{{\nu - 2}} {\nu } \hfill \\ \end{gathered} \end{split}\]

Substituting \(\frac{(\nu- 2)}{\nu }\) we obtain the standardized Student’s distribution:

\[ \begin{equation} f\left( {\frac{{x - \mu }}{\sigma }} \right) = \frac{1} {\sigma }f\left( z \right) = \frac{1} {\sigma }\frac{{\Gamma \left( {\frac{{\nu + 1}} {2}} \right)}}{{\sqrt {\left( {\nu - 2} \right)\pi } \Gamma \left( {\frac{\nu } {2}} \right)}}{\left( {1 + \frac{{{z^2}}} {{\left( {\nu - 2} \right)}}} \right)^{ - \left( {\frac{{\nu + 1}} {2}} \right)}}. \end{equation} \]

The Student distribution has zero skewness and excess kurtosis equal to \(6/(\nu - 4)\) for \(\nu > 4\).

pdqr.dstd(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, log: bool = False) → NDArray[float64]

(Student) Probability Density Function (dstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• log – whether to return the log density

Return type:

a numpy array

pdqr.pstd(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, lower_tail: bool = True) → NDArray[float64]

(Student) Cumulative Probability Function (pstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qstd(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, lower_tail: bool = True) → NDArray[float64]

(Student) Quantile Function (qstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rstd(n: int = 1, mu: float = 0, sigma: float = 1, shape: float = 5, seed: int | None = None) → NDArray[float64]

(Student) Random Number Function (rstd)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Generalized Error distribution

The Generalized Error Distribution is a 3 parameter distribution belonging to the exponential family with conditional density given by,

\[ f\left( x \right) = \frac{{\kappa {e^{ - 0.5{{\left| {\frac{{x - \alpha }} {\beta }} \right|}^\kappa }}}}} {{{2^{1 + {\kappa ^{ - 1}}}}\beta \Gamma \left( {{\kappa ^{ - 1}}} \right)}} \]

with \(\alpha\), \(\beta\) and \(\kappa\) representing the location, scale and shape parameters. Since the distribution is symmetric and unimodal the location parameter is also the mode, median and mean of the distribution (i.e. \(\mu\)). By symmetry, all odd moments beyond the mean are zero. The variance and kurtosis are given by,

\[\begin{split} \begin{gathered} Var\left( x \right) = {\beta ^2}{2^{2/\kappa }}\frac{{\Gamma \left( {3{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {{\kappa ^{ - 1}}} \right)}} \hfill \\ Ku\left( x \right) = \frac{{\Gamma \left( {5{\kappa ^{ - 1}}} \right)\Gamma \left( {{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {3{\kappa ^{ - 1}}} \right)\Gamma \left( {3{\kappa ^{ - 1}}} \right)}} \hfill \\ \end{gathered} \end{split}\]

As \(\kappa\) decreases the density gets flatter and flatter while in the limit as \(\kappa \to \infty\), the distribution tends towards the uniform. Special cases are the Normal when \(\kappa=2\), the Laplace when \(\kappa=1\). Standardization is simple and involves re-scaling the density to have unit standard deviation:

\[\begin{split} \begin{gathered} Var\left( x \right) = {\beta ^2}{2^{2/\kappa }}\frac{{\Gamma \left( {3{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {{\kappa ^{ - 1}}} \right)}} = 1 \hfill \\ \therefore \beta = \sqrt {{2^{ - 2/\kappa }}\frac{{\Gamma \left( {{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {3{\kappa ^{ - 1}}} \right)}}} \hfill \\ \end{gathered} \end{split}\]

Finally, substituting into the scaled density of \(z\):

\[ f\left( {\frac{{x - \mu }} {\sigma }} \right) = \frac{1} {\sigma }f\left( z \right) = \frac{1} {\sigma }\frac{{\kappa {e^{ - 0.5{{\left| {\sqrt {{2^{ - 2/\kappa }}\frac{{\Gamma \left( {{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {3{\kappa ^{ - 1}}} \right)}}} z} \right|}^\kappa }}}}} {{\sqrt {{2^{ - 2/\kappa }}\frac{{\Gamma \left( {{\kappa ^{ - 1}}} \right)}} {{\Gamma \left( {3{\kappa ^{ - 1}}} \right)}}} {2^{1 + {\kappa ^{ - 1}}}}\Gamma \left( {{\kappa ^{ - 1}}} \right)}} \]

pdqr.dged(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, log: bool = False) → NDArray[float64]

(GED) Probability Density Function (dged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.pged(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, lower_tail: bool = True) → NDArray[float64]

(GED) Cumulative Probability Function (pged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qged(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, shape: ArrayLike = 5, lower_tail: bool = True) → NDArray[float64]

(GED) Quantile Function (qged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rged(n: int = 1, mu: float = 0, sigma: float = 1, shape: float = 5, seed: int | None = None) → NDArray[float64]

(GED) Random Number Function (rged)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• shape – the shape parameter (degrees of freedom)

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Skewed Distributions by Inverse Scale Factors

[Fernandez and Steel, 1998] proposed introducing skewness into unimodal and symmetric distributions by introducing inverse scale factors in the positive and negative real half lines.

Given a skew parameter, \(\xi\) (when \(\xi=1\), the distribution is symmetric), the density of a random variable z can be represented as:

\[ \begin{equation} f\left( {z|\xi } \right) = \frac{2} {{\xi + {\xi ^{ - 1}}}}\left[ {f\left( {\xi z} \right)H\left( { - z} \right) + f\left( {{\xi ^{ - 1}}z} \right)H\left( z \right)} \right] \end{equation} \]

where \(\xi \in {\mathbb{R}^ + }\) and \(H(.)\) is the Heaviside function. The absolute moments, required for deriving the central moments, are generated from the following function:

\[ \begin{equation} {M_r} = 2\int_0^\infty {{z^r}f\left( z \right)dz}. \end{equation} \]

The mean and variance are then defined as:

\[\begin{split} \begin{gathered} E\left( z \right) = {M_1}\left( {\xi - {\xi ^{ - 1}}} \right) \hfill \\ Var\left( z \right) = \left( {{M_2} - M_1^2} \right)\left( {{\xi ^2} + {\xi ^{ - 2}}} \right) + 2M_1^2 - {M_2} \hfill \\ \end{gathered} \end{split}\]

The Normal, Student and GED distributions have skew variants which have been standardized to zero mean, unit variance by making use of the moment conditions given above.

Skew Normal distribution

pdqr.dsnorm(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, log: bool = False) → NDArray[float64]

(Skew-Normal) Probability Density Function (dsnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.psnorm(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, lower_tail: bool = True) → NDArray[float64]

(Skew-Normal) Cumulative Probability Function (psnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qsnorm(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, lower_tail: bool = True) → NDArray[float64]

(Skew-Normal) Quantile Function (qsnorm)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rsnorm(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 2, seed: int | None = None) → NDArray[float64]

(Skew-Normal) Random Number Function (rsnorm)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Skew Generalized Error distribution

pdqr.dsged(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, log: bool = False) → NDArray[float64]

(Skew-GED) Probability Density Function (dsged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.psged(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Skew-GED) Cumulative Probability Function (psged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qsged(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Skew-GED) Quantile Function (qsged)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rsged(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 2, shape: float = 4, seed: int | None = None) → NDArray[float64]

(Skew-GED) Random Number Function (rsged)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Skew Student’s T distribution

pdqr.dsstd(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, log: bool = False) → NDArray[float64]

(Skew-Student) Probability Density Function (dsstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.psstd(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Skew-Student) Cumulative Probability Function (psstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qsstd(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Skew-Student) Quantile Function (qsstd)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rsstd(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 2, shape: float = 4, seed: int | None = None) → NDArray[float64]

(Skew-Student) Random Number Function (rsstd)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Johnson’s SU distribution

In the original parameterization of Johnson’s SU distribution, the pdf is given by:

\[ f\left(x, \xi, \lambda, \nu^*, \tau^*\right) = \frac{\tau^*}{\sigma\left(s^2 + 1\right)^{0.5}\sqrt{2\pi}}\exp\left[-\frac{1}{2}z^2\right] \]

for \(-\infty<x<\infty\), \(-\infty<\xi<\infty\), \(\lambda>0\), \(-\infty<\nu^*<\infty\) and \(\tau^*>0\), where

\[ z = \nu^* + \tau^* \sinh^{-1}\left(s\right) = \nu^* + \tau^*\log\left[s + \left(s^2 + 1\right)^{0.5}\right] \]

with \(s = \left(x - \xi)/\lambda\right)\).

Reparameterizing the pdf in terms of the mean and standard deviation we set \(\mu = \xi - \lambda\omega^{0.5}\sinh\left(\nu^*/\tau^*\right)\), \(\sigma = \lambda/c\), \(\nu = -\nu^*\) and \(\tau = \tau^*\). Therefore:

\[ f\left(x, \mu, \sigma, \nu, \tau\right) = \frac{\tau}{c\sigma\left(s^2+1\right)^{0.5}\sqrt(2\pi)}\exp\left[-\frac{1}{2}z^2\right] \]

for \(-\infty<x<\infty\), \(-\infty<\mu<\infty\), \(\sigma>0\), \(-\infty<\nu<\infty\) and \(\tau>0\), where

\[\begin{split} \begin{aligned} z &= -\nu + \tau\sinh^{-1}\left(s\right) = -\nu + \tau\log\left(s + \left(s^2 + 1\right)^{0.5}\right) \hfill \\ s &= \frac{x - \mu + c\sigma\omega^{0.5}\sinh\left(\nu/\tau\right)}{c\sigma} \hfill \\ c &= \left\{0.5\left(\omega - 1\right)\left[\omega\cosh\left(2\nu/\tau\right) + 1\right]\right\}^{-0.5} \hfill \\ \omega &= \exp\left(1/\tau^2\right) \end{aligned} \end{split}\]

The reparameterization is taken from [Rigby et al., 2019] Section 18.4.3 which also contains additional information on the properties of the \(\nu\) and \(\tau\) parameters.

pdqr.djsu(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, log: bool = False) → NDArray[float64]

(Johnson’s SU) Probability Density Function (djsu)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.pjsu(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Johnson’s SU) Cumulative Probability Function (pjsu)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qjsu(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Johnson’s SU) Quantile Function (qjsu)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rjsu(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 2, shape: float = 4, seed: int | None = None) → NDArray[float64]

(Johnson’s SU) Random Number Function (rjsu)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Generalized Hyperbolic distribution

In distributions where the expected moments are functions of all the parameters, it is not immediately obvious how to perform such a transformation. In the case of the Generalized Hyperbolic (GH) distribution, because of the existence of location and scale invariant parameterizations and the possibility of expressing the variance in terms of one of those, namely the \((\zeta, \rho)\), the task of standardizing and estimating the density can be broken down to one of estimating those 2 parameters, representing a combination of shape and skewness, followed by a series of transformation steps to demean, scale and then translate the parameters into the \((\alpha, \beta, \delta, \mu)\) parameterization for which standard formula exist for the likelihood function. The \((\xi, \chi)\) parameterization, which is a simple transformation of the \((\zeta, \rho)\), could also be used in the first step and then transformed into the latter before proceeding further. The only difference is the kind of ‘immediate’ inference one can make from the different parameterizations, each providing a different direct insight into the kind of dynamics produced and their place in the overall GH family particularly with regards to the limit cases.

The package performs estimation using the \((\zeta, \rho)\) parameterization, after which a series of steps transform those parameters into the \((\alpha, \beta, \delta, \mu)\) while at the same time including the necessary recursive substitution of parameters in order to standardize the resulting distribution.

Consider the standardized Generalized Hyperbolic Distribution. Let \(\varepsilon_t\) be a r.v. with mean \((0)\) and variance \(({\sigma}^2)\) distributed as \(\textrm{GH}(\zeta, \rho)\), and let \(z\) be a scaled version of the r.v. \(\varepsilon\) with variance \((1)\) and also distributed as \(\textrm{GH}(\zeta, \rho)\) (the parameters \(\zeta\) and \(\rho\) do not change as a result of being location and scale invariant). The density \(f(.)\) of \(z\) can be expressed as

\[ f(\frac{\varepsilon_t}{\sigma}; \zeta ,\rho ) = \frac{1}{\sigma}f_t(z;\zeta ,\rho ) = \frac{1}{\sigma}f_t(z;\tilde \alpha, \tilde \beta, \tilde \delta ,\tilde \mu ), \]

where we make use of the \((\alpha, \beta, \delta, \mu)\) parameterization since we can only naturally express the density in that parameterization. The steps to transforming from the \((\zeta, \rho)\) to the \((\alpha, \beta, \delta, \mu)\) parameterization, while at the same time standardizing for zero mean and unit variance are given henceforth. Let

\[\begin{split} \begin{eqnarray} \zeta & = & \delta \sqrt {{\alpha ^2} - {\beta ^2}} \hfill \\ \rho & = & \frac{\beta }{\alpha }, \hfill \end{eqnarray} \end{split}\]

which after some substitution may be also written in terms of \(\alpha\) and \(\beta\) as,

\[\begin{split} \begin{eqnarray} \alpha & = & \frac{\zeta }{{\delta \sqrt {(1 - {\rho ^2})} }},\hfill\\ \beta & = &\alpha \rho.\hfill \end{eqnarray} \end{split}\]

For standardization we require that,

\[\begin{split} \begin{eqnarray} E\left(X\right) & = & \mu + \frac{{\beta \delta }}{{\sqrt {{\alpha ^2} - {\beta ^2}} }}\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}} = \mu + \frac{{\beta {\delta ^2}}}{\zeta }\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}} = 0 \hfill \\ \therefore \mu & = & - \frac{{\beta {\delta ^2}}}{\zeta }\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}}\hfill \\ Var\left(X\right) & = & {\delta ^2}\left(\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{\zeta {K_\lambda }\left(\zeta \right)}} + \frac{{{\beta ^2}}}{{{\alpha ^2} - {\beta ^2}}}\left(\frac{{{K_{\lambda + 2}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}} - {\left(\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}}\right)^2}\right)\right) = 1 \hfill\nonumber \\ \therefore \delta & = & {\left(\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{\zeta {K_\lambda }\left(\zeta \right)}} + \frac{{{\beta ^2}}}{{{\alpha ^2} - {\beta ^2}}}\left(\frac{{{K_{\lambda + 2}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}} - {\left(\frac{{{K_{\lambda + 1}}\left(\zeta \right)}}{{{K_\lambda }\left(\zeta \right)}}\right)^2}\right)\right)^{ - 0.5}} \hfill \end{eqnarray} \end{split}\]

Since we can express, \(\beta^2/\left(\alpha^2 - \beta^2\right)\) as,

\[ \frac{{{\beta ^2}}}{{{\alpha ^2} - {\beta ^2}}} = \frac{{{\alpha ^2}{\rho ^2}}}{{{a^2} - {\alpha ^2}{\rho ^2}}} = \frac{{{\alpha ^2}{\rho ^2}}}{{{a^2}\left(1 - {\rho ^2}\right)}} = \frac{{{\rho ^2}}}{{\left(1 - {\rho ^2}\right)}}, \]

then we can re-write the formula for \(\delta\) in terms of the estimated parameters \(\hat\zeta\) and \(\hat\rho\) as,

\[ \delta = {\left(\frac{{{K_{\lambda + 1}}\left(\hat \zeta \right)}}{{\hat \zeta {K_\lambda }\left(\hat \zeta \right)}} + \frac{{{{\hat \rho }^2}}}{{\left(1 - {{\hat \rho }^2}\right)}}\left(\frac{{{K_{\lambda + 2}}\left(\hat \zeta \right)}}{{{K_\lambda }\left(\hat \zeta \right)}} - {\left(\frac{{{K_{\lambda + 1}}\left(\hat \zeta \right)}}{{{K_\lambda }\left(\hat \zeta \right)}}\right)^2}\right)\right)^{ - 0.5}} \]

Transforming into the \((\tilde \alpha ,\tilde \beta ,\tilde \delta ,\tilde \mu )\) parameterization proceeds by first substituting the above value of \(\delta\) into the equation for \(\alpha\) and simplifying:

\[\begin{split} \begin{eqnarray} \tilde \alpha & = & \,{\frac{{\hat \zeta \left( {\frac{{{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{\hat \zeta {{\text{K}}_\lambda }\left( {\hat \zeta } \right)}} + \frac{{{{\hat \rho }^2}\left( {\frac{{{{\text{K}}_{\lambda + 2}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}} - \frac{{{{\left( {{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)} \right)}^2}}}{{{{\left( {{{\text{K}}_\lambda }\left( {\hat \zeta } \right)} \right)}^2}}}} \right)}}{{\left( {1 - {{\hat \rho }^2}} \right)}}} \right)}}{{\sqrt {(1 - {{\hat \rho }^2})} }}^{0.5}}, \hfill\nonumber \\ & = &\,{\frac{{\left( {\frac{{\hat \zeta {{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}} + \frac{{{{\hat \zeta }^2}{{\hat \rho }^2}\left( {\frac{{{{\text{K}}_{\lambda + 2}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}} - \frac{{{{\left( {{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)} \right)}^2}}}{{{{\left( {{{\text{K}}_\lambda }\left( {\hat \zeta } \right)} \right)}^2}}}} \right)}}{{\left( {1 - {{\hat \rho }^2}} \right)}}} \right)}}{{\sqrt {(1 - {\hat \rho ^2})} }}^{0.5}}, \hfill\nonumber \\ & = & {\left( {\left. {\frac{{\frac{{\hat \zeta {{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}}}}{{(1 - {{\hat \rho }^2})}} + \frac{{{\hat \zeta ^2}{\hat \rho ^2}\left( {\frac{{{{\text{K}}_{\lambda + 2}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}\frac{{{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}} - \frac{{{{\left( {{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)} \right)}^2}}}{{{{\left( {{{\text{K}}_\lambda }\left( {\hat \zeta } \right)} \right)}^2}}}} \right)}}{{{{\left( {1 - {{\hat \rho }^2}} \right)}^2}}}} \right)} \right.^{0.5}}, \hfill\nonumber \\ & = & {\left( {\left. {\frac{{\frac{{\hat \zeta {{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}}}}{{(1 - {{\hat \rho }^2})}}\left(1 + \frac{{\hat \zeta {{\hat \rho }^2}\left( {\frac{{{{\text{K}}_{\lambda + 2}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}} - \frac{{{{\text{K}}_{\lambda + 1}}\left( {\hat \zeta } \right)}}{{{{\text{K}}_\lambda }\left( {\hat \zeta } \right)}}} \right)}}{{\left( {1 - {{\hat \rho }^2}} \right)}}\right)} \right)} \right.^{0.5}}. \hfill \end{eqnarray} \end{split}\]

Finally, the rest of the parameters are derived recursively from \(\tilde\alpha\) and the previous results,

\[\begin{split} \begin{eqnarray} \tilde \beta & = & \tilde \alpha \hat \rho,\hfill\\ \tilde \delta & = & \frac{{\hat \zeta }}{{\tilde \alpha \sqrt {1 - {{\hat \rho }^2}} }}, \hfill\\ \tilde \mu & = & \frac{{ - \tilde \beta {{\tilde \delta }^2}{K_{\lambda + 1}}\left(\hat \zeta \right)}}{{\hat \zeta {K_\lambda }\left(\hat \zeta \right)}}.\hfill \end{eqnarray} \end{split}\]

For the use of the \((\xi, \chi)\) parameterization in estimation, the additional preliminary steps of converting to the \((\zeta, \rho)\) are,

\[\begin{split} \begin{eqnarray} \zeta & = & \frac{1}{{{{\hat \xi }^2}}} - 1, \hfill\\ \rho & = & \frac{{\hat \chi }}{{\hat \xi }}. \hfill \end{eqnarray} \end{split}\]

pdqr.dsgh(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 0.5, shape: ArrayLike = 4, lamda: ArrayLike = 0, log: bool = False) → NDArray[float64]

(Standardized Generalized Hyperbolic) Probability Density Function (dsgh)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lamda – the second shape parameter of the GH distribution related to the GIG distribution.

• log – whether to return the log density

Return type:

a numpy array

pdqr.psgh(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 0.5, shape: ArrayLike = 4, lamda: ArrayLike = 0, lower_tail: bool = True) → NDArray[float64]

(Standardized Generalized Hyperbolic) Cumulative Probability Function (psgh)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lamda – the second shape parameter of the GH distribution related to the GIG distribution.

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qsgh(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 0.5, shape: ArrayLike = 4, lamda: ArrayLike = 0, lower_tail: bool = True) → NDArray[float64]

(Standardized Generalized Hyperbolic) Quantile Function (qsgh)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lamda – the second shape parameter of the GH distribution related to the GIG distribution.

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rsgh(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 0.5, shape: float = 4, lamda: float = 0, seed: int | None = None) → NDArray[float64]

(Standardized Generalized Hyperbolic) Random Number Function (rsgh)

Parameters:

• n – the number of draws

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lamda – the second shape parameter of the GH distribution related to the GIG distribution.

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array

Generalized Hyperbolic Skew Student’s T distribution

The GH Skew-Student distribution was popularized by [Aas and Haff, 2006] because of its uniqueness in the GH family in having one tail with polynomial and one with exponential behavior. This distribution is a limiting case of the GH when \(\alpha \to \left| \beta \right|\) and \(\lambda=-\nu/2\), where \(\nu\) is the shape parameter of the Student distribution. The domain of variation of the parameters is \(\beta \in \mathbb{R}\) and \(\nu>0\), but for the variance to be finite \(\nu>4\), while for the existence of skewness and kurtosis, \(\nu>6\) and \(\nu>8\) respectively. The density of the random variable \(x\) is then given by:

\[ f\left( x \right) = \frac{{{2^{\left( {1 - \nu } \right)/2}}{\delta ^\nu }{{\left| \beta \right|}^{\left( {\nu + 1} \right)/2}}{K_{\left( {\nu + 1} \right)/2}}\left( {\sqrt {{\beta ^2}\left( {{\delta ^2} + {{\left( {x - \mu } \right)}^2}} \right)} } \right)\exp \left( {\beta \left( {x - \mu } \right)} \right)}} {{\Gamma \left( {\nu /2} \right)\sqrt \pi {{\left( {\sqrt {{\delta ^2} + {{\left( {x - \mu } \right)}^2}} } \right)}^{\left( {\nu + 1} \right)/2}}}} \]

To standardize the distribution to have zero mean and unit variance, I make use of the first two moment conditions for the distribution which are:

\[\begin{split} \begin{gathered} E\left( x \right) = \mu + \frac{{\beta {\delta ^2}}} {{\nu - 2}} \hfill \\ Var\left( x \right) = \frac{{2{\beta ^2}{\delta ^4}}} {{{{\left( {\nu - 2} \right)}^2}\left( {\nu - 4} \right)}} + \frac{{{\delta ^2}}} {{\nu - 2}} \hfill \\ \end{gathered} \end{split}\]

We require that \(Var(x)=1\), thus:

\[ \delta = {\left( {\frac{{2{{\bar \beta }^2}}}{{{{\left( {\nu - 2} \right)}^2}\left( {\nu - 4} \right)}} + \frac{1}{{\nu - 2}}} \right)^{ - 1/2}} \]

where I have made use of the \(4^{th}\) parameterization of the GH distribution given in [Prause, 1999] where \(\hat \beta = \beta \delta\). The location parameter is then rescaled by substituting into the first moment formula \(\delta\) so that it has zero mean:

\[ \bar \mu = - \frac{{\beta {\delta ^2}}}{{\nu - 2}} \]

Therefore, we model the GH Skew-Student using the location-scale invariant parameterization \((\bar \beta, \nu)\) and then translate the parameters into the usual GH distribution’s \((\alpha, \beta, \delta, \mu)\), setting \(\alpha = abs(\beta)+1e-12\).

pdqr.dsghst(x: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, log: bool = False) → NDArray[float64]

(Standardized Generalized Hyperbolic Skew Student) Probability Density Function (dsghst)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to x (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• x – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• log – whether to return the log density

Return type:

a numpy array

pdqr.psghst(q: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Standardized Generalized Hyperbolic Skew Student) Cumulative Probability Function (psghst)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to q (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• q – a vector of quantiles

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.qsghst(p: ArrayLike, mu: ArrayLike = 0, sigma: ArrayLike = 1, skew: ArrayLike = 2, shape: ArrayLike = 4, lower_tail: bool = True) → NDArray[float64]

(Standardized Generalized Hyperbolic Skew Student) Quantile Function (qsghst)

In cases when the parameters are vectors of size n (>1), then they must be equal in size to p (unless they are of size 1 in which case they are recycled). This is achieved using the numpy vectorize function.

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• lower_tail – if TRUE (default), probabilities are P[X ≤ x] otherwise, P[X > x]

Return type:

a numpy array

pdqr.rsghst(n: int = 1, mu: float = 0, sigma: float = 1, skew: float = 2, shape: float = 4, seed: int | None = None) → NDArray[float64]

(Standardized Generalized Hyperbolic Skew Student) Random Number Function (rsghst)

Parameters:

• p – a vector of probabilities

• mu – the mean

• sigma – the standard deviation

• skew – the skew parameter

• shape – the shape parameter

• seed – an optional value to initialize the random seed generator

Return type:

a numpy array