In [2]:

    

import numpy as np
import pandas as pd
from scipy import stats

import pyreto
from pyreto.tests import utilities


    

In [3]:

    

p = lambda alpha, xmin: stats.pareto(alpha-1, loc=0, scale=xmin)


    

In [4]:

    

prng = np.random.RandomState(42)
alpha, xmin = 3.0, 1.0
pareto_data = pd.Series(pyreto.distributions.Pareto._rvs(prng, 2500, xmin, alpha), name='samples')


    

In [5]:

    

# check that I get same estimate for alpha given reported xmin...
result1 = pyreto.distributions.Pareto.fit(pareto_data, xmin=xmin)


    

In [6]:

    

result1.params['alpha']


    

    
Out[6]:





2.9998119111215154

In [7]:

    

result1.xmin


    

    
Out[7]:





1.0

In [8]:

    

f = p(result1.params['alpha'], 1.0)


    

In [9]:

    

f.pdf(5.0)


    

    
Out[9]:





0.016003339044528487

In [10]:

    

pyreto.distributions.Pareto._pdf(5.0, 1.0, result1.params['alpha'])


    

    
Out[10]:





0.016003339044528487

In [12]:

    

f = stats.pareto(2.0, 1.0)


    

In [13]:

    

stats.pareto.


    

In [28]:

    

pyreto.distributions.Pareto._pdf??


    

In [22]:

    

# check that I get same estimate for alpha given reported xmin...
result2 = pyreto.distributions.Pareto.fit(pareto_data, xmin=None, quantile=0.9, method='brute')


    

In [23]:

    

result2.params['alpha']


    

    
Out[23]:





1.9909069323487545

In [24]:

    

result2.xmin


    

    
Out[24]:





1.1331781782153696

In [25]:

    

# check that I get same estimate for alpha given reported xmin...
result3 = pyreto.distributions.Pareto.fit(pareto_data, xmin=None, quantile=0.9, method='bounded')


    

In [26]:

    

result3.params['alpha']


    

    
Out[26]:





1.9924312888829019

In [27]:

    

result3.xmin


    

    
Out[27]:





1.3693896782872266

In [7]:

    

prng = np.random.RandomState(1234)
cauchy_rvs = stats.cauchy.rvs(size=2500, random_state=prng)
cauchy_data = pd.Series(cauchy_rvs, name='samples')


    

In [8]:

    

# check that I get same estimate for alpha given reported xmin...
result = pyreto.distributions.Pareto.fit(cauchy_data, xmin=None, quantile=0.99, method='brute')


    

In [12]:

    

# will raise assertion error if test fails!
desired_alpha = 2.0
utilities.test_scaling_exponent_estimation(desired_alpha, result)


    

In [2]:

    

prng = np.random.RandomState(4321)
desired_alpha, beta = 1.5, 1.0 
stable_rvs = stats.levy_stable.rvs(desired_alpha - 1, beta, size=1000, random_state=prng)
stable_data = pd.Series(stable_rvs, name='samples')


    

In [3]:

    

result = pyreto.distributions.Pareto.fit(stable_data, xmin=None, quantile=0.99, method='bounded')


    

In [4]:

    

# will raise assertion error if test fails!
utilities.test_scaling_exponent_estimation(desired_alpha, result)


    

In [84]:

    

# goodness of fit test should fail to reject Pareto tails!
test_kwargs = {'seed': None, 'result': result, 'data': stable_data,
               'xmin': None, 'quantile': 0.99, 'method': 'bounded'}
pvalue, _ = pyreto.distributions.Pareto.test_goodness_of_fit(**test_kwargs)
assert pvalue > 0.1, "Goodness of fit test Type I error!"


    

In [ ]: