无法使用因子变量在 svyglm 中指定起始 glm 参数，R 调查包

Question

我正在使用加权分析并使用svyglm 来分析来自复杂加权方案的无回复数据。我想通过将binomial(link=log) 指定为家庭来拟合对数二项式模型来估计在大多数情况下适合的流行率。但是，在默认拟合器无法找到一组起始系数的情况下，我发现在大多数情况下都可以使用的方便设置是设置Start <- c(log(mean(response.var)), rep(0, ncov))。

当我将start 提供给survey 包中的svyglm 函数时，R 抛出了一个我似乎无法解析的错误。似乎只要协变量之一是一个因素。

例子：

library(survey)
data(api)
apistrat$qmeal <- with(apistrat, cut(meals, quantile(meals)))
dstrat<-svydesign(id=~1,strata=~stype, weights=~pw, data=apistrat, fpc=~fpc)

还有一个有问题的 GLM 示例，建模一些荒谬的东西来重现错误：

> svyglm(awards ~ qmeal +emer, family=quasibinomial(link=log), design=dstrat)
Error: no valid set of coefficients has been found: please supply starting values

好的...所以我指定：Start <- c(log(mean(api$awards, na.rm=T)), 0, 0, 0, 0)

> svyglm(awards ~ cut(meals, quantile(meals)) +emer, family=quasibinomial, design=dstrat, start=start)

 > svyglm(awards ~ qmeal +emer, family=quasibinomial(link=log), design=dstrat, start=start)
Error in glm.fit(x = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,  : 
  length of 'start' should equal 5 and correspond to initial coefs for c("(Intercept)", "qmeal(20.8,39.5]", "qmeal(39.5,69]", "qmeal(69,100]", , "emer")

有趣的是，start 的长度为 5。我进一步注意到，svyglm 始终会产生一个额外的,（在最后一个 qmeal 变量和“emer”之间查找），但缺少条目。这在提供给标准glm时没有这样的问题：

glm(awards ~ qmeal +emer, family=quasibinomial(link=log), data=apistrat, start=start)

产生正确的输出：

Call:  glm(formula = awards ~ qmeal + emer, family = quasibinomial(link = log), 
    data = apistrat, start = start)

Coefficients:
     (Intercept)  qmeal(20.8,39.5]    qmeal(39.5,69]     qmeal(69,100]              emer  
        -0.59276           0.13058           0.31311           0.24698          -0.01389  

Degrees of Freedom: 198 Total (i.e. Null);  194 Residual
  (1 observation deleted due to missingness)
Null Deviance:      272.7 
Residual Deviance: 265.7    AIC: NA

Answer 1

您的代码实际上并没有运行，但如果我这样做了 start <- c(log(mean(apistrat$awards=="Yes", na.rm=T)), 0, 0, 0, 0)

我确实收到了您所询问的错误。发生这种情况是因为对glm 的调用（故意）在设计对象中查找它的参数，然后在svyglm（不是故意）内部查找它的参数。形式参数start 在那里不可见。但是start 是stats 包中的一个函数的名称，由于复杂的原因[1]，这就是你得到的参数。它没有长度 5。杂散的逗号是虚假的[2]

解决方法是明确指定 start 参数而不是变量，因此不需要查找 svyglm(awards ~ qmeal +emer, family=quasibinomial(link=log), design=dstrat, start=c(log(mean(apistrat$awards=="Yes", na.rm=T)), 0, 0, 0, 0))

另一种解决方法是指定不是现有函数名称的内容，以便查找进入下一个级别并找到您的变量。例如，这些都对我有用： initial <- c(log(mean(apistrat$awards=="Yes", na.rm=T)), 0, 0, 0, 0) svyglm(awards ~ qmeal +emer, family=quasibinomial(link=log), design=dstrat, start=initial) rose <- c(log(mean(apistrat$awards=="Yes", na.rm=T)), 0, 0, 0, 0) svyglm(awards ~ qmeal +emer, family=quasibinomial(link=log), design=dstrat, start=rose) 我会尝试在下一个版本中解决这个问题。

[1] 不，如果不进行大量实验，我再精确不过了。它们很复杂。

[2] glm.fit 使用deparse 将名称向量转换为字符串，并且向量足够长，可以跨越两行，这就是逗号的来源。如果变量被称为 m 而不是 qmeal 你不会得到逗号。您可能会争辩说这是 glm.fit 中的一个错误，但您可能会付出更多的努力。

Answer 2

glm 在svyglm 中的调用方式似乎存在问题。将向量名称start 替换为与svyglm 的参数名称不匹配的任何内容（例如x）即可解决问题。