scipy.stats.johnsonsu中的a和b参数是什么?

What are the a and b parameters in scipy.stats.johnsonsu?

我正在尝试将 Johnson SU 分布拟合到标准普尔 500 指数 return 的经验集。我的理解(免责声明:不是数学家)是这个分布包含第三和第四时刻(偏斜和峰度)。除了 loc(均值)和 scale(标准偏差)之外,johnsonsu 还采用两个附加参数,ab。但是这些参数的顺序和规格比较混乱。

这就是我的困惑所在:如果我将 return 引入 SPDR S&P 500 ETF Trust (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 应该 return locscale 参数。 (正态分布所需的一切。)检查出来:

import scipy.stats as scs

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

True
True

但不太清楚return由scs.johnsonsu.fit编辑的是什么:

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

这些应该是 four parameters 分布:xi, gamma, delta, lam.

但我无法让他们回到经验平均值,应该是:

即:

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')

得到一堆 NaNs。

它们是Johnson SU的参数。请记住,您所获得的样本均值与分布均值不同。这是平均值的表达式

这里是方差表达式:

在您的代码中,ξ 为 loc,λ 为 scale,γ 为 a,δ 为 b。 sinh-1(x) 等于 log(x + sqrt(1 + x2)).

因此检查拟合的return,为所有四个参数赋值,然后计算分布均值并与样本均值进行比较。如果有效,重复方差练习

更新

我试过你的代码并建议检查均值和方差,效果很好,请在下面检查

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