R语言高级数据管理

Published on 2017 - 04 - 18

一个数据处理难题

要讨论数值和字符处理函数,让我们首先考虑一个数据处理问题。一组学生参加了数学、科学和英语考试。为了给所有学生确定一个单一的成绩衡量指标,需要将这些科目的成绩组合起来。另外,你还想将前20%的学生评定为A,接下来20%的学生评定为B,依次类推。最后,你希望按字母顺序对学生排序。数据如表1所示。

学生姓名 数  学 科  学 英  语
John Davis 502 95 25
Angela Williams 600 99 22
Bullwinkle Moose 412 80 18
David Jones 358 82 15
Janice Markhammer 495 75 20
Cheryl Cushing 512 85 28
Reuven Ytzrhak 410 80 15
Greg Knox 625 95 30
Joel England 573 89 27
Mary Rayburn 522 86 18

观察此数据集,马上可以发现一些明显的障碍。首先,三科考试的成绩是无法比较的。由于它们的均值和标准差相去甚远,所以对它们求平均值是没有意义的。你在组合这些考试成绩之前,必须将其变换为可比较的单元。其次,为了评定等级,你需要一种方法来确定某个学生在前述得分上百分比排名。再次,表示姓名的字段只有一个,这让排序任务复杂化了。为了正确地将其排序,需要将姓和名拆开。

以上每一个任务都可以巧妙地利用R中的数值和字符处理函数完成。在讲解完下一节中的各种函数之后,我们将考虑一套可行的解决方案,以解决这项数据处理难题。

数值和字符处理函数

本节我们将综述R中作为数据处理基石的函数,它们可分为数值(数学、统计、概率)函数和字符处理函数。在阐述过每一类函数以后,我将为你展示如何将函数应用到矩阵和数据框的列(变量)和行(观测)上。

数学函数

表2列出了常用的数学函数和简短的用例。

函  数 描  述
abs(x) 绝对值
abs(-4) 返回值为4
sqrt(x) 平方根sqrt(25)返回值为5和25^(0.5)等价
ceiling(x) 不小于x的最小整数 ceiling(3.475) 返回值为4
floor(x) 不大于x的最大整数floor(3.475)返回值为3
trunc(x) 向 0 的方向截取的x中的整数部分
trunc(5.99) 返回值为5
round(x, digits=n) 将x舍入为指定位的小数 round(3.475, digits=2)返回值为3.48
signif(x, digits=n) 将x舍入为指定的有效数字位数 signif(3.475, digits=2) 返回值为3.5
cos(x)、sin(x) 、tan(x) 余弦、正弦和正切 cos(2)返回值为–0.416
acos(x) 、asin(x) 、atan(x) 反余弦、反正弦和反正切 acos(-0.416)返回值为2
cosh(x) 、sinh(x) 、tanh(x) 双曲余弦、双曲正弦和双曲正切 sinh(2)返回值为3.627
acosh(x) 、asinh(x) 、atanh(x) 反双曲余弦、反双曲正弦和反双曲正切 asinh(3.627)返回值为2
log(x,base=n)、log(x)、log10(x) 对x取以n为底的对数,
为了方便起见log(x)为自然对数,log10(x)为常用对数,
log(10)返回值为2.3026,log10(10)返回值为1
exp(x) 指数函数 exp(2.3026)返回值为10

对数据做变换是这些函数的一个主要用途。例如,你经常会在进一步分析之前将收入这种存在明显偏倚的变量取对数。数学函数也被用作公式中的一部分,用于绘图函数(例如x 对 sin(x))和在输出结果之前对数值做格式化。

表2中的示例将数学函数应用到了标量(单独的数值)上。当这些函数被应用于数值向量、矩阵或数据框时,它们会作用于每一个独立的值。例如,sqrt(c(4, 16, 25))的返回值为c(2, 4, 5)。

统计函数

常用的统计函数如表3所示,其中许多函数都拥有可以影响输出结果的可选参数。

函  数 描  述
mean(x) 平均数,mean(c(1,2,3,4)) 返回值为2.5
median(x) 中位数,median(c(1,2,3,4))返回值为2.5
sd(x) 标准差,sd(c(1,2,3,4))返回值为1.29
var(x) 方差,var(c(1,2,3,4))返回值为1.67
mad(x) 绝对中位差(median absolute deviation)
mad(c(1,2,3,4))返回值为1.48
quantile(x,probs) 求分位数。其中x为待求分位数的数值型向量,
probs为一个由[0,1]之间的概率值组成的数值向量
# 求x的30%和84%分位点
y <- quantile(x, c(.3,.84))
range(x) 求值域,x <- c(1,2,3,4),
range(x)返回值为c(1,4)
diff(range(x))返回值为3
sum(x) 求和,sum(c(1,2,3,4))返回值为10
diff(x, lag=n) 滞后差分,lag用以指定滞后几项。默认的lag值为1,
x<- c(1, 5, 23, 29)
diff(x)返回值为c(4, 18, 6)
min(x) 求最小值,min(c(1,2,3,4))返回值为1
max(x) 求最大值,max(c(1,2,3,4))返回值为4
scale(x,center=TRUE, scale=TRUE) 为数据对象x按列进行中心化(center=TRUE)或标准化(center=TRUE,scale=TRUE);

