数组

数组(array)是一个带有多个下标且型态相同的元素集合，例如数值所构成的数组。在

R 中有一些简单的函数可以建立与处理数组，特别是针对矩阵的处理(矩阵在 R 中是数组的一种)。

数组有一个特别的属性称为维度向量(dimension vector)，此向量是一个由正整数所构成的向量，如果它的长度为

k，那麽该数组就是

k 维的，例如一个矩阵是 2 维数组。数组中元素的下标可以从 1 一直标到维度向量中对应元素的值。若要查看一个向量的维度，可使用 dim() 函数。

事实上数组就是定义了维度属性(dim)之后的向量。假设

z 是一个包含 20 个元素的向量：

1 z

我们可以使用 dim() 函数指定其维度属性，让此向量变成一个 4 × 5 的数组，而 R 会将原本向量中的资料依照 column major order 的方式排列：

1 dim(z)

现在的 z 为(请留意各元素排列的顺序)：

[,1] [,2] [,3] [,4] [,5]

[1,]

1

5

9

13

17

[2,]

2

6

10

14

18

[3,]

3

7

11

15

19

[4,]

4

8

12

16

20

也可已将 z 指定为三维数组：

1 dim(z)

现在的 z 则为：

, , 1

[,1] [,2] [,3] [,4] [,5]

[1,]

1

3

5

7

9

[2,]

2

4

6

8

10

, , 2

[,1] [,2] [,3] [,4] [,5]

[1,]

11

13

15

17

19

[2,]

12

14

16

18

20

在指定维度向量时，对于数组内部所储存的资料其实没有影响，任何维度的数组资料实际上都是储存在一个资料向量中，改变维度只是影响数组中每个资料元素对应到

资料向量中的储存位置而已，R 在存取数组元素时，是依照维度属性计算指定的元素在资料向量中的位置，再去存取资料向量中储存的元素。

另外还有其他数组相关函数，如 matrix() 与 array()，详细介绍请参阅

array() 函数章节。

当数组的维度是一维时，这种数组的处理和向量完全一致(包括输出的方式)，而这有时会导致使用者溷淆。

数组索引

数组的元素可以透过中括号指定下标来存取，各个下标用逗点隔开，例如：

1 arr

2 arr[2, 3]

除此之外，也可以利用索引向量取得指定的子数组：

1 arr

2 subarr

数组 arr 为

[,1] [,2] [,3] [,4] [,5]

[1,]

1

6

11

16

21

[2,]

2

7

12

17

22

[3,]

3

8

13

18

23

[4,]

4

9

14

19

24

[5,]

5

10

15

20

25

而子数组 subarr 为

[,1] [,2] [,3]

[1,]

6

11

16

[2,]

7 12

17

[3,]

8

13

18

若是不指定下标，则 R 会选取对应下标的所有元素，例如：

1 arr

2 subarr2

子数组 subarr2 为

[,1] [,2] [,3] [,4] [,5]

[1,]

1

6

11

16

21

[2,]

2

7

12

17

22

[3,]

3

8

13

18

23

若所有下标都不指定 arr[,,] 则表示整个数组，这与忽略下标直接使用 arr 效果是一样的。

如果一个多维数组只给了一个下标或索引向量，在这种情况下，R 会把数组当作向量来使用，而忽略维度属性。

1 arr

2 arr[3]

此时 arr 为

[,1] [,2] [,3] [,4]

[1,]

1

10

5

7

[2,]

2

4

6

8

但是，如果单个索引不是一个向量而是一个数组，结果可能就不是这样了，具体的实例可以看下一章节的讨论。

索引数组

数组也可以使用索引数组的方式存取任意不规则元素集合，这种方式我们直接使用范例来说明比较清楚，首先假设我们有一个 4

× 5 的二维数组：

1 x

向量 x 为

[,1] [,2] [,3] [,4] [,5]

[1,]

1

5

9

13

17

[2,]

2

6

10

14

18

[3,]

3

7

11

15

19

[4,]

4

8

12

16

20

我们可以使用索引数组的方式将 x[1,3]、x[2,2]

