矩陣與陣列 在

關於 矩陣與陣列

在處理同質的原子資料集合時（例如全部是數值，或全部是字元），R 提供以下這些結構：

• 用於一維資料的vector，以單一索引存取。
• 用於二維、矩形（而非不規則）資料的matrix。
• 用於一般 N 維情況的array。matrix其實只是 N=2 的array。

矩陣

矩陣本質上就是一個向量，再加上定義形狀所需的行數與列數。

可以用matrix()函式來建立矩陣：

# data already available 
matrix(1:6, nrow=2, ncol=3)
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6

# empty matrix to be filled in later
matrix(nrow=2, ncol=3)
#>      [,1] [,2] [,3]
#> [1,]   NA   NA   NA
#> [2,]   NA   NA   NA

請注意，資料是以列主序儲存的（也就是先往下、再往右）。這和 Fortran、Matlab、Julia 相同，但和 C/C++ 與 NumPy 不同。

用dim()可以取得矩陣的行數與列數，用nrow()或ncol()則可以取得單一維度。使用length()只會回傳底層向量的長度，也就是元素的總數（上面例子中是 6）。

其他函式也會產生矩陣輸出，例如str_split_fixed()。

library(stringr)

str_split_fixed(c("split some words", "and some more"), " ", 3)
#>      [,1]    [,2]   [,3]   
#> [1,] "split" "some" "words"
#> [2,] "and"   "some" "more" 

矩陣索引

矩陣索引和向量索引的用法完全相同，只是必須依序指定行與列。把某個維度留空，就能取得整行或整列。

使用上面同一個矩陣：

m <- matrix(1:6, nrow=2, ncol=3)
m
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6

m[2, 3]    # row 2, col 3
#> [1] 6

m[1:2, 3]  # rows 1:2, col 3
#> [1] 5 6

 m[2, ]    # all of row 2
#> [1] 2 4 6

請注意，只取單一行或單一列的子集時，回傳的會是向量，而不是矩陣。

class(m)
#> [1] "matrix" "array" 
class(m[2,])
#> [1] "integer"

這種dimension reduction並不總是我們想要的，特別是在預期輸入為矩陣的函式管線中。

指定drop=FALSE可以保證矩陣輸入會回傳矩陣。請注意，只有一行或一列的矩陣，型別_並不_等同於向量。

m[2, , drop = FALSE]
#>      [,1] [,2] [,3]
#> [1,]    2    4    6

as.matrix()函式會把向量轉成只有單一列的矩陣。

as.matrix(1:3)
#>      [,1]
#> [1,]    1
#> [2,]    2
#> [3,]    3

和向量一樣，只要指定新值，就能修改矩陣的局部內容。

m[, 2] <- c(10, 20)
m
#>      [,1] [,2] [,3]
#> [1,]    1   10    5
#> [2,]    2   20    6

新增與刪除行或列

處理向量時，用c()函式就足以插入與附加資料。 矩陣的結構比較複雜，因此在指定所需內容時要更加小心。

rbind()與cbind()函式會把新的行或列綁到矩陣上。

m
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6
m2 <- matrix(7:9, nrow = 1, ncol = 3)

rbind(m, m2)
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6
#> [3,]    7    8    9

刪除行或列時，和向量一樣使用負索引。

m
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6

m[, -2]
     #> [,1] [,2]
#> [1,]    1    5
#> [2,]    2    6

篩選與函式

所有向量函式也都能用在矩陣上，並套用到所有元素。不過，它們可能會失去矩陣的結構，並回傳向量。

sqrt(m)
#>          [,1]     [,2]     [,3]
#> [1,] 1.000000 1.732051 2.236068
#> [2,] 1.414214 2.000000 2.449490

sum(m)
#> [1] 21

m[m > 3]
#> [1] 4 5 6

對矩陣來說，colSums()與rowSums()函式也可能很有用。

這種情況下維度縮減是預期中的結果，所以回傳值會是向量。

# m as defined above
rowSums(m)
#> [1]  9 12
colSums(m)
#> [1]  3  7 11

同樣地，也有rowMeans()與colMeans()函式可以計算平均值。

操作矩陣的通用函式

在其他概念中，我們看過 Tidyverse 函式庫可以簡化向量與陣列（purrr）或資料框（dplyr）的運算。

不過，矩陣與陣列大多不在 Tidyverse 的輸入範圍內（至少目前如此）。

apply()函式是 Base R 的一部分，而且仍然是更精細控制矩陣程式設計的最佳工具。

它的參數依序是：

1. 輸入的矩陣。
2. 要操作的維度（1 代表行，2 代表列）。
3. 要套用的函式。
4. 選擇性的，函式所需的額外參數。