举例来说:

y <- mean(x)

提供了对象x中元素的算术平均数,而:

z <- mean(x, trim = 0.05, na.rm=TRUE)

则提供了截尾平均数,即丢弃了最大5%和最小5%的数据和所有缺失值后的算术平均数。请使用help()了解以上每个函数和其参数的用法。

以下这段代码演示了计算某个数值向量的均值和标准差的两种方式。

> x <- c(1,2,3,4,5,6,7,8)

> mean(x)       # 简洁的方式
[1] 4.5
> sd(x)
[1] 2.449490

> n <- length(x)     # 冗长的方式
> meanx <- sum(x)/n
> css <- sum((x - meanx)^2)
> sdx <- sqrt(css / (n-1))
> meanx
[1] 4.5
> sdx
[1] 2.449490

数据的标准化

默认情况下,函数scale()对矩阵或数据框的指定列进行均值为0、标准差为1的标准化:

newdata <- scale(mydata)

要对每一列进行任意均值和标准差的标准化,可以使用如下的代码:

newdata <- scale(mydata)*SD + M

其中的M是想要的均值,SD为想要的标准差。在非数值型的列上使用scale()函数将会报错。要对指定列而不是整个矩阵或数据框进行标准化,你可以使用这样的代码:

newdata <- transform(mydata, myvar = scale(myvar)*10+50)

此句将变量myvar标准化为均值50、标准差为10的变量。

数据的中心化
所谓数据的中心化是指数据集中的各项数据减去数据集的均值。例如有数据集1, 2, 3, 6, 3,其均值为3,那么中心化之后的数据集为1-3,2-3,3-3,6-3,3-3,即:-2,-1,0,3,0
数据的标准化
所谓数据的标准化是指中心化之后的数据在除以数据集的标准差,即数据集中的各项数据减去数据集的均值再除以数据集的标准差。
例如有数据集1, 2, 3, 6, 3,其均值为3,其标准差为1.87,那么标准化之后的数据集为(1-3)/1.87,(2-3)/1.87,(3-3)/1.87,(6-3)/1.87,(3-3)/1.87,即:-1.069,-0.535,0,1.604,0
数据中心化和标准化的意义是一样的,为了消除量纲对数据结构的影响。

概率函数

你可能在疑惑为何概率函数未和统计函数列在一起。(你真的对此有些困惑,对吧?)虽然根据定义,概率函数也属于统计类,但是它们非常独特,应独立设一节进行讲解。概率函数通常用来生成特征已知的模拟数据,以及在用户编写的统计函数中计算概率值。

在R中,概率函数形如:

[dpqr]distribution_abbreviation()

其中第一个字母表示其所指分布的某一方面:

  • d = 密度函数(density)
  • p = 分布函数(distribution function)
  • q = 分位数函数(quantile function)
  • r = 生成随机数(随机偏差)

常用的概率函数列于表4中。

分布名称 缩  写
Beta分布 beta
二项分布 binom
柯西分布 cauchy
(非中心)卡方分布 chisq
指数分布 exp
F分布 f
Gamma分布 gamma
几何分布 geom
超几何分布 hyper
对数正态分布 lnorm
Logistic分布 logis
多项分布 multinom
负二项分布 nbinom
正态分布 norm
泊松分布 pois
Wilcoxon符号秩分布 signrank
t分布 t
均匀分布 unif
Weibull分布 weibull
Wilcoxon秩和分布 wilcox

我们不妨先看看正态分布的有关函数,以了解这些函数的使用方法。如果不指定一个均值和一个标准差,则函数将假定其为标准正态分布(均值为0,标准差为1)。密度函数(dnorm)、分布函数(pnorm)、分位数函数(qnorm)和随机数生成函数(rnorm)的使用示例见如下。

在区间[-3, 3]上绘制标准正态曲线

x <- pretty(c(-3,3), 30) 
y <-dnorm(x) 
plot(x, y, 
  type="l",
  xlab="NormalDeviate",
  ylab="Density",
  yaxs="i"
)

位于 z=1.96 左侧的标准正态曲线下方面积是多少?
pnorm(1.96)等于0.975

均值为500,标准差为100的正态分布的0.9分位点值为多少?

qnorm(.9, mean=500, sd=100)等于628.16

生成50个均值为50,标准差为10的正态随机数
rnorm(50, mean=50, sd=10)

设定随机数种子

在每次生成伪随机数的时候,函数都会使用一个不同的种子,因此也会产生不同的结果。你可以通过函数set.seed()显式指定这个种子,让结果可以重现(reproducible)。代码清单2给出了一个示例。这里的函数runif()用来生成0到1区间上服从均匀分布的伪随机数。

