r - Stata vs. R：Delta 方法为 logit 模型的相对风险 SE 提供了不同的结果

Question

我一直在尝试估计 logit 回归（使用相对风险）中协变量的条件平均处理效果以及用于推理目的的标准误差。delta 方法对于计算这些处理效果的标准误差是必要的。我一直在尝试使用 Stata 用户编写的命令重新创建这些结果adjrr，以计算 R 中标准误差和置信区间的相对风险。

Stata 中的adjrr命令计算调整后的相对风险（我的项目感兴趣的条件平均处理效果），它是使用 delta 方法的 SE。R 中的deltamethod命令应该会产生相同的结果，但事实并非如此。

如何在 R 中复制来自 Stata 的结果？

我使用了以下自行生成的数据：（https://migariane.github.io/DeltaMethodEpiTutorial.nb.html）。

R代码如下：

generateData <- function(n, seed){
    set.seed(seed)
    age <- rnorm(n, 65, 5)
    age65p <- ifelse(age>=65, T, F)
    cmbd <- rbinom(n, size=1, prob = plogis(1 - 0.05 * age))
    Y <- rbinom(n, size=1, prob = plogis(1 - 0.02* age - 0.02 * cmbd))
    data.frame(Y, cmbd, age, age65p)
}

# Describing the data
data <- generateData(n = 1000, seed = 777) 
str(data)

logfit <- glm(Y ~ age65p + cmbd, data = data, family = binomial)
summary(logfit)

p1 <- predict(logfit, newdata = data.frame(age65p = T, cmbd = 0), type="response")
p0 <- predict(logfit, newdata = data.frame(age65p = F, cmbd = 0), type="response")
rr <- p1 / p0
rr
0.8123348 #result

library(msm)
se_rr_delta <- deltamethod( ~(1 + exp(-x1)) / (1 + exp(-x1 -x2)), coef(logfit), vcov(logfit))
se_rr_delta
0.6314798 #result

Stata代码（使用相同的数据）：

logit Y i.age65p i.cmbd

adjrr age65p

//results below
R1  = 0.3685 (0.0218)    95% CI  (0.3259, 0.4112)
R0  = 0.4524 (0.0222)    95% CI  (0.4090, 0.4958)
ARR = 0.8146 (0.0626)    95% CI  (0.7006, 0.9471)
ARD = -0.0839 (0.0311)    95% CI  (-0.1449, -0.0229)
p-value (R0 = R1):  0.0071
p-value (ln(R1/R0) = 0):  0.0077

标签： rstata