<jupyter_start><jupyter_text>Introduction  P.H.A.T. - **P**areto **H**ybrids with **A**symmetric **T**ails The **Phat** distribution is an attempt to address the issues of fat tails in two-tailed data. It is a two-tailed, fully-continuous, well-defined asymmetric power law probability distribution. It is a mixture model of two Pareto hybrid distributions, as described in [2009 by <NAME> and <NAME>](https://www.researchgate.net/publication/226293435_A_hybrid_Pareto_model_for_asymmetric_fat-tailed_data_The_univariate_case) with:+ Gaussian body + distinct Pareto power laws in either tail.The distribution requires only 4 parameters:+ $\mu, \sigma$ in the Gaussian body+ $\xi_{\text{left}}, \xi_{\text{right}}$, being the inverse tail index (1/$\alpha$) for either Paretian tail.The **phat-tails** package makes available several methods to fit a given time-series dataset to the parameters of the Phat distribution and produce a forecast with the results.  Installation Installation available via `pip````console$ pip install phat-tails```  Quickstart All pertinent classes and functions are imported via the module `phat`.<jupyter_code>%load_ext autoreload
%autoreload 2
import phat as ph<jupyter_output><empty_output><jupyter_text>The probability distribution is found in the main class, `Phat`, which mimics the structure of the [continuous distributions found in scipy.stats](https://docs.scipy.org/doc/scipy/tutorial/stats/continuous.htmlcontinuous-distributions-in-scipy-stats). We pass the four parameters to instantiate the distribution. For simplicity, we will show the distribution with equal tail indices initially.<jupyter_code>mean, sig, shape = 0, 1, 1/5
phat_dist = ph.Phat(mean, sig, shape, shape)<jupyter_output><empty_output><jupyter_text>Below is a complete rendering of the distribution, with breakdown among the component Gaussian and Pareto tails.<jupyter_code>import numpy as np
import scipy.stats as scist
import matplotlib.pyplot as plt
import seaborn as sns; sns.set(style = 'whitegrid')

x = np.linspace(-5+mean, 5+mean, 1000)

fig, ax1 = plt.subplots(1,1,figsize=(10,6))

ax1.plot(x, phat_dist.pdf(x), alpha=.5)
x_body = np.linspace(phat_dist.left.a, phat_dist.right.a, 100)
x_left = np.linspace(x[0], phat_dist.left.a, 100)
x_right = np.linspace(phat_dist.right.a, x[-1], 100)
ax1.plot(x_body, phat_dist.pdf(x_body), c='C1', ls='--', label='Gaussian Body')
ax1.plot(x_left, phat_dist.pdf(x_left), c='C2', ls='--', label='Left Paretian')
ax1.plot(x_right, phat_dist.pdf(x_right), c='C4', ls='--', label='Right Paretian')

ax1.axvline(
    phat_dist.left.a, .825, .925, 
    c='r', lw=1, label=r'Left Junction; $x = a_l$')
ax1.axvline(phat_dist.right.a, .825, .925, c='r', lw=1, label=r'Right Junction; $x = a_r$')

paramtxt = 'Body'
paramtxt += '\n'
paramtxt += r'$\mu, \sigma$ = ' + f'{phat_dist.mu:.0f}, {phat_dist.sig:.0f}'

ax1.text(
    .5, .5, paramtxt, ha='center',
    transform=ax1.transAxes,
    bbox=dict(boxstyle='round', ec='.8', fc='w')
)

paramtxt = r'Left Tail$_{}$'
paramtxt += '\n'
paramtxt += r'$\epsilon_{l}$ = ' + f'{phat_dist.left.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$a_l$ = ' + f'{phat_dist.left.a:.2f}'
paramtxt += '\n'
paramtxt += r'$b_l$ = ' + f'{1 / phat_dist.left.b:.2f}'
paramtxt += '\n'
paramtxt += r'$\alpha_l$ = ' + f'{1 / phat_dist.left.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$\gamma_l = $' + f'{phat_dist.right.gamma:.2f}'

ax1.text(
    .02,.3, paramtxt,
    transform=ax1.transAxes,
    bbox=dict(boxstyle='round', ec='.8', fc='w')
)
paramtxt = r'Right Tail$_{}$'
paramtxt += '\n'
paramtxt += r'$\epsilon_r$ = ' + f'{phat_dist.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$a_r$ = ' + f'{phat_dist.right.a:.2f}'
paramtxt += '\n'
paramtxt += r'$b_r$ = ' + f'{1 / phat_dist.right.b:.2f}'
paramtxt += '\n'
paramtxt += r'$\alpha_r$ = ' + f'{1 / phat_dist.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$\gamma_r = $' + f'{phat_dist.right.gamma:.2f}'

ax1.text(
    .98,.3, paramtxt, ha='right',
    transform=ax1.transAxes,
    bbox=dict(boxstyle='round', ec='.8', fc='w')
)
ax1.set_xlabel('x')
ax1.set_ylabel('f(x)', loc='top', rotation='horizontal')

ax1.legend()
ax1.set_title('PDF - Phat')

plt.show()<jupyter_output><empty_output><jupyter_text>Below we demonstrate the ability to generate asymmetric tails. We overlay two different Phat distributions, one with symmetric tail indices of $\alpha=2$ and the other with *asymmetric* tail indices, $\alpha_{\text{left}}=2$ and $\alpha_{\text{right}}=20$.We can see that the left tails are identical. In the right tails, the distributions appear to differ only modestly, however, [this difference leads to dramatically different effects](demo.ipynbCompare-Fit-with-Gaussian-and-T).<jupyter_code>mean, sig = 0, 1
shape_l1, shape_r = 1/2, 1/2
dist1 = ph.Phat(mean, sig, shape_l1, shape_r)
shape_l2, shape_r = 1/2, 1/20
dist2 = ph.Phat(mean, sig, shape_l2, shape_r,)
mean, sig = 0, 1
x = np.linspace(-4+mean, 7+mean, 1000)
shape_l1, shape_r = 1/2, 1/2
dist1 = ph.Phat(mean, sig, shape_l1, shape_r)
shape_l2, shape_r = 1/2, 1/20
dist2 = ph.Phat(mean, sig, shape_l2, shape_r,)

fig, ax1 = plt.subplots(1,1,figsize=(10,6))

ax1.plot(x, dist1.pdf(x), label='Fatter Right')
ax1.plot(x, dist2.pdf(x), label='Thinner Right')

paramtxt = r'Fatter Right$_{}$'
paramtxt += '\n'
paramtxt += r'$\epsilon_{l}$ = ' + f'{dist1.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$a_l$ = ' + f'{dist1.right.a:.2f}'
paramtxt += '\n'
paramtxt += r'$b_l$ = ' + f'{1 / dist1.right.b:.2f}'
paramtxt += '\n'
paramtxt += r'$\alpha_l$ = ' + f'{1 / dist1.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$\gamma_l = $' + f'{dist1.right.gamma:.2f}'

ax1.text(
    .02,.95, paramtxt, va='top',
    transform=ax1.transAxes,
    bbox=dict(boxstyle='round', ec='.8', fc='w')
)
paramtxt = r'Thinner Right$_{}$'
paramtxt += '\n'
paramtxt += r'$\epsilon_{l}$ = ' + f'{dist2.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$a_l$ = ' + f'{dist2.right.a:.2f}'
paramtxt += '\n'
paramtxt += r'$b_l$ = ' + f'{1 / dist2.right.b:.2f}'
paramtxt += '\n'
paramtxt += r'$\alpha_l$ = ' + f'{1 / dist2.right.xi:.2f}'
paramtxt += '\n'
paramtxt += r'$\gamma_l = $' + f'{dist2.right.gamma:.2f}'

ax1.text(
    .02,.2, paramtxt, va='bottom',
    transform=ax1.transAxes,
    bbox=dict(boxstyle='round', ec='.8', fc='w')
)
ax1.set_xlabel('x')
ax1.set_ylabel('f(x)', loc='top', rotation='horizontal')

ax1.legend()

ax1.set_title('PDF')

plt.show()<jupyter_output><empty_output><jupyter_text>The `Phat` class has common methods such as `pdf`, `cdf`, `sf`, `ppf`. It can also calculate negative log-likelihood and first and second moments. Derivations are [found here](moments.ipynb).<jupyter_code>mean, sig, shape_l, shape_r = 0,1, 1/5, 1/4
phat_dist = ph.Phat(mean, sig, shape_l, shape_r)
phat_dist.pdf(10)
phat_dist.cdf([.05,1,-0.1])
phat_dist.sf([.05])
assert phat_dist.sf([.05]) == 1 - phat_dist.cdf([.05])
import numpy as np
phat_dist.ppf(np.linspace(0,1,5))
phat_dist.nll(1) # Negative Log-Likelihood
phat_dist.mean()
phat_dist.std()<jupyter_output><empty_output><jupyter_text>It can also generate random variables (and standardized random variables).<jupyter_code>phat_dist.rvs(20)
phat_dist.std_rvs(20)<jupyter_output><empty_output><jupyter_text>Importantly, Phat captures the undefined moments that result when $\alpha < 2$.<jupyter_code>shape_l, shape_r, mean, sig = 1, 1, 0, 1
phat_dist = ph.Phat(mean, sig, shape_l, shape_r)
phat_dist.mean()
phat_dist.var()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18375838926174495, 'fraction_numerical': 0.031140939597315436, 'mean_word_length': 5.087418300653595, 'pattern_counts': {'":': 0, '<': 32, '<?xml version=': 0, '>': 31, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '252004', 'license': 'MIT', 'n_tokens_mistral': 3107, 'n_tokens_neox': 2951, 'n_words': 669, 'path': 'docs/notebooks/intro.ipynb', 'repo_name': 'rskene/phat', 'size': '7464'}