> runif(5)
[1] 0.8725344 0.3962501 0.6826534 0.3667821 0.9255909
> runif(5)
[1] 0.4273903 0.2641101 0.3550058 0.3233044 0.6584988
> set.seed(1234)
> runif(5)
[1] 0.1137034 0.6222994 0.6092747 0.6233794 0.8609154
> set.seed(1234)
> runif(5)
[1] 0.1137034 0.6222994 0.6092747 0.6233794 0.8609154

通过手动设定种子,就可以重现你的结果了。这种能力有助于我们创建会在未来取用的,以及可与他人分享的示例。

生成多元正态数据

在模拟研究和蒙特卡洛方法中,你经常需要获取来自给定均值向量和协方差阵的多元正态分布的数据。MASS包中的mvrnorm()函数可以让这个问题变得很容易。其调用格式为:

mvrnorm(n, mean, sigma)

其中n是你想要的样本大小,mean为均值向量,而sigma是方差-协方差矩阵(或相关矩阵)。在代码清单3中,你将从一个参数如下所示的三元正态分布中抽取500个观测。

均值向量 230.7 146.7 3.6
协方差阵 15360.8
6721.2
-47.1
6721.2
4700.9
-16.5
-47.1
-16.5
0.3
> library(MASS)
> options(digits=3)
> set.seed(1234)             # ❶ 设定随机数种子

> mean <- c(230.7, 146.7, 3.6)             # ❷ 指定均值向量、协方差阵
> sigma <- matrix(c(15360.8, 6721.2, -47.1,
                     6721.2, 4700.9, -16.5,
                      -47.1,  -16.5,   0.3), nrow=3, ncol=3)

> mydata <- mvrnorm(500, mean, sigma)         # ❸ 生成数据
> mydata <- as.data.frame(mydata)
> names(mydata) <- c("y","x1","x2")

> dim(mydata)          # ❹ 查看结果
[1] 500 3
> head(mydata, n=10)
       y    x1   x2
1   98.8  41.3 4.35
2  244.5 205.2 3.57
3  375.7 186.7 3.69
....

代码清单3中设定了一个随机数种子,这样就可以在之后重现结果❶。你指定了想要的均值向量和方差—协方差阵❷,并生成了500个伪随机观测❸。为了方便,结果从矩阵转换为数据框,并为变量指定了名称。最后,你确认了拥有500个观测和3个变量,并输出了前10个观测❹。请注意,由于相关矩阵同时也是协方差阵,所以其实可以直接指定相关关系的结构。

R中的概率函数允许生成模拟数据,这些数据是从服从已知特征的概率分布中抽样而得的。

字符处理函数

数学和统计函数是用来处理数值型数据的,而字符处理函数可以从文本型数据中抽取信息,或者为打印输出和生成报告重设文本的格式。举例来说,你可能希望将某人的姓和名连接在一起,并保证姓和名的首字母大写,抑或想统计可自由回答的调查反馈信息中含有秽语的实例(instance)数量。一些最有用的字符处理函数见表6。

函  数 描  述
nchar(x) 计算x中的字符数量
x <- c("ab", "cde","fghij")
length(x)返回值为 3 (参见表7)
nchar(x[3])返回值为5
substr(x, start, stop) 提取或替换一个字符向量中的子串
x <- "abcdef"
substr(x, 2, 4)返回值为"bcd"
substr(x, 2, 4) <- "22222"(x将变成"a222ef")
grep(pattern, x, ignore. case=FALSE, fixed=FALSE) 在x中搜索某种模式。若fixed=FALSE,
则pattern为一个正则表达式。
若fixed=TRUE,则pattern为一个文本字符串。
返回值为匹配的下标
grep("A",c("b","A","c"),fixed=TRUE) 返回值为2
sub(pattern, replacement, x, ignore.case=FALSE, fixed=FALSE) 在x中搜索pattern,并以文本replacement将其替换。若fixed=FALSE,则pattern为一个正则表达式。若fixed=TRUE,则pattern为一个文本字符串
sub("\s",".","Hello There")返回值为Hello.There。
注意,"\s"是一个用来查找空白的正则表达式;使用"\s"而不用"\"的原因是,后者是R中的转义字符
strsplit(x, split, fixed=FALSE) 在split处分割字符向量x中的元素。若fixed=FALSE,则pattern为一个正则表达式。若fixed=TRUE,则pattern为一个文本字符串
y <- strsplit("abc", "")
将返回一个含有1个成分、3个元素的列表,包含的内容为"a" "b" "c"
unlist(y)[2]和sapply(y, "[", 2)均会返回"b"
paste(…, sep="") 连接字符串,分隔符为sep
paste("x", 1:3,sep="")返回值为c("x1", "x2", "x3")
paste("x",1:3,sep="M")返回值为c("xM1","xM2" "xM3")
paste("Today is", date())返回值为Today is Thu Jun 25 14:17:32 2011
toupper(x) 大写转换
toupper("abc") 返回值为"ABC"
tolower(x) 小写转换
tolower("ABC")返回值为"abc"

