base R编写分组计算均值及置信区间函数的问题排查

错误原因 #

• 循环内使用result[,列号]赋值，是对整列所有行赋值，每次循环都会覆盖全列已有值，最终仅保留最后一次循环的分组计算结果
• 没有提取t.test返回置信区间的上下限元素：tva$conf.int是长度为2的向量，第一位为下限、第二位为上限，直接赋值会触发R的向量循环填充规则，导致行数值交替重复
• 未传入置信水平参数，默认使用95%置信区间，后续封装函数时需要关联传入的alpha参数，设置conf.level = 1 - alpha

修正后的非函数版本代码 #

data(chickwts)

groups <- as.factor(chickwts[,2])
nr_groups <- levels(groups)
alpha <- 0.05 # 对应95%置信区间，可按需调整

result <- matrix(0, nrow = length(nr_groups), ncol = 4)
colnames(result) <- c("Lower CI", "Mean", "Upper CI", "Obs.")
rownames(result) <- levels(groups)

result[,4] <- table(chickwts[,2])
# 给t.test传入置信水平参数
group_test <- by(chickwts[, 1], groups, t.test, conf.level = 1 - alpha)

for (i in levels(groups)) {
  tva <- group_test[[i]]
  row_idx <- which(rownames(result) == i)
  # 分别赋值对应行的三个字段，提取置信区间的上下限
  result[row_idx, 1] <- tva$conf.int[1]
  result[row_idx, 2] <- tva$estimate
  result[row_idx, 3] <- tva$conf.int[2]
}

# 查看结果
print(result)

封装为接收三个参数的函数版本 #

group_ci_cal <- function(X, Y, Z) {
  # 参数校验
  if(!is.numeric(X)) stop("X必须为数值向量")
  if(!is.factor(Y)) Y <- as.factor(Y)
  if(length(X) != length(Y)) stop("X和Y长度必须一致")
  if(Z <=0 | Z >=1) stop("Z必须是0到1之间的alpha值")
  
  groups <- levels(Y)
  n_group <- length(groups)
  
  res <- matrix(NA, nrow = n_group, ncol = 4)
  colnames(res) <- c("Lower CI", "Mean", "Upper CI", "Obs.")
  rownames(res) <- groups
  
  res[,4] <- table(Y)
  test_res <- by(X, Y, t.test, conf.level = 1 - Z)
  
  for (i in groups) {
    tmp <- test_res[[i]]
    idx <- which(rownames(res) == i)
    res[idx,1] <- tmp$conf.int[1]
    res[idx,2] <- tmp$estimate
    res[idx,3] <- tmp$conf.int[2]
  }
  return(res)
}

# 调用测试
data(chickwts)
test_output <- group_ci_cal(X = chickwts[,1], Y = chickwts[,2], Z = 0.05)
print(test_output)

运行修正后的代码，输出的各组均值、置信区间将和实际t.test计算结果一致，不会再出现全列重复的问题。

内容的提问来源于stack exchange，提问作者Helan.ge