scipy.stats.johnsonsuのaとbのパラメータは何ですか？

Question

私はJohnson SUの分布を経験豊富なS & P 500に返します。私の理解（免責事項：数学者ではない）は、この分布は3番目と4番目の瞬間（スキュー&尖度）を取り入れているということです。 loc（平均）およびscale（標準偏差）に加えて、johnsonsuは、2つの追加パラメータ、aおよびbを取ります。しかし、これらのパラメータの順序と仕様は混乱します。私の混乱が由来のはここ

です：私はSPDR S & P 500 ETF信託（SPY）に復帰して引っ張った場合、私は以下の経験的な統計情報を取得する：

今

from pandas_datareader.data import DataReader as dr 
r = dr('SPY', 'google', start='2000')['Close'].pct_change().dropna() 
mean, var, std, skew, kurt = r.mean(), r.var(0), r.std(0), r.skew(), r.kurt() # ddof = 0 
# mean: 0.00027732907268771364 
# var: 0.00014416720067485022 
# std: 0.012006964673673785 

私は正規分布に合わせた場合この経験的データ、.fitは、locとscaleのパラメータを返すことになっています。 （正規分布のために必要とされているすべてのこと。）チェックアウトすること：

import scipy.stats as scs 

normmean, normstd = scs.norm.fit(r) 
print(np.allclose(normmean, mean)) 
print(np.allclose(normstd, std)) 

True 
True 

しかし、scs.johnsonsu.fitによって返されるものそれほど明確ではない：

print(scs.johnsonsu.fit(r)) 
(0.098009661042083682, 1.022060362199493, 0.0013471690867203458, 0.0072653444313926403) 

これらの分布にfour parametersする必要があります：xi、ガンマ、デルタ、ラム。

しかし、私は彼らがあるべき経験的平均、にバックアップすることができません。

すなわち：

def johnsonmean(gamma, xi, delta, lam): 
    mean = xi - lam * np.exp(delta ** -2/2) * np.sinh(gamma/delta) 
    return mean 
gamma, xi, delta, lam = scs.johnsonsu.fit(r) # correct order? 
print(johnsonmean(gamma, xi, delta, lam)) 
-inf 

そして

mean, var, skew, kurt = scs.johnsonsu.stats(loc=xi, scale=lam, 
              a=gamma, b=delta, moments='msvk') 

はNaN Sの束を取得します。

Answer 1

これらはJohnson SUのパラメータです。サンプルの平均値が分布の平均と同じではないことを覚えておいてください。ここでは平均値

そして、ここのための式は分散の式である：あなたのコードで

、ξはloc、λだろうscaleなり、γはaだろうδはbとなります。 sinh ^-1（x）はlog（x + sqrt（1 + x ））に等しい。

したがって、フィットの戻りをチェックし、4つのパラメータすべてに値を割り当て、次に分布平均を計算してサンプル平均と比較します。それが動作する場合、分散

UPDATE

私はあなたのコードを試してみましたが、平均と分散のための小切手を提案し、それが

import sys 
import math 

from pandas_datareader.data import DataReader as dr 
import scipy.stats as scs 

def read_data(): 
    return dr('SPY', 'google', start='2000')['Close'].pct_change().dropna() 

def johnsonsu_mean(a, b, loc, scale): 
    """ 
    Johnson SU mean according to https://en.wikipedia.org/wiki/Johnson%27s_SU-distribution 
    """ 
    v = loc - scale * math.exp(0.5/b**2) * math.sinh(a/b) 
    return v 

def johnsonsu_var(a, b, loc, scale): 
    """ 
    Johnson SU variance according to https://en.wikipedia.org/wiki/Johnson%27s_SU-distribution 
    """ 
    t = math.exp(1.0/b**2) 
    v = 0.5*scale**2 * (t - 1.0) * (t * math.cosh(2.0*a/b) + 1.0) 
    return v 

def johnsonsu_median(a, b, loc, scale): 
    """ 
    Johnson SU median according to https://en.wikipedia.org/wiki/Johnson%27s_SU-distribution 
    """ 
    v = loc + scale * math.sinh(-a/b) 
    return v 

def main(r): 
    sample_mean, sample_med, sample_var, sample_std, sample_skew, sample_kurt = r.mean(), r.median(), r.var(0), r.std(0), r.skew(), r.kurt() 

    a, b, loc, scale = scs.johnsonsu.fit(r) # fit the data and get distribution parameters back 

    # distribution mean and variance according to SciPy 
    dist_mean = scs.johnsonsu.mean(a, b, loc, scale) 
    dist_med = scs.johnsonsu.median(a, b, loc, scale) 
    dist_var = scs.johnsonsu.var(a, b, loc, scale) 

    # distribution mean, var vs sample ones 
    print("{0} {1}".format(sample_mean, dist_mean)) 
    print("{0} {1}".format(sample_med, dist_med)) 
    print("{0} {1}".format(sample_var, dist_var)) 
    print("") 

    # distribution mean and variance according to Wiki vs SciPy 
    print("{0} {1}".format(dist_mean, johnsonsu_mean(a, b, loc, scale))) 
    print("{0} {1}".format(dist_var, johnsonsu_var(a, b, loc, scale))) 
    print("{0} {1}".format(dist_med, johnsonsu_median(a, b, loc, scale))) 

if __name__ == "__main__": 
    r = read_data() 
    main(r) 

    sys.exit(0) 

下に確認してください、うまく機能するために練習を繰り返し、出力生成：

0.00028012130615107805 0.00021391570000183283 
0.0005697194131890626 0.0006458197694718355  
0.00014415554662672425 0.00015479059187195545 

0.00021391570000183283 0.00021391570000541633 
0.00015479059187195545 0.00015479059186527505 
0.0006458197694718355 0.0006458197694718355