请注意,函数grep()、sub()和strsplit()能够搜索某个文本字符串(fixed=TRUE)或某个正则表达式(fixed=FALSE,默认值为FALSE)。正则表达式为文本模式的匹配提供了一套清晰而简练的语法。例如,正则表达式:

^[hc]?at

可匹配任意以0个或1个h或c开头、后接at的字符串。因此,此表达式可以匹配hat、cat和at,但不会匹配bat。要了解更多,请参考维基百科的regular expression(正则表达式)条目。

其他实用函数

表7中的函数对于数据管理和处理同样非常实用,只是它们无法清楚地划入其他分类中。

函  数 描  述
length(x) 对象x的长度
x <- c(2, 5, 6, 9)
length(x)返回值为4
seq(from, to, by) 生成一个序列
indices <- seq(1,10,2)
indices的值为c(1, 3, 5, 7, 9)
rep(x, n) 将x重复n次
y <- rep(1:3, 2)
y的值为c(1, 2, 3, 1, 2, 3)
cut(x, n) 将连续型变量x分割为有着n个水平的因子
使用选项ordered_result = TRUE以创建一个有序型因子
pretty(x, n) 创建美观的分割点。通过选取n+1个等间距的取整值,将一个连续型变量x分割为n个区间。绘图中常用
cat(... , file ="myfile", append =FALSE) 连接...中的对象,并将其输出到屏幕上或文件中
firstname <- c("Jane")
cat("Hello" ,firstname, "\n")

表中的最后一个例子演示了在输出时转义字符的使用方法。\n表示新行,\t为制表符,\'为单引号,\b为退格,等等。(键入?Quotes以了解更多。)例如,代码:

name <- "Bob"
cat( "Hello", name, "\b.\n", "Isn\'t R", "\t", "GREAT?\n")
可生成:
Hello Bob.
 Isn’t R         GREAT?

请注意第二行缩进了一个空格。当cat输出连接后的对象时,它会将每一个对象都用空格分开。这就是在句号之前使用退格转义字符(\b)的原因。不然,生成的结果将是“Hello Bob .”。

在数值、字符串和向量上使用我们最近学习的函数是直观而明确的,但是如何将它们应用到矩阵和数据框上呢?这就是下一节的主题。

将函数应用于矩阵和数据框

R函数的诸多有趣特性之一,就是它们可以应用到一系列的数据对象上,包括标量、向量、矩阵、数组和数据框。代码清单4提供了一个示例。

> a <- 5
> sqrt(a)
[1] 2.236068
> b <- c(1.243, 5.654, 2.99)
> round(b)
[1] 1 6 3
> c <- matrix(runif(12), nrow=3)
> c
       [,1]  [,2]  [,3]  [,4]
[1,] 0.4205 0.355 0.699 0.323
[2,] 0.0270 0.601 0.181 0.926
[3,] 0.6682 0.319 0.599 0.215
> log(c)
       [,1]   [,2]   [,3]   [,4]
[1,] -0.866 -1.036 -0.358 -1.130
[2,] -3.614 -0.508 -1.711 -0.077
[3,] -0.403 -1.144 -0.513 -1.538
> mean(c)
[1] 0.444

请注意,在代码清单4中对矩阵c求均值的结果为一个标量(0.444)。函数mean()求得的是矩阵中全部12个元素的均值。但如果希望求的是各行的均值或各列的均值呢?

R中提供了一个apply()函数,可将一个任意函数“应用”到矩阵、数组、数据框的任何维度上。apply函数的使用格式为:

apply(x, MARGIN, FUN, ...)

其中,x为数据对象,MARGIN是维度的下标,FUN是由你指定的函数,而...则包括了任何想传递给FUN的参数。在矩阵或数据框中,MARGIN=1表示行,MARGIN=2表示列。请看代码清单5中的例子。

> mydata <- matrix(rnorm(30), nrow=6)        # ❶ 生成数据
> mydata
         [,1]   [,2]    [,3]   [,4]   [,5]
[1,]  0.71298  1.368 -0.8320 -1.234 -0.790
[2,] -0.15096 -1.149 -1.0001 -0.725  0.506
[3,] -1.77770  0.519 -0.6675  0.721 -1.350
[4,] -0.00132 -0.308  0.9117 -1.391  1.558
[5,] -0.00543  0.378 -0.0906 -1.485 -0.350
[6,] -0.52178 -0.539 -1.7347  2.050  1.569
> apply(mydata, 1, mean)             # ❷ 计算每行的均值
[1] -0.155 -0.504 -0.511  0.154 -0.310  0.165
> apply(mydata, 2, mean)           # ❸ 计算每列的均值
[1] -0.2907  0.0449 -0.5688 -0.3442  0.1906
> apply(mydata, 2, mean, trim=0.2)          # ❹ 计算每列的截尾均值
[1] -0.1699  0.0127 -0.6475 -0.6575  0.2312

