r - glm后如何返回分数向量？

Question

我正在尝试从 R 中的概率中获取分数向量(对数似然的一阶导数)。这是我的示例代码:

library(AER) # Affairs data
data(Affairs)
mydata<-Affairs
mydata$affairs<-with(mydata,ifelse(affairs>0,1,affairs)) # convert to 1 and 0 

glm(affairs ~ gender+ age + yearsmarried + children + 
  religiousness+education + rating,
 family = binomial(link = "probit"),data = mydata)

Call:  glm(formula = affairs ~ gender + age + yearsmarried + children + 
    religiousness + education + rating, family = binomial(link = "probit"), 
    data = mydata)

Coefficients:
  (Intercept)     gendermale            age   yearsmarried    childrenyes  religiousness      education         rating  
      0.76416        0.18882       -0.02440        0.05461        0.20807       -0.18609        0.01551       -0.27271  

Degrees of Freedom: 600 Total (i.e. Null);  593 Residual
Null Deviance:      675.4 
Residual Deviance: 610.5    AIC: 626.5

在 Stata 中，这将是 predict anything,score 在 glm 之后。

Answer 1

可以使用 sandwich 包中的 estfun() 方法提取观察到的数据的分数(也称为估计函数)。所以你首先拟合模型

m <- glm(affairs ~ gender + age + yearsmarried + children + religiousness + education + rating,
  family = binomial(link = "probit"), data = mydata)

然后您可以提取观察分数的 n x k 矩阵:

library("sandwich")
s <- estfun(m)
dim(s)
## [1] 601   8

观察值之和基本上为零:

colSums(s)
##   (Intercept)    gendermale           age  yearsmarried   childrenyes 
## -0.0006048446 -0.0001081596  0.0005041310  0.0098581660  0.0005940238 
## religiousness     education        rating 
## -0.0006689870 -0.0147258776 -0.0016111599 

检查前六个分数给出:

head(s)
##    (Intercept) gendermale        age yearsmarried childrenyes
## 4   -0.3789593 -0.3789593 -14.021496   -3.7895934   0.0000000
## 5   -0.1913665  0.0000000  -5.166896   -0.7654660   0.0000000
## 11  -0.7474727  0.0000000 -23.919126  -11.2120903  -0.7474727
## 16  -0.1564706 -0.1564706  -8.918825   -2.3470591  -0.1564706
## 23  -0.5249512 -0.5249512 -11.548926   -0.3937134   0.0000000
## 29  -0.1611754  0.0000000  -5.157613   -0.2417631   0.0000000
##    religiousness education    rating
## 4     -1.1368780 -6.821268 -1.515837
## 5     -0.7654660 -2.679131 -0.765466
## 11    -0.7474727 -8.969672 -2.989891
## 16    -0.7823530 -2.816471 -0.782353
## 23    -1.0499023 -8.924170 -1.574853
## 29    -0.3223508 -2.739982 -0.805877