代码之家 › 专栏 › 技术社区 › Jake

更新函数中的模型因子级别时工作不正常

function r

Jake · 技术社区 · 7 年前

我在R中运行了一些GLM模型,这些模型与我正在进行的喂养试验有关。我用两个预测因子回归我感兴趣的变量:一个因素有三个水平,一个连续变量。我想将每一级的拦截因素相互比较,以确定它们是否不同。为此,我编写了一个函数(在下面的可复制代码中称为interceptcarbire),它释放因子并更新模型,然后保存每个模型的结果。这是我对截获的所有成对比较的快速方法。

问题是,当我运行函数时,它似乎没有正确地更新模型。 返回的列表中的每个项都是相同的,当它们应该改变时,每个项都有一个不同的因子级别,与其他级别正在比较的“(截距)”不同。我怀疑这与函数的环境有关,但我不确定。我还没有在StackOverflow或Google上找到类似的例子。

下面是一个可重复的例子:

food <- as.factor(rep(c("a", "b", "c"), each = 20))
variable <- rbinom(60, 1, 0.7)
movement <- rgamma(60, 10, 2)
binomial.model <- glm(variable ~ food,
         family = "binomial")
gamma.model <- glm(movement ~ food,
               family = Gamma)
interceptCompare <- function(model, factor) {
  results <- list() # empty list to store results
  for (i in unique(factor)) {
    factor <- relevel(factor, ref = i)
    model <- update(model)
    results[[i]] <- summary(model)$coefficients[1:3, ]
  }
  results <- lapply(results, function(x) round(x, 4))
  return(results)
}

interceptCompare(binomial.model, food)
interceptCompare(gamma.model, food)

2 回复 | 直到 7 年前

Onyambu 7 年前

您需要添加一行,以便更改数据,并在更新中使用它:

interceptCompare <- function(model, factor) {
  results <- list() # empty list to store results

   s <- deparse(substitute(factor))#ADD THIS LINE

  for (i in unique(factor)) {
    factor <- relevel(factor, ref = i)
    model[["model"]][[s]] <- factor #CHANGE THE DATA IN THE MODEL
    model <- update(model,data=model[["model"]])# UPDATE THE MODEL
    results[[i]] <- summary(model)$coefficients[1:3, ]
  }
  results <- lapply(results, function(x) round(x, 4))
  return(results)
}

interceptCompare(binomial.model, food)
$a
            Estimate Std. Error z value Pr(>|z|)
(Intercept)   1.3863     0.5590  2.4799   0.0131
foodb        -0.7673     0.7296 -1.0516   0.2930
foodc        -0.2877     0.7610 -0.3780   0.7054

$b
            Estimate Std. Error z value Pr(>|z|)
(Intercept)   0.6190     0.4688  1.3205   0.1867
fooda         0.7673     0.7296  1.0516   0.2930
foodc         0.4796     0.6975  0.6876   0.4917

$c
            Estimate Std. Error z value Pr(>|z|)
(Intercept)   1.0986     0.5164  2.1275   0.0334
foodb        -0.4796     0.6975 -0.6876   0.4917
fooda         0.2877     0.7610  0.3780   0.7054


interceptCompare(gamma.model, food)
$a
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   0.2246     0.0156 14.3919   0.0000
foodb        -0.0170     0.0213 -0.8022   0.4257
foodc        -0.0057     0.0218 -0.2608   0.7952

$b
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   0.2076     0.0144 14.3919   0.0000
fooda         0.0170     0.0213  0.8022   0.4257
foodc         0.0114     0.0210  0.5421   0.5898

$c
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   0.2189     0.0152 14.3919   0.0000
foodb        -0.0114     0.0210 -0.5421   0.5898
fooda         0.0057     0.0218  0.2608   0.7952

MrFlick 7 年前

当你试图从你正在做的公式中替换符号时,你需要更加小心。你需要用R语言能理解的术语来表达它。你要以“食物”的名义传递,而不是像现在这样存储在“食物”向量中的值。这是一个更新,它似乎可以完成你想做的事情

interceptCompare <- function(model, factor) {
  sym <- substitute(factor)
  results <- list() # empty list to store results
  for (i in unique(factor)) {
    change <- eval(bquote(~.-.(sym)+relevel(.(sym), ref=.(i))))
    new_model <- update(model, change)
    results[[i]] <- summary(new_model)$coefficients[1:3, ]
  }
  results <- lapply(results, function(x) round(x, 4))
  return(results)
}

在这里,我们用 substitute .然后我们使用 bquote() 建立一个新的公式,该公式将删除原始值,并使用特定的引用重新调整因子变量。然后我们把这个保存到一个新的对象中,这样我们就不会一直更新同一个模型。对于 binomial.model ,这个返回

$`a`
                          Estimate Std. Error z value Pr(>|z|)
(Intercept)                 0.8473     0.4879  1.7364   0.0825
relevel(food, ref = "a")b   0.0000     0.6901  0.0000   1.0000
relevel(food, ref = "a")c  -0.8473     0.6619 -1.2801   0.2005

$b
                          Estimate Std. Error z value Pr(>|z|)
(Intercept)                 0.8473     0.4879  1.7364   0.0825
relevel(food, ref = "b")a   0.0000     0.6901  0.0000   1.0000
relevel(food, ref = "b")c  -0.8473     0.6619 -1.2801   0.2005

$c
                          Estimate Std. Error z value Pr(>|z|)
(Intercept)                 0.0000     0.4472  0.0000   1.0000
relevel(food, ref = "c")a   0.8473     0.6619  1.2801   0.2005
relevel(food, ref = "c")b   0.8473     0.6619  1.2801   0.2005

你可以看到它是如何改变 ref= 每次迭代时