首先生成了一个包含正态随机数的6×5矩阵❶。然后你计算了6行的均值❷,以及5列的均值❸。最后,你计算了每列的截尾均值(在本例中,截尾均值基于中间60%的数据,最高和最低20%的值均被忽略)❹。

FUN可为任意R函数,这也包括你自行编写的函数,所以apply()是一种很强大的机制。apply()可把函数应用到数组的某个维度上,而lapply()和sapply()则可将函数应用到列表(list)上。

数据处理难题的一套解决方案

第一节中提出的问题是:将学生的各科考试成绩组合为单一的成绩衡量指标、基于相对名次(前20%,下20%,等等)给出从A到F的评分、根据学生姓氏和名字的首字母对花名册进行排序。代码清单6给出了一种解决方案。

> options(digits=2)

> Student <- c("John Davis", "Angela Williams", "Bullwinkle Moose",
               "David Jones", "Janice Markhammer", "Cheryl Cushing",
               "Reuven Ytzrhak", "Greg Knox", "Joel England",
               "Mary Rayburn")
> Math <- c(502, 600, 412, 358, 495, 512, 410, 625, 573, 522)
> Science <- c(95, 99, 80, 82, 75, 85, 80, 95, 89, 86)
> English <- c(25, 22, 18, 15, 20, 28, 15, 30, 27, 18)
> roster <- data.frame(Student, Math, Science, English,
                       stringsAsFactors=FALSE)

> z <- scale(roster[,2:4])             # 计算综合得分
> score <- apply(z, 1, mean)
> roster <- cbind(roster, score)
> y <- quantile(score, c(.8,.6,.4,.2))         # 对学生评分
> roster$grade[score >= y[1]] <- "A"
> roster$grade[score < y[1] & score >= y[2]] <- "B"
> roster$grade[score < y[2] & score >= y[3]] <- "C"
> roster$grade[score < y[3] & score >= y[4]] <- "D"
> roster$grade[score < y[4]] <- "F"

> name <- strsplit((roster$Student), " ")           # 抽取姓氏和名字
> lastname <- sapply(name, "[", 2)
> firstname <- sapply(name, "[", 1)
> roster <- cbind(firstname,lastname, roster[,-1])

> roster <- roster[order(lastname,firstname),]          # 根据姓氏和名字排序

> roster
    Firstname   Lastname Math Science English score grade
6      Cheryl    Cushing  512      85      28  0.35     C
1        John      Davis  502      95      25  0.56     B
9        Joel    England  573      89      27  0.70     B
.......

以上代码写得比较紧凑,逐步分解如下。

步骤1 原始的学生花名册已经给出了。options(digits=2)限定了输出小数点后数字的位数,并且让输出更容易阅读。

> options(digits=2)
> roster
             Student Math Science English
1         John Davis  502      95      25
2    Angela Williams  600      99      22
3   Bullwinkle Moose  412      80      18
......

步骤2 由于数学、科学和英语考试的分值不同(均值和标准差相去甚远),在组合之前需 要先让它们变得可以比较。一种方法是将变量进行标准化,这样每科考试的成绩就都是用单位标准差来表示,而不是以原始的尺度来表示了。这个过程可以使用scale()函数来实现:

> z <- scale(roster[,2:4])
> z
        Math Science English
 [1,]  0.013   1.078   0.587
 [2,]  1.143   1.591   0.037
 [3,] -1.026  -0.847  -0.697
 ......

步骤3 然后,可以通过函数mean()来计算各行的均值以获得综合得分,并使用函数 cbind()将其添加到花名册中:

> score <- apply(z, 1, mean)
> roster <- cbind(roster, score)
> roster
            Student   Math  Science English  score
1         John Davis  502      95      25    0.559
2    Angela Williams  600      99      22    0.924
3   Bullwinkle Moose  412      80      18   -0.857
....

步骤4 函数quantile()给出了学生综合得分的百分位数。可以看到,成绩为A的分界点 为0.74,B的分界点为0.44,等等。

> y <- quantile(roster$score, c(.8,.6,.4,.2))
> y
  80%   60%   40%   20%
 0.74  0.44 -0.36 -0.89

步骤5 通过使用逻辑运算符,你可以将学生的百分位数排名重编码为一个新的类别型成绩变量。下面在数据框roster中创建了变量grade。

> roster$grade[score >= y[1]] <- "A"
> roster$grade[score < y[1] & score >= y[2]] <- "B"
> roster$grade[score < y[2] & score >= y[3]] <- "C"
> roster$grade[score < y[3] & score >= y[4]] <- "D"
> roster$grade[score < y[4]] <- "F"
> roster
             Student Math Science English  score grade
