r - GLM 模型在交互式代码中运行，但在我使用 knitr 时不运行

Question

我在使用 knitr 时遇到问题。具体来说，我有一个模型在控制台中运行得非常好，但是当我尝试编织文档时，R 会抛出一个错误。

加载数据集(available here to facilitate replication)

scabies <- read.csv(file = "S1-Dataset_CSV.csv", header = TRUE, sep = ",")
scabies$agegroups <- as.factor(cut(scabies$age, c(0,10,20,Inf), labels = c("0-10","11-20","21+"), include.lowest = TRUE)) 
scabies$agegroups <-relevel(scabies$agegroups, ref = "21+")
scabies$house_cat <- as.factor(cut(scabies$house_inhabitants, c(0,5,10,Inf), labels = c("0-5","6-10","10+"), include.lowest = TRUE))
scabies$house_cat <- relevel(scabies$house_cat, ref = "0-5")
scabies <- scabies %>% mutate(scabies = case_when(scabies_infestation=="yes"~1,
                                                  scabies_infestation=="no"~0)) %>%
                      mutate(impetigo = case_when(impetigo_active=="yes" ~1,
                                                  impetigo_active=="no" ~0))

拟合模型

scabiesrisk <- glm(scabies~agegroups+gender+house_cat,data=scabies,family=binomial())
scabiesrisk_OR <- exp(cbind(OR= coef(scabiesrisk), confint(scabiesrisk)))
scabiesrisk_summary <- summary(scabiesrisk)
scabiesrisk_summary <- cbind(scabiesrisk_OR, scabiesrisk_summary$coefficients)
scabiesrisk_summary

此代码在控制台中运行得非常好。但是当我尝试 knitr 我得到:

Error in model.frame.default(formula = scabies ~ agegroups + gender +: invalid type(list) for variable 'scabies Calls: ... glm-> eval -> eval -> -> model.frame.default

Answer 1

我能够重现您描述的问题，但还没有完全理解幕后发生的事情。
这个 Markdown block 很有趣:

```{r}
scabiesrisk_OR <- exp(cbind(OR= coef(scabiesrisk), confint((scabiesrisk))))
scabiesrisk_summary <- summary(scabiesrisk)
scabiesrisk_summary <- cbind(scabiesrisk_OR, scabiesrisk_summary$coefficients)
scabiesrisk_summary
```

如果我手动快速依次执行 block 中的行(ctrl+Enter x 4)，有时我会收到两条分析消息:

Waiting for profiling to be done...
Waiting for profiling to be done...

在这种情况下，summary(scabiesrisk) 是一个矩阵:

> class(scabiesrisk_summary)
[1] "matrix" "array" 

如果我手动慢慢执行 block 中的行，我只会收到一条分析消息:

Waiting for profiling to be done...

summary(scabiesrisk) 是一个 summary.glm :

> class(scabiesrisk_summary)
[1] "summary.glm"

看起来分析是在单独的线程上启动的，并且根据它是否完成，摘要函数没有相同的行为。如果分析完成，它会返回预期的 summary.glm 对象，但如果不是这样，它会启动另一个分析并返回一个矩阵。
特别是，对于矩阵 scabiesrisk_summary$coefficients 不可用，我在这种情况下收到以下错误消息:

Error in scabiesrisk_summary$coefficients : 
  $ operator is invalid for atomic vectors

编织时也可能发生这种情况:编织开销是否会使分析速度变慢从而出现问题？

找到解决方法 here (使用 confint.default 而不是 confint)，我无法重现上述问题:

scabiesrisk_OR <- exp(cbind(OR= coef(scabiesrisk), confint.default((scabiesrisk))))
scabiesrisk_summary <- summary(scabiesrisk)
scabiesrisk_summary <- cbind(scabiesrisk_OR, scabiesrisk_summary$coefficients)
scabiesrisk_summary
                       OR      2.5 %    97.5 %   Estimate Std. Error
(Intercept)    0.09357141 0.06984512 0.1253575 -2.3690303  0.1492092
agegroups0-10  2.20016940 1.60953741 3.0075383  0.7885344  0.1594864
agegroups11-20 2.53291768 1.79985894 3.5645415  0.9293719  0.1743214
gendermale     1.44749159 1.13922803 1.8391682  0.3698321  0.1221866
house_cat6-10  1.30521927 1.02586104 1.6606512  0.2663710  0.1228792
house_cat10+   1.17003712 0.67405594 2.0309692  0.1570355  0.2813713
                   z value     Pr(>|z|)
(Intercept)    -15.8772359 9.110557e-57
agegroups0-10    4.9442116 7.645264e-07
agegroups11-20   5.3313714 9.747386e-08
gendermale       3.0267824 2.471718e-03
house_cat6-10    2.1677478 3.017788e-02
house_cat10+     0.5581076 5.767709e-01

所以你也可以试试这个。

与直接可读的 R 函数 confint.defaut 不同，confint 是 S3 调度方法(感谢 @Ben Bolker 注释中的内部引用)，以及我还没有进一步调查什么可以解释这种令人惊讶的行为。

另一个选项似乎将 scabiesrisk_summary 保存在另一个变量中。
我很努力，但这样做后永远无法重现该问题:

```{r}
scabiesrisk_OR <- exp(cbind(OR= coef(scabiesrisk), confint((scabiesrisk))))
scabiesrisk_summary <- summary(scabiesrisk)
scabiesrisk_final <- cbind(scabiesrisk_OR, scabiesrisk_summary$coefficients)
scabiesrisk_final
```