用于GLM伽玛分布的SAS和R之间的标准误差

Question

我将模型的R和SAS的GLM输出与Gamma分布进行比较.点估计是相同的,但它们具有不同的标准误差估计,因此具有不同的p值.

有谁知道为什么？我想知道R和SAS是否使用不同的方法来估计标准误差？也许MLE与时刻的方法？

R示例代码

set.seed(2)
test = data.table(y = rnorm(100, 1000, 100), x1 = rnorm(100, 50, 20), x2 = rgamma(100, 0.01))
model = summary(glm(formula = y ~ x1+x2 , family = Gamma(link = "log"), data = test))

使用此处生成的相同数据,我使用以下代码在SAS中运行模型:

proc genmod data= test_data;
                model y =  x1 x2 /link= log dist= gamma;
    run;

R的输出如下:

Call:
glm(formula = y ~ x1 + x2, family = Gamma(link = "log"), data = test)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-0.26213  -0.08456  -0.01033   0.08364   0.20878  

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  6.9210757  0.0324674 213.170   <2e-16 ***
x1          -0.0003371  0.0005985  -0.563    0.575    
x2           0.0234097  0.0627251   0.373    0.710    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Gamma family taken to be 0.01376099)

    Null deviance: 1.3498  on 99  degrees of freedom
Residual deviance: 1.3436  on 97  degrees of freedom
AIC: 1240.6

Number of Fisher Scoring iterations: 4

SAS的输出:

Answer 1

默认情况下，R 计算scale=1/分散参数的方式与sas/genmod/model 选项scale=pearson 相同。尺度参数的选择会影响 SE。请参阅此处的文档： https://support.sas.com/documentation/cdl/en/statug/63033/HTML/default/viewer.htm#statug_genmod_sect022.htm

默认情况下，SAS/genmod 给出形状参数的 MLE。假设拟合的伽玛模型存储在列表“fit”中。要在 R 中获取此内容，请加载 MASS 库，然后键入

gamma.shape(fit)

这给出了形状参数 alpha 的 MLE。如果您随后输入

summary(fit, dispersion=1/gamma.shape(fit)$alpha))

摘要函数在计算 SE 时将使用 alpha 的 MLE，并且它们将与 SAS/genmod 完全匹配。

我将就此单独发表一篇文章。虽然summary.glm给出了正确的SE（使用指定的离差值），但它不会打印正确的AIC值（它不使用指定的离差，而是使用通过皮尔逊残差计算的值）。差异很小，但我将其称为错误。