1         John Davis  502      95      25  0.559     B
2    Angela Williams  600      99      22  0.924     A
3   Bullwinkle Moose  412      80      18 -0.857     D
.......

步骤6 你将使用函数strsplit()以空格为界把学生姓名拆分为姓氏和名字。把strsplit()应用到一个字符串组成的向量上会返回一个列表:

> name <- strsplit((roster$Student), " ")
> name

[[1]]
[1] "John"  "Davis"

[[2]]
[1] "Angela"   "Williams"

[[3]]
[1] "Bullwinkle" "Moose"

.......

步骤7 你可以使用函数sapply()提取列表中每个成分的第一个元素,放入一个储存名字 的向量,并提取每个成分的第二个元素,放入一个储存姓氏的向量。"["是一个可以提取某个对象的一部分的函数——在这里它是用来提取列表name各成分中的第一个或第二个元素的。你将使用cbind()把它们添加到花名册中。由于已经不再需要student变量,可以将其丢弃(在下标中使用?)。

> Firstname <- sapply(name, "[", 1)
> Lastname <- sapply(name, "[", 2)
> roster <- cbind(Firstname, Lastname, roster[,-1])
> roster
    Firstname   Lastname Math Science English  score grade
1        John      Davis  502      95      25  0.559     B
2      Angela   Williams  600      99      22  0.924     A
3  Bullwinkle      Moose  412      80      18 -0.857     D
.....

步骤8 最后,可以使用函数order()依姓氏和名字对数据集进行排序:

> roster[order(Lastname,Firstname),]
    Firstname   Lastname  Math  Science  English  score   grade
6      Cheryl    Cushing  512      85       28     0.35     C
1        John      Davis  502      95       25     0.56     B
9        Joel    England  573      89       27     0.70     B
......

控制流

在正常情况下,R程序中的语句是从上至下顺序执行的。但有时你可能希望重复执行某些语句,仅在满足特定条件的情况下执行另外的语句。这就是控制流结构发挥作用的地方了。

R拥有一般现代编程语言中都有的标准控制结构。首先你将看到用于条件执行的结构,接下来是用于循环执行的结构。

为了理解贯穿本节的语法示例,请牢记以下概念:

  • 语句(statement)是一条单独的R语句或一组复合语句(包含在花括号{ } 中的一组R语句,使用分号分隔);
  • 条件(cond)是一条最终被解析为真(TRUE)或假(FALSE)的表达式;
  • 表达式(expr)是一条数值或字符串的求值语句;
  • 序列(seq)是一个数值或字符串序列。

在讨论过控制流的构造后,我们将学习如何编写函数。

重复和循环

循环结构重复地执行一个或一系列语句,直到某个条件不为真为止。循环结构包括for和while结构。

for结构

for循环重复地执行一个语句,直到某个变量的值不再包含在序列seq中为止。语法为:

for (var in seq) statement

在下例中:

for (i in 1:10) print("Hello")

单词Hello被输出了10次。

while结构

while循环重复地执行一个语句,直到条件不为真为止。语法为:

while (cond) statement

作为第二个例子,代码:

i <- 10
while (i > 0) {print("Hello"); i <- i - 1}

又将单词Hello输出了10次。请确保括号内while的条件语句能够改变,即让它在某个时刻不再为真——否则循环将永不停止!在上例中,语句:

i <- i - 1

在每步循环中为对象i减去1,这样在十次循环过后,它就不再大于0了。反之,如果在每步循环都加1的话,R将不停地输出Hello。这也是while循环可能较其他循环结构更危险的原因。

在处理大数据集中的行和列时,R中的循环可能比较低效费时。只要可能,最好联用R中的内建数值/字符处理函数和apply族函数。

条件执行

在条件执行结构中,一条或一组语句仅在满足一个指定条件时执行。条件执行结构包括if-else、ifelse和switch。

if-else结构

控制结构if-else在某个给定条件为真时执行语句。也可以同时在条件为假时执行另外的语句。语法为:

if (cond) statement
if (cond) statement1 else statement2

示例如下:

if (is.character(grade)) grade <- as.factor(grade)
if (!is.factor(grade)) grade <- as.factor(grade) else print("Grade already is a factor")

在第一个实例中,如果grade是一个字符向量,它就会被转换为一个因子。在第二个实例中,两个语句择其一执行。如果grade不是一个因子(注意符号!),它就会被转换为一个因子。如果它是一个因子,就会输出一段信息。

ifelse结构

ifelse结构是if-else结构比较紧凑的向量化版本,其语法为:

ifelse(cond, statement1, statement2)

若cond为TRUE,则执行第一个语句;若cond为FALSE,则执行第二个语句。示例如下:

ifelse(score > 0.5, print("Passed"), print("Failed"))
outcome <- ifelse (score > 0.5, "Passed", "Failed")