与 x[3,1] 取出来构成一组新的向量，首先产生一个索引数组 i：

1 i

数组 i 为

[,1] [,2]

[1,]

1

3

[2,]

2

2

[3,]

3

1

索引数组的行(column)数必须等于被存取数组的维度，其列数则没有限制，每一列(row)指定一个元素，一列中的各个数字就是对应元素的下标。以索引数组取出指定的元素：

1 y

新的向量 y 为

[1] 9 6 3

除了取出元素，索引数组也可以改变数组中指定元素的值：

1 x[i]

此时 x 变为

[,1] [,2] [,3] [,4] [,5]

[1,]

1

5

30

13

17

[2,]

2

31

10

14

18

[3,]

32

7

11

15

19

[4,]

4

8

12

16

20

array() 函数

要建立一个数组，除了将向量指定 dim 属性的方法之外，也可直接透过 array() 函数将向量转换而得到，此函数第一个参数指定资料向量，而第二个参数指定数组维度：

1 h

2 arr

其作用等同于

1 h

2 arr

3 dim(arr)

如果指定的资料向量长度不够时，那麽 R 会将资料向量重複使用直到长度符合维度向量所指定的大小为止，我们看一个特别的例子：

1 Z

这会产生一个所有元素都是 0 的三维数组。

数组可以用于算术运算式中，其运算方式与向量的运算方式相似，对每一个数组中的元素个别做运算，最后得到的结果也是一个数组，所有参与计算的数组其为度必须相同，而这个维度也就是最后计算结果的维度，例如：

1 arr1

2 arr2

3 arr3

4 arr4

数组 arr4 为

[,1] [,2] [,3] [,4]

[1,]

12

38

80 138

[2,]

23

57

107 173

数组的外积

外积运算(outer product)是数组非常重要的运算之一，透过外积运算子

%o% 可以得到两个向量的外积：

1 a

2 b

3 ab

得到的矩阵 ab 为

[,1] [,2] [,3] [,4]

[1,]

1

2

3

4

[2,]

2

4

6

8

[3,]

3

6

9

12

另一更一般化的种作法是使用 outer() 函数：

1 ab

outer() 函数的第三个参数是指定运算函数，若指定为乘法(*)就是一般的外积，我们可以自订这个运算函数，例如假设我们想看函数

f(x,y) = cos(y)/(1+x^2) 在 [-5,5] x [-5,5] 平面上的图形，首先我们先指定 x 与 y 方向的格点与 f 函数：

1 x

2 y

3 f

再以 outer() 计算每个点的数值：

1 z

最后再把图形画出来

zlim

zlen

# height color lookup table

colorlut

# assign colors to heights for each point

col

persp3d(x, y, z, aspect = c(1, 1, 0.5),

col = "lightblue", xlab = "X", ylab = "Y",

zlab = "cos(y)/(1 + x^2)")

得到的图形为

关于 3D 绘图函数 persp3d() 的用法请参考 RGL 套件。

关于数学上的外积可参考 Wiki 上的说明。

广义转置

aperm(A, perm) 函数可以用来重新排列 A 数组。第一个参数指定被重新排列的数组，而第二个参数则是指定排列顺序，此参数是一个由 {1, ..., k}

所组成的排列(permutation)，其中

k 为数组的维度，这个函数会产生一与 A 相同道小的数组，但原本 A 的第 perm[j] 个维度将会变成新数组的第 j 个维度，这种操作实际上是对矩阵的一种广义转置。例如若 A 为一个二维数组，也就是矩阵：

1 A

则 A 的转置矩阵 B 为

1 B

而若是只是简单矩阵转置，有个更简单的函数 t() 可以使用：

1 B

所得到的结果是一样的。

矩阵工具矩阵就是二维的数组，然而他在各个领域被广泛的使用，因此 R 中也很许多特别针对矩阵运算的运算子与函数，例如前面提到的矩阵转置函数 t()，以及

nrow() 与 ncol() 函数(这两个函数分别传回矩阵的列(row)数与行(column)数)。

以下介绍一些矩阵的运算：