Base R 中已經有rowMeans()函式，但我們也可以用apply()搭配mean()來模擬它。

randm <- sample(1:10, 16, replace = TRUE) |> matrix(4, 4)
randm
#>      [,1] [,2] [,3] [,4]
#> [1,]    4    9    8    3
#> [2,]    6    8    5    5
#> [3,]    4    2    1    2
#> [4,]    4    2    5    4

apply(randm, 1, mean)
#> [1] 6.00 6.00 2.25 3.75

rowMeans(randm)
#> [1] 6.00 6.00 2.25 3.75

線性代數

現代 R「顯然」是一門資料科學語言：主要是為了處理表格式資料而設計，也受到統計學家、流行病學家和社會科學家等使用者的歡迎。 這個領域的主要競爭對手是 Python 的科學運算生態系，包括 Pandas 或 Polars 這類套件。

不過，在更早的幾十年裡，R 也曾被視為可能的工程語言，直接與 Matlab 競爭。

或許出乎意料的是，Base R 至今仍然支援完整的線性代數運算，可以作用在向量與矩陣上。 這些功能如今在終端使用者的程式碼中占多少分量並不清楚，但這類技術在_內部_實作統計函式時仍然非常重要。

R 的線性代數至今仍是一些大學課程的教學內容，例如這門計量經濟學課程。 關於數學（而非程式設計）的優質教材，可以看3Blue1Brown。

對角線

方陣（N x N)）A的主對角線，從左上角（A[1,1]）延伸到右下角（A[N,N]）。

diag()函式的行為會因為傳入的引數不同而有所不同。

傳入向量v時，diag(v)會建立一個對角矩陣。

diag(1:3)
#>      [,1] [,2] [,3]
#> [1,]    1    0    0
#> [2,]    0    2    0
#> [3,]    0    0    3

傳入方陣A時，diag(A)會以向量的形式回傳主對角線上的元素。

A
#>      [,1] [,2] [,3]
#> [1,]   10    5    2
#> [2,]    9    3    7
#> [3,]    8   10    4
diag(A)
#> [1] 10  3  4

傳入單一整數k時，diag(k)會建立一個k × k的單位矩陣：這或許和你想的不一樣。

即使是嚴謹的 R 教科書，講到這裡時也常常會加上一句_「真是想不到」_這類的註解。

diag(3)
#>      [,1] [,2] [,3]
#> [1,]    1    0    0
#> [2,]    0    1    0
#> [3,]    0    0    1

轉置與反矩陣

t()函式可以轉置任何矩陣，用來交換行與列。

m
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6
t(m)
#>      [,1] [,2]
#> [1,]    1    2
#> [2,]    3    4
#> [3,]    5    6

對複數矩陣來說，結果_並不_是共軛複數。 使用Conj(t(A))可以取得矩陣A的共軛轉置（伴隨矩陣）。

對更高維度的陣列，aperm()可以置換維度。

要求方陣A的反矩陣，可以使用solve(A)。 我們有些人覺得這個語法……有點出乎意料。

把矩陣乘上它的反矩陣會得到單位矩陣（可能會有極小的捨入誤差）。

A <- matrix(sample(1:10, 9, replace = TRUE), 3, 3)
A
#>      [,1] [,2] [,3]
#> [1,]   10    5    2
#> [2,]    9    3    7
#> [3,]    8   10    4

solve(A)
#>             [,1]        [,2]        [,3]
#> [1,]  0.16666667  0.00000000 -0.08333333
#> [2,] -0.05747126 -0.06896552  0.14942529
#> [3,] -0.18965517  0.17241379  0.04310345

A %*% solve(A) |> round(digits = 10)
     [,1] [,2] [,3]
#> [1,]    1    0    0
#> [2,]    0    1    0
#> [3,]    0    0    1

乘法

一般的算術運算子，例如+、-、*、^，會以逐元素的方式運算矩陣。如果維度不符，就會套用循環的一般規則。

m
#>      [,1] [,2] [,3]
#> [1,]    1    3    5
#> [2,]    2    4    6

m * m
#>      [,1] [,2] [,3]
#> [1,]    1    9   25
#> [2,]    4   16   36

# same calculation, different syntax
m ^ 2
#>      [,1] [,2] [,3]
#> [1,]    1    9   25
#> [2,]    4   16   36

說明其他各種向量與矩陣乘法的細節超出了我們的範圍，不過簡單摘要如下：

• %*%是matmult()的中綴運算子形式，用來計算矩陣乘積。
• %o%是outer()的中綴運算子形式，用來計算外積。
• crossprod()用來計算叉積。