R中for循环按分组应用自定义汇总函数时分组失效问题排查

编写for循环将自定义的summarise汇总函数应用到数据框所有数值列时，循环输出结果忽略分组因子，两组返回的统计值完全相同；不使用for循环、单独对单个列执行该函数时，可以得到正确的分组统计结果。

复现代码 #

# 示例数据框
dexadf <- structure(list(participant = c("pt04", "pt75", "pt21", "pt73", 
"pt27", "pt39", "pt43", "pt52", "pt69", "pt49", "pt50", "pt56", 
"pt62", "pt68", "pt22", "pt64", "pt54", "pt79", "pt36", "pt26", 
"pt65", "pt38"), group = structure(c(1L, 2L, 2L, 1L, 1L, 2L, 
1L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 1L
), .Label = c("c", "e"), class = "factor"), sex = structure(c(2L, 
1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 
2L, 1L, 1L, 2L, 1L), .Label = c("m", "f"), class = "factor"), 
    fm_bdc3 = c(18.535199635968, 23.52996574649, 17.276246451976, 
    11.526088555461, 23.805048656112, 23.08597823716, 28.691020942436, 
    28.968097858499, 23.378093165331, 22.491725344661, 14.609015054932, 
    19.734914019306, 31.947412973684, 25.152298171274, 12.007356801787, 
    20.836128108938, 22.322230884349, 14.777652101515, 21.389572717608, 
    16.992853675086, 14.138189878472, 17.777235203826), fm_rec3 = c(18.545007190636, 
    23.017181869742, 17.031403417007, 11.227201061887, 23.581434653208, 
    21.571120542136, 28.919246372213, 28.138632765662, 22.990408911436, 
    22.274932676852, 14.012586350504, 19.066675709151, 30.897705534847, 
    24.491614222412, 11.670939246332, 20.306494543464, 22.052263684182, 
    14.252973638341, 21.028701096846, 17.207104923059, 13.172159777361, 
    17.610831079442), fm_chg = c(0.00980755466799721, -0.512783876747999, 
    -0.244843034968998, -0.298887493573998, -0.223614002904, 
    -1.514857695024, 0.228225429777002, -0.829465092836998, -0.387684253894999, 
    -0.216792667809003, -0.596428704428, -0.668238310155001, 
    -1.049707438837, -0.660683948862001, -0.336417555455, -0.529633565474001, 
    -0.269967200167002, -0.524678463173998, -0.360871620761998, 
    0.214251247972999, -0.966030101111, -0.166404124383998), 
    fm_percchg = c(0.00052913132097943, -0.0217928016671462, 
    -0.0141722355981437, -0.0259313896588437, -0.00939355370091154, 
    -0.0656180855522784, 0.00795459423472242, -0.0286337438132355, 
    -0.0165832282022865, -0.00963877445980213, -0.0408260722701251, 
    -0.0338607155572751, -0.0328573534170568, -0.0262673392452288, 
    -0.028017619615079, -0.025419001203338, -0.0120940958619099, 
    -0.0355048596062299, -0.0168713805332318, 0.0126083147698213, 
    -0.0683277073949869, -0.00936051767758492)), row.names = c(NA, 
-22L), class = "data.frame")

# 自定义汇总函数
summbygrp <- function(x) {
  group_by(dexadf, group) %>%
        summarise(
          count = n(),
          mean = mean({{x}}, na.rm = TRUE),
          sd = sd({{x}}, na.rm = TRUE)
        ) %>%
        mutate(se = sd / sqrt(11),
               lower.ci = mean - qt(1 - (0.05 / 2), 11 - 1) * se,
               upper.ci = mean + qt(1 - (0.05 / 2), 11 - 1) * se
              ) 
} 

# 循环应用函数
coln = 1
for (col in dexadf) {
  print(colnames(dexadf)[coln])
  coln = coln + 1
  if(is.numeric(col)) {
    print(summbygrp(col))
  } else {next}
}

错误输出示例 #

[1] "fm_bdc3"
# A tibble: 2 × 7
  group count  mean    sd    se lower.ci upper.ci
  <fct> <int> <dbl> <dbl> <dbl>    <dbl>    <dbl>
1 c        11  20.6  5.48  1.65     16.9     24.3
2 e        11  20.6  5.48  1.65     16.9     24.3

两个分组所有统计量完全一致，和实际数据分布不符。

• 循环遍历方式错误：for (col in dexadf) 拿到的col是从数据框中抽离的独立原子向量，已经脱离了和原数据框、分组列的绑定关系，不是dplyr非标准求值能识别的列引用。
• 运算符使用场景错误：{{}}是dplyr提供的、用于在函数内透传未求值列名的拥抱运算符，仅能识别数据掩码环境内的列对象。传入外部独立向量时，它不会按分组切分数据计算，只会对整个向量求统计值，因此两个分组返回的都是全量样本的计算结果，自然完全相同。
• 函数写法存在硬编码问题：函数内部直接调用全局环境的dexadf，没有把数据集作为参数传入，扩展性差，也容易引发环境相关的计算错误。

方案1：修正循环逻辑和函数写法 #

调整函数参数，传入数据集和列名，用.data[[colname]]的方式引用列，循环时遍历列名而非列向量：

library(dplyr)

# 修正后的函数
summbygrp <- function(data, colname) {
  data %>%
    group_by(group) %>%
    summarise(
      count = n(),
      mean = mean(.data[[colname]], na.rm = TRUE),
      sd = sd(.data[[colname]], na.rm = TRUE)
    ) %>%
    mutate(se = sd / sqrt(11),
           lower.ci = mean - qt(1 - (0.05 / 2), 11 - 1) * se,
           upper.ci = mean + qt(1 - (0.05 / 2), 11 - 1) * se
    ) 
}

# 提取所有数值列的列名
num_cols <- names(dexadf)[sapply(dexadf, is.numeric)]

# 遍历列名执行计算
for (col in num_cols) {
  print(col)
  print(summbygrp(dexadf, col))
}

方案2：用tidyverse原生语法替代for循环（推荐） #

不需要写显式循环，用across一次性批量计算所有数值列的统计量，输出结果更规整：

library(dplyr)
library(tidyr)

dexadf %>%
  group_by(group) %>%
  summarise(
    across(where(is.numeric), 
           list(mean = ~mean(.x, na.rm = TRUE),
                sd = ~sd(.x, na.rm = TRUE))
    ),
    count = n()
  ) %>%
  pivot_longer(-c(group, count), names_to = c("var", ".value"), names_sep = "_") %>%
  mutate(se = sd / sqrt(count),
         lower.ci = mean - qt(1 - (0.05 / 2), count - 1) * se,
         upper.ci = mean + qt(1 - (0.05 / 2), count - 1) * se
  )

注：方案2里把硬编码的样本量11改成了动态读取每组count值，后续分组样本量变化时不需要手动改代码，容错性更高。

内容的提问来源于stack exchange，提问作者aguan3