R数据分析大数据当中的化整为零（Split-Apply-Combine）策略
R数据分析大数据当中的化整为零（Split-Apply-Combine）策略
引子： 我们常常会遇到这样的问题，数据量很大，并不需要依顺序来依次处理。合理分块处理，并最终整合起来是一个不错的选择。这也就是所谓的Split-Apply-Combine Strategy策略。这在速度上会有比做一个loop有优势，因为它可以并行处理数据。
什么时候我们需要使用到化整为零的策略呢？有以下三种情况：
- 数据需要分组处理
- 数据需要按照每行或者每列来处理
- 数据需要分级处理，和分组很类似，但是分级时需要考虑分级之间的关系。
化整为零策略有点类似于由Google推广的map-reduce策略。当然map-reduce策略的基础是网格，而这里的Split-Apply-Combine的基础完全可以是单机，甚至不支持并行处理的单机都可以。
然而，化整为零并不是一个很直观的编程过程。最直观的过程是使用Loop循环。这里使用一个例子来讲解一下如何实现化整为零策略。在plyr包中有数据ozone,它是一个三维矩阵(24X24X72)，其中最后一维72是指的6年12个月每个月的结果。也就是ozone是一个包括了连续72个月24Ｘ24的三维矩阵数据。三维分别是lat,long,time。我们需要由对时间robust linear model之后的残基residuals。
> library(plyr) # need for dataset ozone> library(MASS) # need for function rlm> month <- ordered(rep(1:12, length=72)) #set time sequence> #try one set> one <- ozone[1,1,]> model <- rlm(one ~ month - 1); model Call:rlm(formula = one ~ month - 1)Converged in 9 iterations Coefficients: month1 month2 month3 month4 month5 month6 month7 month8 month9 month10 month11 month12 264.3964 259.2036 255.0000 252.0052 258.5089 265.3387 274.0000 276.6724 277.0000 285.0000 283.6036 273.1964 Degrees of freedom: 72 total; 60 residual Scale estimate: 4.45 > deseas <- resid(model) |
|---|
现在我们对每一组数据都做相同的处理。首先使用for loop，这样比较直观
> deseasf <- function(value) rlm(value ~ month -1) #the function> models <- as.list(rep(NA, 24*24)) #prepare the variable> dim(models) <- c(24, 24)> deseas <- array(NA, c(24,24,72)) #prepare the variable> dimnames(deseas) <- dimnames(ozone)> for (i in seq_len(24)) { #for loop for first dimension+ for(j in seq_len(24)) { #for loop for second dimension+ mod <- deseasf(ozone[i, j, ]) #apply function+ models[[i, j]] <- mod #save data+ deseas[i, j, ] <- resid(mod) #get residure+ }+ } |
|---|
接下来我们使用化整为零的策略。因为数据可以分成24X24块来处理，每一块都是单独运算，可以并行处理。而使用for loop，只能一块接一块的处理，在速度上可能没有并行处理来得快。而在R当中，有一系列相关的函数，apply, lapply, sapply, tapply, mapply, sweep。我们先得了解这些函数，然后再来应用它们。最简单的是apply。
其形式是apply(array, margin, function, …)。首先,apply的对象是矩阵array或者matrix。它的第二个参数是指的维度，如果你的array是一个二维矩阵，需要按横排的方式计算每一排的平均值，那么你的第二个参数就应该是1。如果需要按纵列的方式计算每一列的平均值，那么第二个参数就应该是2。当然还可以使用c(1,2)这样的方式来设置第二个参数，就是并行计算每个值。第三个参数是需要应用的函数。之后的…是需要传入函数的其它参数。而apply的返回值就是由function来确定的，它可能是vector, array or list。下面举个例子比较容易理解。
> x<-cbind(x1=3,x2=c(4:1,2:5))> dimnames(x)[[1]]<-letters[1:8]> x x1 x2 a 3 4b 3 3c 3 2d 3 1e 3 2f 3 3g 3 4h 3 5> apply(x,2,mean,trim =.2) #在这里，trim =.2就是mean(x, trim = 0, na.rm = FALSE, ...)函数当中的一个参数。x1 x2 3 3 > apply(x,1,mean,trim =.2) a b c d e f g h 3.5 3.0 2.5 2.0 2.5 3.0 3.5 4.0 > col.sums <- apply(x, 2, sum)> row.sums <- apply(x, 1, sum)> rbind(cbind(x, Rtot = row.sums), Ctot = c(col.sums, sum(col.sums))) x1 x2 Rtot a 3 4 7b 3 3 6c 3 2 5d 3 1 4e 3 2 5f 3 3 6g 3 4 7h 3 5 8Ctot 24 24 48> sum.plus.y <- function(x,y){+ sum(x) + y+ }> apply(x, 1, sum.plus.y, 3) #使用自定义函数 a b c d e f g h 10 9 8 7 8 9 10 11 |
|---|
理解了apply,就可比较容易地理解lapply, sapply, vapply了。这三者针对的对象是list或者Vector。其形式为
lapply(X, FUN, ...)sapply(X, FUN, ..., simplify = TRUE, USE.NAMES = TRUE)vapply(X, FUN, FUN.VALUE, ..., USE.NAMES = TRUE) |
|---|
我们看到，它没有了apply当中所需要的第二个参数margin，其原因就是输入对象不是array或者matrix，而是list或者Vector。先举个最简单的应用例子。
> x <- list(a = 1:10, beta = exp(-3:3), logic = c(TRUE,FALSE,FALSE,TRUE))> x$a [1] 1 2 3 4 5 6 7 8 9 10 $beta[1] 0.04978707 0.13533528 0.36787944 1.00000000 2.71828183 7.38905610 20.08553692 $logic[1] TRUE FALSE FALSE TRUE > lapply(x,mean)$a[1] 5.5 $beta[1] 4.535125 $logic[1] 0.5 |
|---|
也就是说，其x的维度应该是一维的，当然，其下面的元素可以是多维的。那么如果x是一个矩阵呢？我们先来看例子。
> x<-cbind(x1=3,x2=c(4:1,2:5))> dimnames(x)[[1]]<-letters[1:8]> x x1 x2 a 3 4b 3 3c 3 2d 3 1e 3 2f 3 3g 3 4h 3 5> x<-as.data.frame(x)> as.list(x)$x1[1] 3 3 3 3 3 3 3 3 $x2[1] 4 3 2 1 2 3 4 5 > lapply(x,function(.ele) mean(.ele))$x1[1] 3 $x2[1] 3 > sapply(x,mean)x1 x2 3 3 > vapply(x,mean,1)x1 x2 3 3 |
|---|
从它们的说明文件我们知道，无论你传入的x是什么，它首先做的一步说是使用as.list来将其转换成一个一维的list。所以，一个data.frame传入lapply之后，它的colnames将会转换成list的names,而rownames可能会丢失。比较可知，lapply和sapply的差别在于，lapply的返回值是一个list，而sapply的返回值是一个矩阵。sapply的返回值其实就是在lapply的基础上再使用了simplify2array(x, higher=TRUE)函数，使用其结果变成一个array。为了更清楚地了解sapply和vapply，我们看下面的例子
> i39 <- sapply(3:9, seq)> i39[[1]][1] 1 2 3 [[2]][1] 1 2 3 4 [[3]][1] 1 2 3 4 5 [[4]][1] 1 2 3 4 5 6 [[5]][1] 1 2 3 4 5 6 7 [[6]][1] 1 2 3 4 5 6 7 8 [[7]][1] 1 2 3 4 5 6 7 8 9 > sapply(i39, fivenum) [,1] [,2] [,3] [,4] [,5] [,6] [,7][1,] 1.0 1.0 1 1.0 1.0 1.0 1[2,] 1.5 1.5 2 2.0 2.5 2.5 3[3,] 2.0 2.5 3 3.5 4.0 4.5 5[4,] 2.5 3.5 4 5.0 5.5 6.5 7[5,] 3.0 4.0 5 6.0 7.0 8.0 9> vapply(i39, fivenum,+ c(Min. = 0, "1st Qu." = 0, Median = 0, "3rd Qu." = 0, Max. = 0)) [,1] [,2] [,3] [,4] [,5] [,6] [,7]Min. 1.0 1.0 1 1.0 1.0 1.0 11st Qu. 1.5 1.5 2 2.0 2.5 2.5 3Median 2.0 2.5 3 3.5 4.0 4.5 53rd Qu. 2.5 3.5 4 5.0 5.5 6.5 7Max. 3.0 4.0 5 6.0 7.0 8.0 9 |
|---|
其中，fivenum函数会返回一组数的最小值，四分位低值(lower-hinge)，中值，四分位高值(upper-hinge)，以及最大值。从上面的比较中，我们很清楚的看到，sapply返回值的排列形式，以list的names为colnames。可以想象，它使用的是按列填充matrix的方式输出的。而vapply是在sapply的基础上，为rownames做出了定义。
除了上面介绍的，还有tapply,mapply,sweep等。它们的定义如下。如果需要了解和掌握它们，需要熟悉上面介绍的apply, lapply, sapply以及vapply。还需要了解split。所以这里就不多加解释了（因为篇幅会很长）。
tapply(X, INDEX, FUN = NULL, ..., simplify = TRUE)mapply(FUN, ..., MoreArgs = NULL, SIMPLIFY = TRUE, USE.NAMES = TRUE)sweep(x, MARGIN, STATS, FUN="-", check.margin=TRUE, ...) |
|---|
有了上面的了解，我们可以看一下最前面的for loop可以怎样改写为apply形式了。
> models <- apply(ozone, 1:2, deseasf) #这里相当于for loop当中的for(i in seq_len(24)){for(j in seq_len(24)){mod<-deseasf(ozone[i,j,]); models[[i,j]]<-mod;}}, 但是运算却是并行处理的。> resids_list <- lapply(models, resid)> resids <- unlist(resids_list)> dim(resids) <- c(72, 24, 24)> deseas <- aperm(resids, c(2, 3, 1))> dimnames(deseas) <- dimnames(ozone) |
|---|
如果，我们知道的并不是ozone这样一个24X24X72的已知确切维度的三维数组，假如我们只有一个名为ozonedf的matrix，它的5列分别为lat, long, time, month,和value。我们如果需要做上述的分析应该怎么办呢？在思路上，我们的想法可能会是先从ozonedf出发生成一个类似ozone这样子的数据，然后再使用apply,lapply这样的函数来完成就可以。第一步生成ozone这样子的数据，就是化整为零策略(Split-Apply-Combine)的第一步了。
R> deseasf_df <- function(df) {+ rlm(value ~ month - 1, data = df)+ }R> pieces <- split(ozonedf, list(ozonedf$lat, ozonedf$long))R> models <- lapply(pieces, deseasf_df)R> results <- mapply(function(model, df) {+ cbind(df[rep(1, 72), c("lat", "long")], resid(model))+ }, models, pieces)R> deseasdf <- do.call("rbind", results) |
|---|
我们再来举一个airquality的例子，我们先来构造一个比较复杂的数据。
> aq<-airquality> aq$source<-"qiuworld"> head(aq) Ozone Solar.R Wind Temp Month Day source1 41 190 7.4 67 5 1 qiuworld2 36 118 8.0 72 5 2 qiuworld3 12 149 12.6 74 5 3 qiuworld4 18 313 11.5 62 5 4 qiuworld5 NA NA 14.3 56 5 5 qiuworld6 28 NA 14.9 66 5 6 qiuworld> aq1<-aq> aq1$source<-c("gaziou")> head(aq1) Ozone Solar.R Wind Temp Month Day source1 41 190 7.4 67 5 1 gaziou2 36 118 8.0 72 5 2 gaziou3 12 149 12.6 74 5 3 gaziou4 18 313 11.5 62 5 4 gaziou5 NA NA 14.3 56 5 5 gaziou6 28 NA 14.9 66 5 6 gaziou> set.seed(123)> x1<-runif(nrow(aq1),-0.5,0.5)> head(x1)[1] -0.21242248 0.28830514 -0.09102308 0.38301740 0.44046728 -0.45444350> aq1$Temp<-aq1$Temp+x1> head(aq1$Temp)[1] 66.78758 72.28831 73.90898 62.38302 56.44047 65.54556> aq<-rbind(aq,aq1)> aq<-aq[order(aq$Wind),]> head(aq) Ozone Solar.R Wind Temp Month Day source53 NA 59 1.7 76.00000 6 22 qiuworld206 NA 59 1.7 76.29892 6 22 gaziou121 118 225 2.3 94.00000 8 29 qiuworld274 118 225 2.3 94.14789 8 29 gaziou126 73 183 2.8 93.00000 9 3 qiuworld279 73 183 2.8 93.48422 9 3 gaziou |
|---|
我们的任务就是比较不同source来源的每个Month的平均Temp。思路应该是先把数据按照source和Month分成小块，计算出来其Temp的平均值，然后输出。这就是一个完整而简单的Split-Apply-Combine的过程了。
> pieces <- split(aq, list(aq$source,aq$Month))> pieces[1:2]$gaziou.5 Ozone Solar.R Wind Temp Month Day source183 115 223 5.7 78.64711 5 30 gaziou164 7 NA 6.9 74.45683 5 11 gaziou154 41 190 7.4 66.78758 5 1 gaziou184 37 279 7.4 76.46302 5 31 gaziou155 36 118 8.0 72.28831 5 2 gaziou180 NA NA 8.0 57.04407 5 27 gaziou160 23 299 8.6 65.02811 5 7 gaziou163 NA 194 8.6 68.95661 5 10 gaziou166 11 290 9.2 66.17757 5 13 gaziou165 16 256 9.7 68.95333 5 12 gaziou173 11 44 9.7 62.45450 5 20 gaziou174 1 8 9.7 59.38954 5 21 gaziou176 4 25 9.7 61.14051 5 23 gaziou167 14 274 10.9 68.07263 5 14 gaziou157 18 313 11.5 62.38302 5 4 gaziou169 14 334 11.5 64.39982 5 16 gaziou172 30 322 11.5 67.82792 5 19 gaziou170 34 307 12.0 65.74609 5 17 gaziou177 32 92 12.0 61.49427 5 24 gaziou181 23 13 12.0 67.09414 5 28 gaziou156 12 149 12.6 73.90898 5 3 gaziou168 18 65 13.2 57.60292 5 15 gaziou161 19 99 13.8 59.39242 5 8 gaziou158 NA NA 14.3 56.44047 5 5 gaziou159 28 NA 14.9 65.54556 5 6 gaziou179 NA 266 14.9 58.20853 5 26 gaziou182 45 252 14.9 80.78916 5 29 gaziou175 11 320 16.6 73.19280 5 22 gaziou178 NA 66 16.6 57.15571 5 25 gaziou171 6 78 18.4 56.54206 5 18 gaziou162 8 19 20.1 61.05144 5 9 gaziou $qiuworld.5 Ozone Solar.R Wind Temp Month Day source30 115 223 5.7 79 5 30 qiuworld11 7 NA 6.9 74 5 11 qiuworld1 41 190 7.4 67 5 1 qiuworld31 37 279 7.4 76 5 31 qiuworld2 36 118 8.0 72 5 2 qiuworld27 NA NA 8.0 57 5 27 qiuworld7 23 299 8.6 65 5 7 qiuworld10 NA 194 8.6 69 5 10 qiuworld13 11 290 9.2 66 5 13 qiuworld12 16 256 9.7 69 5 12 qiuworld20 11 44 9.7 62 5 20 qiuworld21 1 8 9.7 59 5 21 qiuworld23 4 25 9.7 61 5 23 qiuworld14 14 274 10.9 68 5 14 qiuworld4 18 313 11.5 62 5 4 qiuworld16 14 334 11.5 64 5 16 qiuworld19 30 322 11.5 68 5 19 qiuworld17 34 307 12.0 66 5 17 qiuworld24 32 92 12.0 61 5 24 qiuworld28 23 13 12.0 67 5 28 qiuworld3 12 149 12.6 74 5 3 qiuworld15 18 65 13.2 58 5 15 qiuworld8 19 99 13.8 59 5 8 qiuworld5 NA NA 14.3 56 5 5 qiuworld6 28 NA 14.9 66 5 6 qiuworld26 NA 266 14.9 58 5 26 qiuworld29 45 252 14.9 81 5 29 qiuworld22 11 320 16.6 73 5 22 qiuworld25 NA 66 16.6 57 5 25 qiuworld18 6 78 18.4 57 5 18 qiuworld9 8 19 20.1 61 5 9 qiuworld > avgTemp<-lapply(pieces,function(.ele) mean(.ele$Temp))> head(avgTemp)$gaziou.5[1] 65.63339 $qiuworld.5[1] 65.54839 $gaziou.6[1] 79.02871 $qiuworld.6[1] 79.1 $gaziou.7[1] 83.90313 $qiuworld.7[1] 83.90323 > avgTemp<-do.call(rbind,avgTemp)> head(avgTemp) [,1]gaziou.5 65.63339qiuworld.5 65.54839gaziou.6 79.02871qiuworld.6 79.10000gaziou.7 83.90313qiuworld.7 83.90323 |
|---|
由上面的过程我们可以看来，其实化整为零策略在实现起来，就是分三步走，使用split将数据化分成小块，使用lapply函数对小块进行计算，最后使用do.call使用函数将其整理成我们需要的形式。这个过程，其实使用plyr包来实现，就更为便洁了。同样是上面的操作，使用plyr的话，只需要一行即可。
> avgTemp1<-ddply(aq,.(source,Month),function(.ele) mean(.ele$Temp))> avgTemp1 source Month V11 gaziou 5 65.633392 gaziou 6 79.028713 gaziou 7 83.903134 gaziou 8 83.986115 gaziou 9 76.893196 qiuworld 5 65.548397 qiuworld 6 79.100008 qiuworld 7 83.903239 qiuworld 8 83.9677410 qiuworld 9 76.90000 |
|---|
plyr包给我们提供了非常简洁的书写方式。我们来看看其主要函数的定义方式。
Inputoutput | Array | Data frame | List | Discarded |
|---|---|---|---|---|
Array | aaply | adply | alply | a_ply |
Data frame | daply | ddply | dlply | d_ply |
List | laply | ldply | llply | l_ply |
a*ply(.data, .margins, .fun, ..., .progress = "none")d*ply(.data, .variables, .fun, ..., .progress = "none")l*ply(.data, .fun, ..., .progress = "none") |
|---|
我们可以看出，函数名称的第一个字母代表输入的形式，它们分别是a->Array, d->Data frame, l->List。而第二字母代表输出的形式，它们的定义同前。对于输入为array和data frame的，函数的第二个参数为data的margins或者variables。对于margins,可以是
- .margins = 1 #以行为单位
- .margins = 2 #以列为单位
- .margins = c(1,2) #以individual cell为单位
需要注意的是，这里的每一个参数都使用了.点号起始。在作者看来，这样可以避免误操作。
腾讯云开发者