矩阵乘法

矩阵的乘法运算子为 %*%，例如：

A

x

b

若 X 与 Y 为相同大小的方阵：

X

Y

则

1 X * Y

为矩阵中对应元素相乘，而

1 X %*% Y

则是数学上的矩阵相乘。

n × 1 或者 1 × n 的矩阵在情况允许的时候可以作为一个长度为 n 的向量使用，反之，向量在矩阵相乘运算式中出现时，可能会被自动转换成与矩阵对应的行向量或列向量。

crossprod(x, y) 函数所作的运算与 t(x) %*% y 相同，但是在运算上效率更高。如果 crossprod() 第二个参数不指定，其预设与第一个参数一样。

diag() 函数的功能会因为其传入的参数不同而有所差异，当传入一个向量时，diag() 传回以该向量元素为对角元素的对角矩阵：

1 X

得到的 X 为

[,1] [,2] [,3]

[1,]

1

0

0

[2,]

0

2

0

[3,]

0

0

3

当传入的参数是一个矩阵时，diag() 会传回此矩阵的对角元素：

X

d

得到的 d 为

[1] 1 5 9

这和 Matlab 中 diag() 的用法完全一致。不过有点溷乱的是，如果传入的参数是单一个值，那么 diag(k) 的结果就是 k × k 的单位矩阵。

线性方程式与反矩阵

线性联立方程式求解的问题可使用 solve() 函数，例如：

A

x

b

几设我们只知道 A 与 b 要求 x 则使用：

1 x.sol

所得到的 x.sol 即为 x 的解。

另外 solve() 函数亦可用来计算反矩阵，用法就是直接将矩阵传入：

1

solve(A)

不过这一般比较少用到，若使用

1 x

方式来求解跟 solve(A,b) 相较之下，其计算效率较低而且潜在数值上的不稳定性。

用于多变数计算的二次方程式

1 x %*% A^{-1} %*% x

可以透过像

1 x %*%

solve(A, x)

的方式计算得到，而不直接计算 A 的反矩阵。

特徵值与特征向量

eigen() 函数可用来计算矩阵的特征值和特征向量，这个函数的传回值是一个含有 values 和 vectors 两个元素的列表(list)，例如：

X

ev

得到的 ev 为一个列表，而 ev$values 为特征值：

[1] 10.8150729 0.1849271

ev$vectors 为特征向量：

[,1]

[,2]

[1,] -0.3773994 -0.9798642

[2,] -0.9260506 0.1996652

如果我们只需要特征值，我们可以使用如下的命令：

eig.vals

所得到的 eig.vals 为特征值，而特征向量则被捨弃。如果以下面的运算式作为一个命令，

eigen(X)

则会输出所有特征值与特征向量。

SVD 分解与行列式

svd() 函数用在矩阵的 SVD 分解(singular value decomposition，请参考

Wiki 上的说明)，svd() 所传回的结果是由 d, u 与 v 构成的一个列表，d 为奇异值(singular value)构成的向量，u

与 v 为矩阵(其 columns 为 orthonormal)，而

u 的行空间(column space)等于

X 的行空间，而 v 的行空间等于 X 的列空间(row space)，我们以实例说明比较清楚：

X

X.svd

得到的 X.svd 为

$d

[1] 9.933252 4.508935

$u

[,1]

[,2]

[1,] -0.8206658 0.5674691

[2,] -0.3271692 -0.5627514

[3,] -0.4684741 -0.6010738

$v

[,1]

[,2]

[1,] -0.7882054 0.6154122

[2,] -0.6154122 -0.7882054

而原本的 X 就会等于

X.svd$u %*% diag(X.svd$d) %*% t(X.svd$v)

如果 M 是一个方阵，我们可以使用以下方式计算 M 行列式的绝对值：

absdetM

若是只要计算行列式，在 R 中有内建的函数 det() 可以使用。

Least Squares Fitting 与 QR 分解

lsfit() 函数可以处理 Least Squares Fitting 的问题：

ans

