使用似然法提取svyglm风险差异置信区间时遇系数有效性错误

我尝试用带恒等链接（identity link）的二项族GLM拟合调查加权数据以计算风险差异，已通过提供初始值成功运行svyglm模型，但使用似然法提取不对称置信区间时，收到错误提示：

Error: no valid set of coefficients has been found: please supply starting values

样本数据（前30条） #

data <- structure(list(binary_outcome= structure(c(1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("0", "1"), class = "factor"), 
    binary_exposure= structure(c(1L, 1L, NA, 1L, 1L, NA, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, NA, 1L, NA, 1L, 
    NA, NA, NA, NA, 1L, NA, 2L, NA), .Label = c("1", "2"), label = "binary exposure", class = "factor"), 
    wt = c(191217, 89001, 16578, 188901, 7787, 7787, 3876, 
    27569, 7208, 7208, 8594, 170173, 3430, 7208, 53034, 53034, 
    5037, 5037, 786465, 177210, 7849, 9038, 6353, 3837, 89733, 
    188973, 176366, 188736, 3456, 27478), block = c("a", 
    "a", "b", "c", "w", "e", 
    "j", "c", "g", "r", "j", 
    "r", "v", "f", "g", "d", 
    "a", "a", "v", "t", "l", 
    "q", "d", "n", 
    "m", "b", "g", "k", "m", 
    "g"), strata_study = c("C R", "K U", 
    "TN R", "WB R", "K R", "K R", 
    "TN U", "WB U", "TN U", "TN U", 
    "O U", "M U", "K U", "TN U", 
    "WB U", "WB U", "K R", "K R", 
    "D U", "WB R", "TN U", "TN U", 
    "TN U", "K U", "D U", "WB R", "C R", 
    "WB R", "K U", "WB U")), class = c("grouped_df", 
"tbl_df", "tbl", "data.frame"), row.names = c(NA, -30L), groups = structure(list(
    binary_outcome = structure(1:2, .Label = c("0", "1"), class = "factor"), 
    .rows = structure(list(1:15, 16:30), ptype = integer(0), class = c("vctrs_list_of", 
    "vctrs_vctr", "list"))), row.names = c(NA, -2L), class = c("tbl_df", 
"tbl", "data.frame"), .drop = TRUE))

拟合模型代码 #

svy_design <- svydesign(ids = ~block, strata = ~strata_study, weights = ~wt, data = data, nest = TRUE, fpc = NULL)

summary(svyglm(I(binary_outcome==1)~binary_exposure, design=svy_design,
               family=quasi(link="identity",variance="mu(1-mu)"),  
               start=c(0.1,0.6)))

模型运行结果 #

Call:
svyglm(formula = I(binary_outcome == 1) ~ binary_exposure, design = svy_design, 
    family = quasi(link = "identity", variance = "mu(1-mu)"), 
    start = c(0.1, 0.1))

Survey design:
svydesign(ids = ~block, strata = ~strata_study, weights = ~wt, 
    data = data, nest = TRUE, fpc = NULL)

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)  
(Intercept)        0.3575     0.2043   1.750   0.1107  
binary_exposure2   0.6425     0.2043   3.145   0.0104 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for quasi family taken to be 0.7467369)

Number of Fisher Scoring iterations: 2

出错的置信区间调用代码 #

confint(svyglm(I(binary_outcome==1)~binary_exposure, design=svy_design,
               family=quasi(link="identity",variance="mu(1-mu)"),  
               start=c(0.1,0.6)), method="likelihood")

解决方案 #

报错原因是confint()在内部重新拟合模型时，没有继承你传入的初始值参数。解决方式是先将拟合好的模型保存为对象，再对该对象调用confint()：

# 先拟合模型并保存
fit <- svyglm(I(binary_outcome==1)~binary_exposure, design=svy_design,
              family=quasi(link="identity",variance="mu(1-mu)"),  
              start=c(0.1,0.6))

# 对已拟合完成的模型调用似然法置信区间
confint(fit, method="likelihood")

另外需要注意：恒等链接的二项模型本身容易出现拟合不稳定（预测值可能超出[0,1]范围），初始值的选择对收敛至关重要，你当前使用的初始值已经能让模型收敛，复用该初始值即可确保confint内部的拟合过程顺利完成。

内容的提问来源于stack exchange，提问作者s.stats