R中的嵌套表和计算带有置信区间的汇总统计

Question

这个问题是关于统计程序R。

数据

我有一个数据框 study_data，它有 100 行，每行代表不同的人，三列 gender、height_category 和 freckles。变量gender 是一个因子，取“男性”或“女性”的值。变量height_category 也是一个因子，取“tall”或“short”的值。变量freckles 是一个连续的数值变量，表示个体有多少雀斑。

以下是一些示例数据（感谢 Roland）：

set.seed(42)
DF <- data.frame(gender=sample(c("m","f"),100,T),
      height_category=sample(c("tall","short"),100,T),
      freckles=runif(100,0,100))

问题 1

我想创建一个嵌套表，将这些患者分为“男性”与“女性”，进一步细分为“高”与“矮”，然后计算每个子组中的患者数量以及具有上下 95% 置信区间的雀斑中位数。

示例

该表应如下所示，其中# 符号被替换为适当的计算结果。

gender height_category n median_freckles LCI UCI

male              tall #               #   #   #
                 short #               #   #   #
female            tall #               #   #   #
                 short #               #   #   #

问题 2

计算完这些结果后，我想创建一个条形图。 y 轴将是雀斑的中位数。 x 轴将分为男性和女性。但是，这些部分将按高度类别细分（因此总共有四个条形，每两个一组）。我想在条形顶部叠加 95% 的置信区间。

我尝试过的

我知道我可以使用MASS 库和xtabs 命令制作嵌套表：

ftable(xtabs(formula = ~ gender + height_category, data = study_data))

但是，我不确定如何将计算雀斑数量的中位数合并到此命令中，然后让它显示在汇总表中。我也知道ggplot2 可用于制作条形图，但我不知道如何做到这一点，因为我无法计算出我首先需要的数据。

Answer 1

您确实应该提供一个可重现的示例。无论如何，您可能会发现 library(plyr) 很有帮助。请注意这些置信区间，因为如果 n

library(plyr)
ddply(df, .(gender, height_category), summarize,
      n=length(freckles), median_freckles=median(freckles), 
      LCI=qt(.025, df=length(freckles) - 1)*sd(freckles)/length(freckles)+mean(freckles), 
      UCI=qt(.975, df=length(freckles) - 1)*sd(freckles)/length(freckles)+mean(freckles))

编辑：我忘了在情节上添加一点。假设我们将之前的结果保存为tab：

library(ggplot2)
library(reshape)
m.tab <- melt(tab, id.vars=c("gender", "height_category"))
dodge <- position_dodge(width=0.9)
ggplot(m.tab, aes(fill=height_category, x=gender, y=median_freckles))+
  geom_bar(position=dodge) + geom_errorbar(aes(ymax=UCI, ymin=LCI), position=dodge, width=0.25)

Answer 2

set.seed(42)
DF <- data.frame(gender=sample(c("m","f"),100,T),
                 height_category=sample(c("tall","short"),100,T),
                 freckles=runif(100,0,100))


library(plyr)
res <- ddply(DF,.(gender,height_category),summarise,
             n=length(na.omit(freckles)),
             median_freckles=quantile(freckles,0.5,na.rm=TRUE),
             LCI=quantile(freckles,0.025,na.rm=TRUE),
             UCI=quantile(freckles,0.975,na.rm=TRUE))

library(ggplot2)
p1 <- ggplot(res,aes(x=gender,y=median_freckles,ymin=LCI,ymax=UCI,
                     group=height_category,fill=height_category)) +
  geom_bar(stat="identity",position="dodge") +
  geom_errorbar(position="dodge")
print(p1)

#a better plot that doesn't require to precalculate the stats
library(hmisc)
p2 <- ggplot(DF,aes(x=gender,y=freckles,colour=height_category)) + 
  stat_summary(fun.data="median_hilow",geom="pointrange",position = position_dodge(width = 0.4))
print(p2)