1. 程式人生 > >R中的多維陣列和矩陣

R中的多維陣列和矩陣

1.生成陣列或矩陣

陣列可以看成是帶多個下標的型別相同的元素的集合,常用的數值型的陣列如矩陣,也可以有其他型別(字元型、邏輯型、複數型)。R可以很容易地生成和處理陣列,特別是矩陣(二維陣列)。

陣列有一個特殊的屬性叫做維數向量,維數向量是一個元素取正整數值得向量,其長度是陣列的維數,比如維數向量有兩個元素時陣列為二維陣列(矩陣)。維數向量的每一個元素指定了該下標的上界,下標的下界總為1.

(()將向量定義成陣列

向量只有定義了維數向量後才能被看作為陣列,比如:

> z<-1:12

> dim(z)<-c(3,4)

> z

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

[1,]   1    4    7  10

[2,]   2    5    8  11

[3,]   3    6    9  12

注意:矩陣的元素時按列存放的。也可以吧向量定義為一維陣列,例如:

> dim(z)<-12

> z

 [1] 1  2  3 4  5  6 7  8  9 10 11 12

(2)用array()函式構造多維陣列

array(data,dim,dimnames)

其中data是一個向量資料,dim()是陣列各維的長度,預設時為原向量的長度。dimnames是陣列維的名字,預設為空。如:

> x<-array(1:20,dim=c(4,5))

> 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

另一種方式:

> z<-array(0,dim=c(3,4,2))

> z

, , 1

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

[1,]   0    0    0   0

[2,]   0    0    0   0

[3,]   0    0    0   0

, , 2

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

[1,]   0    0    0   0

[2,]   0    0    0   0

[3,]   0    0    0   0

如此定義了一個3*4*2的三維陣列,元素均為0。這種方法常用來對陣列作初始化。

(3)用matrix()函式構造矩陣

matrix(data=NA,nrow=1,ncol=1,byrow=FALSE,dimnames=NULL)

其中data是一個向量資料,nrow是矩陣的行數,ncol為列數,nrow與ncol的乘積是矩陣元素個數,byrow項控制排列元素是否按行進行,當為T時按行放置,當為F時按列放置,dimnames給定行和列的名稱,預設時為空。如:

> A<-matrix(1:15,nrow=3,ncol=5,byrow=TRUE)

> A

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

[1,]   1    2    3   4    5

[2,]   6    7    8   9   10

[3,]  11   12   13  14   15

以下兩種格式與前面的格式是等價的:

> A<-matrix(1:15,nrow=3,byrow=TRUE)

> A<-matrix(1:15, ncol=5,byrow=TRUE)

如果把語句中的byrow=TRUE去掉或者改為byrow=FALSE,則資料按列放置:

> A<-matrix(1:15,nrow=3,ncol=5)

> A

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

[1,]   1    4    7  10   13

[2,]   2    5    8  11   14

[3,]   3    6    9  12   15

(4)用合併命令構造矩陣

用rbind()、cbind()將兩個或兩個以上的向量或矩陣合併起來,rbind()表示按行合併,cbind()表示按列合併。例:

> x<-c(1:5)

> y<-c(6:10)

> rbind(x,y)

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

x   1    2    3   4    5

y   6    7    8   9   10

> cbind(x,y)

    x  y

[1,] 1 6

[2,] 2 7

[3,] 3 8

[4,] 4 9

[5,] 5 10

2.陣列下標

陣列和向量一樣,可以對陣列中的某些元素進行訪問,或進行運算。

(1)陣列下標

要訪問陣列的某個元素,只要寫出陣列名和方括號內的用逗號分開的下標即可。如:

> a<-1:24

> dim(a)<-c(2,3,4)

> a[2,1,2]

[1] 8

更進一步還可以在每一個下標位置寫一個下標向量,表示這一維取出的所有指定下標的元素,如:

> a[1,2:3,2:3]

    [,1] [,2]

[1,]   9   15

[2,]  11   17

取出了第一下標為1,第二下標為2或3.第三下標為2或3的元素。

注意,因為第一維只有一個下標,所以退化成一個2*2的陣列。

另外,如果略寫某一維的下標,則表示該維全選,如:

> a[1,,]

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

[1,]   1    7   13  19

[2,]   3    9   15  21

[3,]   5   11   17  23

取出所有第一下標為1的元素,得到一個3*4的陣列。

> a[,2,]

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

[1,]   3    9   15  21

[2,]   4   10   16  22

取出所有第二下標為2的元素得到一個2*4的陣列。

> a[1,1,]

[1] 1  7 13 19

得到一個長度為4的向量,不再是陣列。

a[,,]或a[]都表示整個陣列,如:

> a[]<-0

可以在不改變陣列維數的條件下把元素都賦為0.

(2)不規則的陣列下標

在R中可以把陣列中的任意位置的元素作為陣列訪問,方法是用一個二維陣列作為陣列的下標,二維陣列的每一行是一個元素的下標,列數為陣列的維數。

> z<-1:24

> dim(z)<-c(2,3,4)

> b<-matrix(c(1,1,1,2,2,3,1,3,4,2,1,4),ncol=3,byrow=T);b

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

[1,]   1    1    1

[2,]   2    2    3

[3,]   1    3    4

[4,]   2    1    4

> z[b]

[1] 1 16 23 20

表示把2*3*4的陣列z的第[1,1,1],[2,2,3],[1,3,4],[2,1,4]這四個元素取出。

注意:取出的是一個向量,我們還可以對這幾個元素賦值,如:

> z[b]<-c(101,102,103,104)

> z[b]<-0

3.陣列的四則運算

可以對陣列之間進行四則運算,這時進行的是陣列對應元素的四則運算,參加運算的陣列一般應該是相同形狀的(dim屬性完全相同)。如:

> A=B=matrix(1:12,nrow=3,ncol=4)

> A+B

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

[1,]   2    8   14  20

[2,]   4   10   16  22

[3,]   6   12   18  24

> A-B

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

[1,]   0    0    0   0

[2,]   0    0    0   0

[3,]   0    0    0   0

> C<-A*B;D=A/B;C;D

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

[1,]   1   16   49 100

[2,]   4   25   64 121

[3,]   9   36   81 144

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

[1,]   1    1   1    1

[2,]   1    1    1   1

[3,]   1    1    1   1

陣列的加、減運算和數乘運算滿足原矩陣運算的性質,但陣列的乘、除運算實際上是陣列中對應位置的元素作運算。

形狀不一致的向量(或陣列)中的資料進行四則運算,一般的規則是將向量(或陣列)中的資料與對應向量(或陣列)中的資料進行運算,把段向量(或陣列)的資料迴圈使用,從而可以與長向量(或陣列)資料進行匹配,並儘可能保留共同的陣列屬性。如:

> x1<-c(100,200)

> x2<-1:6

> x1+x2

[1] 101 202 103 204 105 206

> x1+x3

    [,1] [,2]

[1,] 101  204

[2,] 202  105

[3,] 103  206

可以看到,當向量與陣列共同運算時,向量按列匹配。當兩個陣列不匹配時,R會提出警告。如:

> x2<-1:5

> x1+x2

[1] 101 202 103 204 105

警告資訊:

In x1 + x2 : 長的物件長度不是短的物件長度的整倍

4.矩陣的運算

(1)轉置

求矩陣A的轉置用函式t()。例如:

> A=matrix(1:12,nrow=3,ncol=4);A

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

[1,]   1    4    7  10

[2,]   2    5    8  11

[3,]   3    6    9  12

> t(A)

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

[1,]   1    2    3

[2,]   4    5    6

[3,]   7    8    9

[4,]  10   11   12

(2)方陣的行列式

det()函式可求方陣行列式的值,例如:

> det(matrix(1:4,ncol=2))

[1] -2

(3)向量的內積

對於n維向量x,可以看成n*1階矩陣或1*n階矩陣。若x與y是相同維數的向量,則x%*%y表示x與y作內積,例如:

> x<-1:5;y<-2*1:5

> x%*%y

    [,1]

[1,] 110

函式crossprod()是內積計算函式(表示交叉乘積),crossprod(x,y)加us那向量x與y的內積,即‘t(x)%*%y‘。crossprod(x)表示x與x的內積。

類似地,tcrossprod(x,y)表示‘x%*% t(y) ‘,即x與y的外積,也稱為叉積,tcrossprod(x)表示x與x的外積。

(4)向量的外積(叉積)

設x,y是n維向量,則x%*o%y表示x與y作外積。如:

> x%o%y

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

[1,]   2    4    6   8   10

[2,]   4    8   12  16   20

[3,]   6   12   18  24   30

[4,]   8   16   24  32   40

[5,]  10   20   30  40   50

函式outer()是外積運算函式,outer(x,y)計算x與y的外積,等價於x%o%y。

函式outer()的一般格式為:outer(x,y,fun)

其中x,y矩陣(或向量),fun是作外積運算函式,預設值為乘法運算。函式outer()在繪製三維曲面時非常有用。

(5)矩陣的乘法

A%*%B表示通常意義下的兩個矩陣的乘積(當然要求矩陣A的列數等於矩陣B的行數),如:

> A=matrix(1:12,nrow=3,ncol=4)

> B=matrix(1:12,nrow=4,ncol=3)

> A%*%B

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

[1,]  70  158  246

[2,]  80  184  288

[3,]  90  210  330

(6)生成對角陣和矩陣取對角運算

若取一個方陣的對角元素,用函式diag()。

對一個向量用函式diag(),將產生以這個向量為對角元素、其餘元素全為0的對角矩陣。

對一個正整數k應用函式diag(),將產生k維單位矩陣。

例:

> A=matrix(1:16,nrow=4,ncol=4)

> A

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

[1,]   1    5    9  13

[2,]   2    6   10  14

[3,]   3    7   11  15

[4,]   4    8   12  16

> diag(A)

[1] 1  6 11 16

> diag(3)

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

[1,]   1    0    0

[2,]   0    1    0

[3,]   0    0    1

(7)解線性方程組和求矩陣的逆矩陣

陣求逆用函式solve()。

用solve(A,b)運算,結果是可解線性方程組Ax=b的解x。若b為單位矩陣,則結果就是A的逆矩陣。solve()中,若b預設,系統預設為單位矩陣,則這個過程便為矩陣求逆。

例:

rnorm()為隨機生成正態

> A=matrix(rnorm(16),4,4);A

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

[1,] 0.67168810  0.77743512 -2.1351143  1.1196029

[2,] 0.70557911 -0.40344244  0.5866727 1.2809767

[3,] 0.57529587 -0.02789363 -1.0624213-1.2417952

[4,] 0.08381951  0.65116106 0.1662728 -0.9488299

> solve(A)

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

[1,] 0.02554831  0.89633961  0.5293070 0.5475198

[2,] 0.34677139 -0.02034261 -0.5381588 1.0860438

[3,] -0.22913613  0.34001906 -0.2186800  0.4748699

[4,] 0.20008471  0.12480673 -0.3608891-0.1770177

(8)矩陣的特徵值和特徵向量

矩陣A的譜分析為A=UDU’,其中D為由A的特徵值組成的對角矩陣,U的列為A的特徵值對應的特徵向量,在R中可用函式eigen()得到U和D。

eigen(x,symmetric,only.value=FALSE)

x為矩陣,symmetric項指定矩陣x是否為對稱矩陣,若不指定,系統將自動檢測x是否為對稱矩陣。only.value項如果是TRUE,只有特徵值計算並返回,否則返回的特徵值和特徵向量。

例:

> A=diag(4)+1

> A

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

[1,]   2    1    1   1

[2,]   1    2    1   1

[3,]   1    1   2    1

[4,]   1    1    1   2

> A.e=eigen(A,symmetric=T)

> A.e

$values

[1] 5 1 1 1

$vectors

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

[1,] -0.5 0.000000e+00  0.0000000  0.8660254

[2,] -0.5 -6.408849e-17  0.8164966 -0.2886751

[3,] -0.5 -7.071068e-01 -0.4082483-0.2886751

[4,] -0.5 7.071068e-01 -0.4082483 -0.2886751

>A.e$vectors%*%diag(A.e$values)%*%t(A.e$vectors)

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

[1,]   2    1    1   1

[2,]   1    2    1   1

[3,]   1    1    2   1

[4,]   1    1    1    2

(9)矩陣的Choleskey分解

對於正定矩陣A,可對其進行Choleskey分解,即A=P’P,其中,P為是上三角矩陣,在R中可用函式chol()來進行Choleskey分解。

正定矩陣是指

例:

> A.c=chol(A)

> A.c

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

[1,] 1.414214 0.7071068 0.7071068 0.7071068

[2,] 0.000000 1.2247449 0.4082483 0.4082483

[3,] 0.000000 0.0000000 1.1547005 0.2886751

[4,] 0.000000 0.0000000 0.0000000 1.1180340

> t(A.c)%*%A.c

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

[1,]   2    1    1   1

[2,]   1    2    1   1

[3,]   1    1    2   1

[4,]   1    1    1   2

(10)矩陣奇異值分解

A為m*n矩陣,A的秩為n,即rank(A)=n,可分解為A=UDV’,其中U’U=V’V=I。在R中可以用函式svd()進行奇異值分解。

例:

> A=matrix(1:18,3,6)

> A

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

[1,]   1    4    7  10   13   16

[2,]   2    5    8  11   14   17

[3,]   3    6    9  12   15   18

> A.s=svd(A)

> A.s

$d

[1] 4.589453e+01 1.640705e+00 3.627301e-16

$u

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

[1,] -0.5290354  0.74394551 0.4082483

[2,] -0.5760715  0.03840487 -0.8164966

[3,] -0.6231077 -0.66713577  0.4082483

$v

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

[1,] -0.07736219 -0.71960032 -0.18918124

[2,] -0.19033085 -0.50893247  0.42405898

[3,] -0.30329950 -0.29826463 -0.45330031

[4,] -0.41626816 -0.08759679 -0.01637004

[5,] -0.52923682  0.12307105 0.64231130

[6,] -0.64220548  0.33373889 -0.40751869

> A.s$u%*%diag(A.s$d)%*%t(A.s$v)

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

[1,]   1    4    7  10   13   16

[2,]   2    5    8  11   14   17

[3,]   3    6    9  12   15   18

5.與矩陣(陣列)運算有關的函式

(1)取矩陣的維數

函式dim()返回一個矩陣的維數,nrow()返回行數,ncol()返回列數。例:

> A=matrix(1:12,3,4)

> A

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

[1,]   1    4    7  10

[2,]   2    5    8  11

[3,]   3    6    9  12

> dim(A)

[1] 3 4

> nrow(A)

[1] 3

> ncol(A)

[1] 4

(2)矩陣的拉直

as.vector()函式可將矩陣轉化為向量,如:

> A<-matrix(1:6,nrow=2);A

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

[1,]   1    3    5

[2,]   2    4    6

> as.vector(A)

[1] 1 2 3 4 5 6

(3)陣列的維名字

陣列可以有一個屬性dimnames儲存各維的各個下標的名字,預設為NULL。如:

>x<-matrix(1:6,ncol=2,dimnames=list(c("A","B","C"),c("a","b")),byrow=T);x

  a b

A 1 2

B 3 4

C 5 6

也可以先定義矩陣x,然後再為dimnames(x)賦值。如:

> x<-matrix(1:6,ncol=2,byrow=T)

>dimnames(x)<-list(c("A","B","C"),c("a","b"))

對於矩陣,還可以使用屬性rownames和colnames來訪問行名和列名,例如:

> x<-matrix(1:6,ncol=2,byrow=T)

>colnames(x)<-c("a","b")

>rownames(x)<-c("A","B","C")

(4)apply()函式

對於向量,可以用sum、mean等函式對其進行計算,對於陣列(矩陣),如果想對其一維(或若干維)進行某種計算,可用apply()函式,其一般形式為

apply(A,margin,fun)

其中A是一個數組(矩陣),margin是固定哪些維不變,fun是用來計算的函式,如:

> A<-matrix(1:6,nrow=2);A

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

[1,]   1    3    5

[2,]   2    4    6

> apply(A,1,sum)

[1] 9 12

> apply(A,2,mean)

[1] 1.5 3.5 5.5

參考:《統計建模與R軟體》 薛毅 陳立萍 編著 清華大學出版社