在程序的行为是二元时,或者希望结构的输入和输出均为向量时,请使用ifelse。

switch结构

switch根据一个表达式的值选择语句执行。语法为:

switch(expr, ...)

其中的...表示与expr的各种可能输出值绑定的语句。通过观察代码清单7中的代码,可以轻松地理解switch的工作原理。

> feelings <- c("sad", "afraid")
> for (i in feelings)
    print(
      switch(i,
        happy  = "I am glad you are happy",
        afraid = "There is nothing to fear",
        sad    = "Cheer up",
        angry  = "Calm down now"
      )
    )

[1] "Cheer up"
[1] "There is nothing to fear"

虽然这个例子比较幼稚,但它展示了switch的主要功能。你将在下一节学习如何使用switch编写自己的函数。

用户自编函数

R的最大优点之一就是用户可以自行添加函数。事实上,R中的许多函数都是由已有函数构成的。一个函数的结构看起来大致如此:

myfunction <- function(arg1, arg2, ... ){
  statements
  return(object)
}

函数中的对象只在函数内部使用。返回对象的数据类型是任意的,从标量到列表皆可。让我们看一个示例。

假设你想编写一个函数,用来计算数据对象的集中趋势和散布情况。此函数应当可以选择性地给出参数统计量(均值和标准差)和非参数统计量(中位数和绝对中位差)。结果应当以一个含名称列表的形式给出。另外,用户应当可以选择是否自动输出结果。除非另外指定,否则此函数的默认行为应当是计算参数统计量并且不输出结果。代码清单8给出了一种解答。

mystats <- function(x, parametric=TRUE, print=FALSE) {
  if (parametric) {
    center <- mean(x); spread <- sd(x)
  } else {
    center <- median(x); spread <- mad(x)
  }
  if (print & parametric) {
    cat("Mean=", center, "\n", "SD=", spread, "\n")
  } else if (print & !parametric) {
    cat("Median=", center, "\n", "MAD=", spread, "\n")
  }
  result <- list(center=center, spread=spread)
  return(result)
}

要看此函数的实战情况,首先需要生成一些数据(服从正态分布的,大小为500的随机样本):

set.seed(1234)
x <- rnorm(500)

在执行语句:

y <- mystats(x)

之后,y$center将包含均值(0.00184),y$spread将包含标准差(1.03),并且没有输出结果。如果执行语句:

y <- mystats(x, parametric=FALSE, print=TRUE)

y$center将包含中位数(–0.0207),y$spread将包含绝对中位差(1.001)。另外,还会输出以下结果:

Median= -0.0207
MAD= 1

下面让我们看一个使用了switch结构的用户自编函数,此函数可让用户选择输出当天日期的格式。在函数声明中为参数指定的值将作为其默认值。在函数mydate()中,如果未指定type,则long将为默认的日期格式:

mydate <- function(type="long") {
  switch(type,
    long  = format(Sys.time(), "%A %B %d %Y"),
    short = format(Sys.time(), "%m-%d-%y"),
    cat(type, "is not a recognized type\n")
  )
}

实战中的函数如下:

> mydate("long")
[1] "Thursday December 02 2010"
> mydate("short")
[1] "12-02-10"
> mydate()
[1] "Thursday December 02 2010"
> mydate("medium")
medium is not a recognized type

请注意,函数cat()仅会在输入的日期格式类型不匹配"long"或"short"时执行。使用一个表达式来捕获错误输入的参数值通常来说是一个好主意。

整合与重构

R中提供了许多用来整合(aggregate)和重塑(reshape)数据的强大方法。在整合数据时,往往将多组观测替换为根据这些观测计算的描述性统计量。在重塑数据时,则会通过修改数据的结构(行和列)来决定数据的组织方式。本节描述了用来完成这些任务的多种方式。

在接下来的两个小节中,我们将使用已包含在R基本安装中的数据框mtcars。这个数据集是从Motor Trend 杂志(1974)提取的,它描述了34种车型的设计和性能特点(汽缸数、排量、马力、每加仑汽油行驶的英里数,等等)。要了解此数据集的更多信息,请参阅help(mtcars)。

转置

转置(反转行和列)也许是重塑数据集的众多方法中最简单的一个了。使用函数t()即可对一个矩阵或数据框进行转置。对于后者,行名将成为变量(列)名。代码清单9展示了一个例子。

> cars <- mtcars[1:5,1:4]
> cars
                   mpg cyl disp  hp
Mazda RX4         21.0   6  160 110
Mazda RX4 Wag     21.0   6  160 110
Datsun 710        22.8   4  108  93
Hornet 4 Drive    21.4   6  258 110
Hornet Sportabout 18.7   8  360 175
> t(cars)
     Mazda RX4 Mazda RX4 Wag Datsun 710 Hornet 4 Drive Hornet Sportabout
mpg         21            21       22.8           21.4              18.7
cyl          6             6        4.0            6.0               8.0
disp       160           160      108.0          258.0             360.0
hp         110           110       93.0          110.0             175.0

