Moving to GB2 Distribution

[davefoote]

I'm sorry if this is a stupid question but when moving to the GB2 Distribution do we use our 'm' estimate directly as 'p' or is there some conversion that needs to take place?

[rickecon]

@davefoote . Good question. The progression is the following, which is also explained in Section 7 of the maximum likelihood notebook:

1. Estimate the GA(x|alpha,beta) using initial guesses of beta_init=Var(x)/E(x) and alpha_init=E(x)/beta_init.
2. Because the GA(x|alpha,beta) distribution is a nested case of the GG(x|alpha,beta,m) distribution such that GA(x|alpha,beta)=GG(x|alpha,beta,m=1), estimate the GG(x|alpha,beta,m) distribution using initial guesses alpha_hat and beta_hat from step (1) estimates of the GA distribution and initial guess m_init=1.
3. Because the GG(x|alpha,beta,m) distribution is a nested case of the GB2(x|a,b,p,q) distribution such that GG(x|alpha,beta,m)=GB2(x|a=m,b=q**(1/m),p=alpha/m,q), estimate the GB2(x|a,b,p,q) distribution using initial guesses a_init=m_hat, p_init=alpha_hat/m_hat from step (2) and b_init and q_init using large q like q=1000. This approximates infinity.

[rickecon]

@davefoote . Here is some good GB2 code.

import numpy as np
import scipy.special as spc

def GB2_pdf(xvals, aa, bb, pp, qq):
    '''
    --------------------------------------------------------------------
    Returns the PDF values from the four-parameter generalized beta 2
    (GB2) distribution. See McDonald and Xu (1995).

    (GB2): f(x; a, b, p, q) = (a * (x ** ((a*p) - 1))) /
        ((b ** (a * p)) * spc.beta(p, q) *
        ((1 + ((x / b) ** a)) ** (p + q)))
    x in [0, infty), alpha, beta, m > 0
    --------------------------------------------------------------------
    INPUTS:
    xvals = (N,) vector, values in the support of generalized beta 2
            (GB2) distribution
    aa    = scalar > 0, generalized beta 2 (GB2) distribution parameter
    bb    = scalar > 0, generalized beta 2 (GB2) distribution parameter
    pp    = scalar > 0, generalized beta 2 (GB2) distribution parameter
    qq    = scalar > 0, generalized beta 2 (GB2) distribution parameter

    OTHER FUNCTIONS AND FILES CALLED BY THIS FUNCTION:
        spc.beta()

    OBJECTS CREATED WITHIN FUNCTION:
    pdf_vals = (N,) vector, pdf values from generalized beta 2 (GB2)
               distribution corresponding to xvals given parameters aa,
               bb, pp, and qq

    FILES CREATED BY THIS FUNCTION: None

    RETURNS: pdf_vals
    --------------------------------------------------------------------
    '''
    pdf_vals = \
        np.float64((aa * (xvals ** (aa * pp - 1))) / ((bb ** (aa * pp)) *
                   spc.beta(pp, qq) *
                   ((1 + ((xvals / bb) ** aa)) ** (pp + qq))))

    return pdf_vals