使用 statsmodels 计算 Gamma GLM 的尺度/分散答案

【问题标题】：Calculating scale/dispersion of Gamma GLM using statsmodels使用 statsmodels 计算 Gamma GLM 的尺度/分散
【发布时间】：2020-05-29 14:09:01
【问题描述】：

我在使用 statsmodels 的 GLM 函数获取模拟数据的分散参数时遇到问题。

import statsmodels.api as sm
import matplotlib.pyplot as plt 
import scipy.stats as stats 
import numpy as np

np.random.seed(1)

# Generate data
x=np.random.uniform(0, 100,50000)
x2 = sm.add_constant(x)
a = 0.5
b = 0.2
y_true = 1/(a+(b*x))
# Add error 
scale = 2 # the scale parameter I'm trying to obtain
shape = y_true/scale # given that, for Gamma, mu = scale*shape
y = np.random.gamma(shape=shape, scale=scale)

# Run model
model = sm.GLM(y, x2, family=sm.families.Gamma()).fit() 

model.summary()

以下是上面的摘要：

请注意，系数估计值是正确的（0.5 和 0.2），但规模 (21.995) 与我设置的规模 (2) 相差甚远。

有人能指出我误解/做错了什么吗？谢谢！

【问题讨论】：

在这种情况下，GLM 中的 scale 只是 pearson_chi2 / df_resid。我不知道这与通常的 Gamma 参数化有何关系。 statmodels GLM 不通过最大似然估计尺度。
据我所知，GLM 参数化对应于y = np.random.gamma(shape=1 / scale, scale=y_true * scale)。
另外，如果将 x 的上限减小到 10，则结果看起来会更好，因为它避免了平均值的小值。
我只是在模式匹配 GLM 和 GAMMA loglike 函数。 Statsmodels 并没有真正使用这种等效性，仅使用 GLM 均值离散参数化，因此记录和支持重新参数化仍然不完整。
（快速搜索）pj.freefaculty.org/guides/stat/Regression-GLM/Gamma/… 中的第 2 部分证实了这一点。在单参数族 GLM 版本中，我们将 Gamma 形状参数作为固定参数，而 Gamma 尺度参数与均值相关，即期望值。这对应于 statsmodels 中的实现。

标签： python statistics statsmodels

【解决方案1】：

正如 Josef 在 cmets 中指出的那样，statsmodels 使用了不同类型的参数化。

【讨论】：