在group_by()中使用mutate()的lm()

Question

我希望找到一种方法，可以将不同级别的分别计算的lm(a~b)函数的残差添加到我的数据表中作为一列。

有人建议我使用sort_by(c)函数，但这似乎无法与lm(a~b)一起使用。

我的工作示例数据如下：

列名为subject、trial和rt的数据框（data.frame）中，我的目标是从一个R函数计算出我最初在SPSS中创建的Zre_SPSS。

我已经尝试过：

data %<>% group_by (subject) %>% 
  mutate(Zre=residuals(lm(log(rt)~trial)))

但它不起作用- Zre得到计算，但不是分别针对每个主题进行计算，而是针对整个数据框。

请问有人能帮我吗？我完全是R（以及编码）的新手，所以如果这个问题很蠢或者是duplicate，请原谅我，可能是我没有理解其他解决方案，或者它们不是我寻找的解决方案。最好的祝福。

按照Ben Bolker的要求，这里是生成来自Excel屏幕截图数据的R代码

#generate data
  subject<-c(1,1,1,1,1,1,2,2,2,2,2,2,3,3,3,3,3,3)
  subject<-factor(subject)
  trial<-c(1,2,3,4,5,6,1,2,3,4,5,6,1,2,3,4,5,6)
  rt<-c(300,305,290,315,320,320,350,355,330,365,370,370,560,565,570,575,560,570)

#Following variable is what I would get after using SPSS code
  ZreSPSS<-c(0.4207,0.44871,-1.7779,0.47787,0.47958,-0.04897,0.45954,0.45487,-1.7962,0.43034,0.41075,0.0407,-0.6037,0.0113,0.61928,1.22038,-1.32533,0.07806)

#make data frame
  sym<-data.frame(subject, trial, rt, ZreSPSS)

- blazej

---

你能否以文本/可复制的形式发布你的示例，而不是作为屏幕截图呢？ - Ben Bolker

你可能想要查看 tidyr::nest 和一个快速博客简介。 - r2evans

这可能会有所帮助。 - Haboryme

一个在 RHS 上带有 trial 而不是作为因子的模型似乎没有太多意义。这真的是您在 SPSS 模型上拟合的数据吗？ - Hong Ooi

@HongOoi 是的，我相信我说得没错 - 一般来说想要从每个参与者的试验次数 ln(log(rt)~trial) 预测 log(rt)（延迟），并保存对应于每个 rt 的残差。 - blazej

Answer 1

在dplyr 0.5版本的mutate中存在一个bug，当在组内使用 lm 时，它仍然会尝试使用整个数据集。您可以改用do来解决：

sym %>% group_by(subject) %>% do(
{
    r <- resid(lm(log(rt) ~ trial, data = .))
    data.frame(., r)
})

这个结果仍然不匹配你的SPSS列，但是它是根据你提供的数据得出的正确结果。你可以通过为每个受试者手动拟合模型并检查残差来验证这一点。

（其他类型的残差包括标准化残差rstandard和学生化残差rstudent。它们仍然不符合你的SPSS数字，但可能是你要寻找的内容。）

- Hong Ooi

---

感谢您的建议（Zre_SPSS实际上是rstandard）。您的代码可以运行，但正如您已经看到的，它与SPSS输出不匹配。我仔细检查了数据，发现完全相同，这太奇怪了！我应该编辑我的问题并发布带有一些语法参考的SPSS数据文件吗？ - blazej

我知道你肯定不喜欢截屏，但这是一个 SPSS 界面的证明[链接]http://imgur.com/a/JHRPP[链接]。 - blazej

我刚刚发现的是，你的公式 rstandard（lm(log(rt) ~ trial)与SPSS输出的学生化残差匹配，但与上面我的截图中的标准化残差不匹配。 - blazej

替代 do(...) 可以尝试使用 mutate(resid = resid(lm(data.frame(log(rt), trial))))。 - G. Grothendieck

@blazej，你是怎么做到让公式的输出结果与SPSS给出的一致的？ - user3553260

我没有匹配成功 - SPSS 计算 ZRE 分数的方式与 R 中的 rstandard 不同。 在这里看一下：http://stackoverflow.com/questions/40062482/standarized-residuals-in-spss-not-maching-r-rstandardlm - blazej

Answer 2

较新版本的dplyr似乎能够处理这个问题（已经测试过dplyr0.7.4）：

sym %>% group_by(subject) %>% do(
{
    r <- resid(lm(log(rt) ~ trial, data = .))
    data.frame(., r)
}) ->a

sym %>% group_by(subject) %>% mutate(

    r =  resid(lm(log(rt) ~ trial))
) ->b

all(a$r==b$r)  #->TRUE

另一个独立测试

# https://dev59.com/2pvga4cB1Zd3GeqP5KdU#40061201
# https://dev59.com/0YHba4cB1Zd3GeqPSozK
# https://github.com/tidyverse/dplyr/issues/2177

# tested with dplyr 0.7.4

# 1) do 
df = group_by(iris,Species) %>% do({
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width, data=.) )
data.frame(., res)
})

# 2) group_by + mutate
# cannot have "data=." in lm
df2 = group_by(iris,Species) %>% mutate(
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width) )
)

# 3) filter + mutate
df3 = filter(iris,Species=='setosa') %>% mutate(
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width, data=.) )
)
df3 = bind_rows(df3,
filter(iris,Species=='versicolor') %>% mutate(
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width, data=.) )
))
df3 = bind_rows(df3,
filter(iris,Species=='virginica') %>% mutate(
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width, data=.) )
))

# 4) across all rows (should not be the same)
df4 = mutate(iris,
res = resid( lm(Sepal.Length~Petal.Length+Petal.Width, data=iris) )
)

# conclusion: all the same, except df4
all(df$res==df2$res)
all(df$res==df3$res)
df$res==df4$res

- Jerry T