Summarizing data 汇总数据
当你想按组对数据进行汇总（包括平均值、标准差等）
这里介绍的三种方法可以根据某些指定变量对数据进行分组
并对每组应用汇总函数（如平均值、标准差等）：
1、函数 ddply() ，来自 plyr 包，最容易使用的
2、函数 summarizeBy() ，来自 doBy 包，比较容易使用
3、函数 aggregate() ，来自R base 包，较难使用
示例数据
不同性别受试者，服用阿斯匹林或安慰剂前后的观测值以及变化差值
由性别和条件的组合来分组：F-安慰剂、F-阿司匹林、M-安慰剂和 M-阿司匹林
并希望计算出每个组的：
1） N（个体数）
2） mean of change（变化均值）
3） standard deviation（标准差）
4） standard error of the mean （均值标准误差）
data <- read.table(header=TRUE, text=' subject sex condition before after change 1 F placebo 10.1 6.9 -3.2 2 F placebo 6.3 4.2 -2.1 3 M aspirin 12.4 6.3 -6.1 4 F placebo 8.1 6.1 -2.0 5 M aspirin 15.2 9.9 -5.3 6 F aspirin 10.9 7.0 -3.9 7 F aspirin 11.6 8.5 -3.1 8 M aspirin 9.5 3.0 -6.5 9 F placebo 11.5 9.0 -2.5 10 M placebo 11.9 11.0 -0.9 11 F aspirin 11.4 8.0 -3.4 12 M aspirin 10.0 4.4 -5.6 13 M aspirin 12.5 5.4 -7.1 14 M placebo 10.6 10.6 0.0 15 M aspirin 9.1 4.3 -4.8 16 F placebo 12.1 10.2 -1.9 17 F placebo 11.0 8.8 -2.2 18 F placebo 11.9 10.2 -1.7 19 M aspirin 9.1 3.6 -5.5 20 M placebo 13.5 12.4 -1.1 21 M aspirin 12.0 7.5 -4.5 22 F placebo 9.1 7.6 -1.5 23 M placebo 9.9 8.0 -1.9 24 F placebo 7.6 5.2 -2.4 25 F placebo 11.8 9.7 -2.1 26 F placebo 11.8 10.7 -1.1 27 F aspirin 10.1 7.9 -2.2 28 M aspirin 11.6 8.3 -3.3 29 F aspirin 11.3 6.8 -4.5 30 F placebo 10.3 8.3 -2.0 ')1、使用 ddply()
library(plyr)# ddply参数ddply(data, variables, fun = NULL, ... , progress = 'none' , drop = TRUE , parallel = FALSE)# - data ： 要处理的数据框# - variables ：分组变量，可以是公式，也可以是字符向量# - fun：处理分组数据的函数cdata <- ddply(data, c("sex", "condition"), summarise, N = length(change), # 长度即总数 mean = mean(change), # 求平均值 sd = sd(change), # 求标准差 se = sd / sqrt(N) # 求标准误)cdata#> sex condition N mean sd se#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230#> 2 F placebo 12 -2.058333 0.5247655 0.1514867#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190#> 4 M placebo 4 -0.975000 0.7804913 0.3902456处理丢失的数据：
如果数据中有 NA 值，就需要设定参数 na.rm=TRUE ；
length() 没有 na.rm 这个选项，因此只能 sum(!is.na(...)) 来计算有多少个非 NA 值
# 先放一点NA值进来dataNA <- datadataNA$change[11:14] <- NAcdata <- ddply(dataNA, c("sex", "condition"), summarise, N = sum(!is.na(change)), mean = mean(change, na.rm=TRUE), sd = sd(change, na.rm=TRUE), se = sd / sqrt(N))cdata#> sex condition N mean sd se#> 1 F aspirin 4 -3.425000 0.9979145 0.4989572#> 2 F placebo 12 -2.058333 0.5247655 0.1514867#> 3 M aspirin 7 -5.142857 1.0674848 0.4034713#> 4 M placebo 3 -1.300000 0.5291503 0.30550502、函数 summarizeBy()
library(doBy)cdata <- summaryBy(change ~ sex + condition, data=data, FUN=c(length,mean,sd))cdata#> sex condition change.length change.mean change.sd#> 1 F aspirin 5 -3.420000 0.8642916#> 2 F placebo 12 -2.058333 0.5247655#> 3 M aspirin 9 -5.411111 1.1307569#> 4 M placebo 4 -0.975000 0.7804913# 重命名 change.length 这一列为 Nnames(cdata)[names(cdata)=="change.length"] <- "N"# 计算标准误cdata$change.se <- cdata$change.sd / sqrt(cdata$N)cdata#> sex condition N change.mean change.sd change.se#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230#> 2 F placebo 12 -2.058333 0.5247655 0.1514867#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190#> 4 M placebo 4 -0.975000 0.7804913 0.3902456处理丢失数据：
na.rm=TRUE 可以传递给 summaryBy() 调用的每个函数，除了 length()
# 写一个新版本的 length 让它可以识别 na.rm=Tlength2 <- function (x, na.rm=FALSE) { if (na.rm) sum(!is.na(x)) else length(x)}dataNA <- datadataNA$change[11:14] <- NAcdataNA <- summaryBy(change ~ sex + condition, data=dataNA, FUN=c(length2, mean, sd), na.rm=TRUE)cdataNA#> sex condition change.length2 change.mean change.sd#> 1 F aspirin 4 -3.425000 0.9979145#> 2 F placebo 12 -2.058333 0.5247655#> 3 M aspirin 7 -5.142857 1.0674848#> 4 M placebo 3 -1.300000 0.5291503扩充：安排一步得出均值、计数、标准差、标准误和置信区间的函数summarySE
写函数：
## 对数据进行汇总## 得到 计数, 均值, 标准差, 标准误 和 置信区间 (95%).## data: 数据框## measurevar: 需要进行汇总的列## groupvars: 提供分组的列## na.rm: 一个布尔值，决定是否忽略NA## conf.interval: 置信区间的百分比范围 (默认是 95%)summarySE <- function(data=NULL, measurevar, groupvars=NULL, na.rm=FALSE, conf.interval=.95, .drop=TRUE) { library(plyr) # 一个新版本的 length 让它可以识别 na.rm=T length2 <- function (x, na.rm=FALSE) { if (na.rm) sum(!is.na(x)) else length(x) } # 这一步开始汇总. datac <- ddply(data, groupvars, .drop=.drop, .fun = function(xx, col) { c(N = length2(xx[[col]], na.rm=na.rm), mean = mean (xx[[col]], na.rm=na.rm), sd = sd (xx[[col]], na.rm=na.rm) ) }, measurevar ) # 重命名 均值 这一列 datac <- rename(datac, c("mean" = measurevar)) datac$se <- datac$sd / sqrt(datac$N) # 计算标准误 # 标准误的置信区间乘数 # 计算置信区间的t统计量: # 例如, 置信区间是.95, 则使用 .975 (高于/低于), 并使用 df=N-1 ciMult <- qt(conf.interval/2 + .5, datac$N-1) datac$ci <- datac$se * ciMult return(datac)}示例计算：一部得到以上所有汇总信息
summarySE(data, measurevar="change", groupvars=c("sex", "condition"))#> sex condition N change sd se ci#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230 1.0731598#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190 0.8691767#> 4 M placebo 4 -0.975000 0.7804913 0.3902456 1.2419358# 当有 NA 值时, 使用 na.rm=TRUEsummarySE(dataNA, measurevar="change", groupvars=c("sex", "condition"), na.rm=TRUE)#> sex condition N change sd se ci#> 1 F aspirin 4 -3.425000 0.9979145 0.4989572 1.5879046#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 7 -5.142857 1.0674848 0.4034713 0.9872588#> 4 M placebo 3 -1.300000 0.5291503 0.3055050 1.3144821用NA填充空白组合：
有时，汇总数据框中会出现空的因子组合，即有可能出现的因素组合，但实际上并没有出现在原始数据框架中。
在汇总数据框中自动用 NA 填入这些组合通常很有用。
为此，在调用 ddply() 或 summarySE() 时设置 .drop=FALSE 。
用法示例：
# 首先把 男性+安慰剂 组合的数据删去dataSub <- subset(data, !(sex=="M" & condition=="placebo"))# 这样我们汇总数据时，男性+安慰剂 组合就会有空缺列summarySE(dataSub, measurevar="change", groupvars=c("sex", "condition"))#> sex condition N change sd se ci#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230 1.0731598#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190 0.8691767# 设置 .drop=FALSE 不放弃这些组合summarySE(dataSub, measurevar="change", groupvars=c("sex", "condition"), .drop=FALSE)#> Warning in qt(conf.interval/2 + 0.5, datac$N - 1): NaNs produced#> sex condition N change sd se ci#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230 1.0731598#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190 0.8691767#> 4 M placebo 0 NaN NA NA NA用0填充空白组合：
原理同上，但改用0填充
写用0填充的函数：
fillMissingCombs <- function(df, factors, measures) { # 列出因子level的组合列表 levelList <- list() for (f in factors) { levelList[[f]] <- levels(df[,f]) } fullFactors <- expand.grid(levelList) dfFull <- merge(fullFactors, df, all.x=TRUE) # 如果测量值中有 NA ，就用 0 代替 for (m in measures) { dfFull[is.na(dfFull[,m]), m] <- 0 } return(dfFull)}用法举例：
# 首先把 男性+安慰剂 组合的数据删去dataSub <- subset(data, !(sex=="M" & condition=="placebo"))cdataSub <- summarySE(dataSub, measurevar="change", groupvars=c("sex", "condition"))cdataSub#> sex condition N change sd se ci#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230 1.0731598#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190 0.8691767# 用0填充空白组合fillMissingCombs(cdataSub, factors=c("sex","condition"), measures=c("N","change","sd","se","ci"))#> sex condition N change sd se ci#> 1 F aspirin 5 -3.420000 0.8642916 0.3865230 1.0731598#> 2 F placebo 12 -2.058333 0.5247655 0.1514867 0.3334201#> 3 M aspirin 9 -5.411111 1.1307569 0.3769190 0.8691767#> 4 M placebo 0 0.000000 0.0000000 0.0000000 0.00000003、使用R自带的aggregate()
缺点是步骤较为繁琐
优点是不用安装新的R包
# 统计每个类别（性别*条件）中的受试者人数cdata <- aggregate(data["subject"], by=data[c("sex","condition")], FUN=length)cdata#> sex condition subject#> 1 F aspirin 5#> 2 M aspirin 9#> 3 F placebo 12#> 4 M placebo 4# 重命名 "subject" 这一列为 "N"names(cdata)[names(cdata)=="subject"] <- "N"cdata#> sex condition N#> 1 F aspirin 5#> 2 M aspirin 9#> 3 F placebo 12#> 4 M placebo 4# 按 sex 排序cdata <- cdata[order(cdata$sex),]cdata#> sex condition N#> 1 F aspirin 5#> 3 F placebo 12#> 2 M aspirin 9#> 4 M placebo 4# 同时保留 __before__ 和 __after__ 这两列:# 按性别和条件获取平均效应大小cdata.means <- aggregate(data[c("before","after","change")], by = data[c("sex","condition")], FUN=mean)cdata.means#> sex condition before after change#> 1 F aspirin 11.06000 7.640000 -3.420000#> 2 M aspirin 11.26667 5.855556 -5.411111#> 3 F placebo 10.13333 8.075000 -2.058333#> 4 M placebo 11.47500 10.500000 -0.975000# 合并数据框cdata <- merge(cdata, cdata.means)cdata#> sex condition N before after change#> 1 F aspirin 5 11.06000 7.640000 -3.420000#> 2 F placebo 12 10.13333 8.075000 -2.058333#> 3 M aspirin 9 11.26667 5.855556 -5.411111#> 4 M placebo 4 11.47500 10.500000 -0.975000# 获取 "change"的样本 (n-1) 标准差cdata.sd <- aggregate(data["change"], by = data[c("sex","condition")], FUN=sd)# 重命名为 change.sdnames(cdata.sd)[names(cdata.sd)=="change"] <- "change.sd"cdata.sd#> sex condition change.sd#> 1 F aspirin 0.8642916#> 2 M aspirin 1.1307569#> 3 F placebo 0.5247655#> 4 M placebo 0.7804913# 合并cdata <- merge(cdata, cdata.sd)cdata#> sex condition N before after change change.sd#> 1 F aspirin 5 11.06000 7.640000 -3.420000 0.8642916#> 2 F placebo 12 10.13333 8.075000 -2.058333 0.5247655#> 3 M aspirin 9 11.26667 5.855556 -5.411111 1.1307569#> 4 M placebo 4 11.47500 10.500000 -0.975000 0.7804913# 最后计算标准误cdata$change.se <- cdata$change.sd / sqrt(cdata$N)cdata#> sex condition N before after change change.sd change.se#> 1 F aspirin 5 11.06000 7.640000 -3.420000 0.8642916 0.3865230#> 2 F placebo 12 10.13333 8.075000 -2.058333 0.5247655 0.1514867#> 3 M aspirin 9 11.26667 5.855556 -5.411111 1.1307569 0.3769190#> 4 M placebo 4 11.47500 10.500000 -0.975000 0.7804913 0.3902456如果数据中有NA值并希望忽略它们，设置 na.rm=TRUE
cdata.means <- aggregate(data[c("before","after","change")], by = data[c("sex","condition")], FUN=mean, na.rm=TRUE)cdata.means#> sex condition before after change#> 1 F aspirin 11.06000 7.640000 -3.420000#> 2 M aspirin 11.26667 5.855556 -5.411111#> 3 F placebo 10.13333 8.075000 -2.058333#> 4 M placebo 11.47500 10.500000 -0.975000