为了节约空间,代码清单9仅使用了mtcars数据集的一个子集。在本节稍后讲解reshape包的时候,你将看到一种更为灵活的数据转置方式。

整合数据

在R中使用一个或多个by变量和一个预先定义好的函数来折叠(collapse)数据是比较容易的。调用格式为:

aggregate(x, by, FUN)

其中x是待折叠的数据对象,by是一个变量名组成的列表,这些变量将被去掉以形成新的观测,而FUN则是用来计算描述性统计量的标量函数,它将被用来计算新观测中的值。

作为一个示例,我们将根据汽缸数和挡位数整合mtcars数据,并返回各个数值型变量的均值(见代码清单10)。

> options(digits=3)
> attach(mtcars)
> aggdata <-aggregate(mtcars, by=list(cyl,gear), FUN=mean, na.rm=TRUE)
> aggdata
  Group.1 Group.2  mpg cyl disp  hp drat   wt qsec  vs   am gear carb
1       4       3 21.5   4  120  97 3.70 2.46 20.0 1.0 0.00    3 1.00
2       6       3 19.8   6  242 108 2.92 3.34 19.8 1.0 0.00    3 1.00
3       8       3 15.1   8  358 194 3.12 4.10 17.1 0.0 0.00    3 3.08
......

在结果中,Group.1表示汽缸数量(4、6或8),Group.2代表挡位数(3、4或5)。举例来说,拥有4个汽缸和3个挡位车型的每加仑汽油行驶英里数(mpg)均值为21.5。

在使用aggregate()函数的时候,by中的变量必须在一个列表中(即使只有一个变量)。你可以在列表中为各组声明自定义的名称,例如by=list(Group.cyl=cyl, Group.gears=gear)。指定的函数可为任意的内建或自编函数,这就为整合命令赋予了强大的力量。但说到力量,没有什么可以比reshape包更强。

reshape包

reshape包是一套重构和整合数据集的绝妙的万能工具。由于它的这种万能特性,可能学起来会有一点难度。我们将慢慢地梳理整个过程,并使用一个小型数据集作为示例,这样每一步发生了什么就很清晰了。由于reshape包并未包含在R的标准安装中,在第一次使用它之前需要使用install.packages("reshape")进行安装。

大致说来,你需要首先将数据“融合”(melt),以使每一行都是一个唯一的标识符—变量组合。然后将数据“重铸”(cast)为你想要的任何形状。在重铸过程中,你可以使用任何函数对数据进行整合。将使用的数据集如表8所示。

ID Time X1 X2
1 1 5 6
1 2 3 5
2 1 6 1
2 2 2 4

在这个数据集中,测量(measurement)是指最后两列中的值(5、6、3、5、6、1、2、4)。每个测量都能够被标识符变量(在本例中,标识符是指ID、Time以及观测属于X1还是X2)唯一地确定。举例来说,在知道ID为1、Time为1,以及属于变量X1之后,即可确定测量值为第一行中的5。

融合

数据集的融合是将它重构为这样一种格式:每个测量变量独占一行,行中带有要唯一确定这个测量所需的标识符变量。要融合表8中的数据,可使用以下代码:

library(reshape)
md <- melt(mydata, id=(c("id", "time")))

你将得到如表9所示的结构。

ID Time 变  量
1 1 X1 5
1 2 X1 3
2 1 X1 6
2 2 X1 2
1 1 X2 6
1 2 X2 5
2 1 X2 1
2 2 X2 4

注意,必须指定要唯一确定每个测量所需的变量(ID和Time),而表示测量变量名的变量(X1或X2)将由程序为你自动创建。

既然已经拥有了融合后的数据,现在就可以使用cast()函数将它重铸为任意形状了。

重铸

cast()函数读取已融合的数据,并使用你提供的公式和一个(可选的)用于整合数据的函数将其重塑。调用格式为:

newdata <- cast(md, formula, FUN)

其中的md为已融合的数据,formula描述了想要的最后结果,而FUN是(可选的)数据整合函数。其接受的公式形如:

rowvar1 + rowvar2 +~ colvar1 + colvar2 +

在这一公式中,rowvar1 + rowvar2 + ...定义了要划掉的变量集合,以确定各行的内容,而colvar1 + colvar2 + ...则定义了要划掉的、确定各列内容的变量集合。参见图1中的示例。

由于右侧(d、e和f)的公式中并未包括某个函数,所以数据仅被重塑了。反之,左侧的示例(a、b和c)中指定了mean作为整合函数,从而就对数据同时进行了重塑与整合。例如,(a)中给出了每个观测所有时刻中在X1和X2上的均值;示例(b)则给出了X1和X2在时刻1和时刻2的均值,对不同的观测进行了平均;在(c)中则是每个观测在时刻1和时刻2的均值,对不同的X1和X2进行了平均。

参考文档