CODEWITHSUNDEEP
pythonscipyscipy.stats
Luna Paluna
Luna Paluna

Reputation: 73

Different results when using the st.gumbel_r.ppf()

I want to calculate the value for a given probability using the st.gumbel_r.ppf(). I'm comparing it with the analytical solution and it is giving me completely different results, anyone knows why? I obtained the values for the scale and location using the moment of methods

P_class = [0.14285714, 0.28571429, 0.42857143, 0.57142857, 0.71428571, 0.85714286, 0.999999]

u = 8.590342451210152
alpha = 0.1841827435642898

h_class1 = st.gumbel_r.ppf(P_class, scale = u, loc = alpha)
h_class = u-np.log(-np.log(P_class))/alpha

Results
h_class1 = [ -5.53466431,  -1.7516637 ,   1.6076281 ,   5.17091797, 9.54112426,  16.24661736, 118.86414528]
h_class = [4.97583548,  7.36682128,  9.49000861, 11.74212971, 14.50424958, 18.74235061, 83.60013865]

I want to get the same h_class results when using the scipy function

Upvotes: 1

Views: 212

Answers (2)

piterbarg
piterbarg

Reputation: 8219

You got your scale and location parameters mixed up, the correct expression for h_class is below.

import scipy.stats as st
P_class = [0.14285714, 0.28571429, 0.42857143, 0.57142857, 0.71428571, 0.85714286, 0.999999]

u = 8.590342451210152
alpha = 0.1841827435642898

h_class1 = st.gumbel_r.ppf(P_class, scale = u, loc = alpha)
h_class = alpha-np.log(-np.log(P_class))*u
print(h_class1)
print(h_class)

and you get a perfect match

[ -5.5346644   -1.7516636    1.60762813   5.17091794   9.54112411
  16.24661755 118.86414528]
[ -5.5346644   -1.7516636    1.60762813   5.17091794   9.54112411
  16.24661755 118.86414528]

Upvotes: 1

Max Pierini
Max Pierini

Reputation: 2269

Given your scale-location values, the right-skewed Gumpel ditribution you defined is

dist = sps.gumbel_r(loc=alpha, scale=u)
x = np.linspace(dist.ppf(.0001), dist.ppf(.99999), 1000)
pdf = dist.pdf(x)
cdf = dist.cdf(x)
fig, ax = plt.subplots(1, 2, figsize=(12, 4))
ax[0].plot(x, pdf)
ax[0].set_title('PDF')
ax[1].plot(x, cdf)
ax[1].set_title('CDF')
for p in P_class:
    ax[1].axhline(p, color='r', ls='--', lw=1)
    ax[1].axvline(dist.ppf(p), color='k', ls=':', lw=1)
    print('CDF:', p, 'X:', dist.ppf(p))
plt.show()

with the ppf (Percent Point Function) you're evaluating the value of x at the value p of the Cumulative Density Function (CDF)

enter image description here

CDF: 0.14285714 X: -5.534664397573627
CDF: 0.28571429 X: -1.7516635973033354
CDF: 0.42857143 X: 1.6076281292815586
CDF: 0.57142857 X: 5.170917935306729
CDF: 0.71428571 X: 9.541124106213145
CDF: 0.85714286 X: 16.246617550548397
CDF: 0.999999 X: 118.86414527936681

So maybe, either the scale-location you got is wrong or what you're searching for is different.

Upvotes: 0

Related Questions