statsmodels-绘制拟合分布

Question

以下代码适合使用以下简化的广义线性模型： statsmodels

model = smf.glm('Y ~ 1', family=sm.families.NegativeBinomial(), data=df)
results = model.fit()

这给出了系数和标准差：

               coef stderr   
Intercept    2.9471  0.120

现在，我想以图形方式将变量Y（直方图）的实际分布与来自模型的分布进行比较。

但是我需要两个参数r并p对其进行评估stats.nbinom(r,p)和绘制。

有没有办法从拟合结果中检索参数？如何绘制PMF？

Answer 1

统计模型中的广义线性模型GLM当前不估计负二项分布的额外参数。负二项式仅对于固定形状参数属于指数分布族。

但是，statsmodels还具有负二项式作为最大似然模型，discrete_model其中估计了所有参数。

用于计数回归的负二项式的参数化是根据平均值或期望值，这与scipy.stats.nbinom中的参数化不同。实际上，负二项式计数回归有两种不同的常用参数设置，通常称为nb1和nb2

下面是恢复scipy.stats.nbinom参数，一个迅速写好的剧本n=size，并p=prob从估计的参数。一旦有了scipy.stats.distribution的参数，就可以使用所有可用的方法rvs，pmf等。

这样的东西应该在statsmodels中可用。

在一些示例运行中，我得到了这样的结果

data generating parameters 50 0.25
estimated params           51.7167511571 0.256814610633
estimated params           50.0985814878 0.249989725917

另外，由于潜在的指数重新参数化，scipy优化器有时会收敛。在这些情况下，提供更好的起始值或使用Nelder-Mead作为优化方法通常会有所帮助。

import numpy as np
from scipy import stats
import statsmodels.api as sm

# generate some data to check
nobs = 1000
n, p = 50, 0.25
dist0 = stats.nbinom(n, p)
y = dist0.rvs(size=nobs)
x = np.ones(nobs)

loglike_method = 'nb1'  # or use 'nb2'
res = sm.NegativeBinomial(y, x, loglike_method=loglike_method).fit(start_params=[0.1, 0.1])

print dist0.mean()
print res.params

mu = res.predict()   # use this for mean if not constant
mu = np.exp(res.params[0])   # shortcut, we just regress on a constant
alpha = res.params[1]

if loglike_method == 'nb1':
    Q = 1
elif loglike_method == 'nb2':    
    Q = 0

size = 1. / alpha * mu**Q
prob = size / (size + mu)

print 'data generating parameters', n, p
print 'estimated params          ', size, prob

#estimated distribution
dist_est = stats.nbinom(size, prob)

顺便说一句：我之前遇到过这个问题，但是没有时间看一下 https://github.com/statsmodels/statsmodels/issues/106