这样就得到 Least Squares Fitting 的结果，其中 y 是观测向量，X 是设计矩阵，更详细的说明可以查看线上使用手册，另外也可以参考回归诊断(regression

diagnostics)函数

ls.diag()。

实际上，你在回归分析中可能已经习惯使用 lm() (参见线性模型)

而不是 lsfit()。

区块矩阵

cbind() 与 rbind() 函数可以把向量或矩阵拼成一个新的矩阵，cbind() 是把向量或矩阵横向拼成一个大矩阵，而 rbind() 则是以纵向的方式，这两个函数传入的所有参数必须是任意长度的向量或是矩阵，如果传入的参数是向量，则 cbind()

会把向量当成行向量(column vector)使用，而

rbind() 会把向量当成列向量使用(row vector)，例如：

a

b

c

d

得到的 c 为

a b

[1,] 1 4

[2,] 2 5

[3,] 3 6

而 d 为

[,1]

[,2] [,3]

a

1

2

3

b

4

5

6

如果传入的向量长度不一，则较短的向量会被重複使用，直到与最长的向量长度相同。例如：

w

x

y

z

得到的 y 为

w x

[1,] 1 5

[2,] 2 6

[3,] 1 7

[4,] 2 8

而 z 为

[,1]

[,2] [,3] [,4]

w

1

2

1

2

x

5

6

7

8

若传入的是矩阵，则 cbind() 只能接受所有矩阵的列(row)数都相同的情形，而

rbind() 则只能接受所有矩阵的行(column)数都相同的情形：

A

B

C

得到的 C 为

[,1] [,2] [,3]

[1,] 1

4

7

[2,]

2

5

8

[3,]

3

6

9

rbind() 与 cbind() 的传回值都是矩阵型态，因此 cbind(x) 和 rbind(x) 是个把向量 x 分别转换成行向量或列向量最简单的方法。

数组与 c() 函数

前面提到的 cbind() 与 rbind() 函数是考虑 dim 属性的连接函数，而 c() 函数则是不考虑这些数值物件的 dim 和 dimnames 属性，这一点在有些时候非常有用。

将一个数组强制转换成一般向量的标准方式是使用函数 as.vector()：

X

Xvec

但也可以使用 c() 函数：

[splua]Xvec

两者功用相同，但习惯上一般建议使用前者。

因子分布表

前面提到过因子(factor)会把资料分群，而

table() 函数可以从同样长度的不同因子中计算出频率表，如果输入 k 个因子参数，会得到 k 维的频率分布数组。

延续因子章节的范例，statef 是由州名所产生的因子，可使用 table() 函数得到频率分布表：

statefr

得到的 statefr 为

statef

act nsw

nt qld sa

tas vic

wa

2

6

2

5

4

2

5

4

这些频率会依照因子的 levels 排序与标记。下面这个指令也有同样的作用，只是比较繁琐：

statefr

延续因子章节 tapply() 函数的例子，我们使用 cut() 函数产生一个将收入(incomes)分成几个群组：

incomef

breaks = 35 + 10 * (0:7)))

然后我们结合 statef 与 incomef 可得到一个二维频率分布表：

table(incomef, statef)

得到的频率分布表为

statef

incomef act nsw nt qld sa tas vic wa

(35,45]

1

1

0

1

0

0

1 0

(45,55]

1

1

1

1

2

0

1 3

(55,65]

0

3

1

3

2

2

2 1

(65,75]

0

1

0

0

0

0

1 0

使用类似的方法也可以得到更高维的频率分布表。

向量与数组混合运算

向量与数组间对应元素混合运算的确切规则有点诡异，并且很难在参考文献中找到很明确的的解释。根据经验，我们列出下面的一些比较可靠的说明：

运算式运算是从左到右进行的。

长度较短的向量运算元会被重复使用直到与较长的运算元长度相同为止。

当只有较短的向量与数组出现时，所有的数组必须有相同的维度，否则会得到错误的结果。

向量运算元比矩阵或数组运算元长时会引起错误。

如果运算式中包含数组结构，且计算时也没有关于向量的错误资讯和强制转换操作，则所得到的计算结果会是一个与数组运算元维度相同的数组。