矩阵的赋值操作

什么是矩阵？

维度为 $(𝑚 \times n)$ 的矩阵 $𝐴$ 是一个包含 $(𝑚 \times n)$ 个元素的集合。这些元素按照由 $𝑚$ 行和 $𝑛$ 列组成的矩形结构进行排列，如下所示：

$$\begin{bmatrix} a_{11} & a_{12} & \dots & a_{1n} \\ a_{21} & a_{22} & \dots & a_{2n} \\ \vdots & \vdots & \vdots & \vdots \\ a_{m1} & a_{m2} & \dots & a_{mn} \end{bmatrix}$$

其中，$𝑎_{𝑖𝑗}$ 满足 $𝑖 = 1,...,𝑚$ 且 $𝑗 = 1,...,𝑛$。

在 R 中进行矩阵赋值

在基础 R 中，matrix() 函数通过传入数据和维度来创建矩阵，并可选择按列（默认）或按行填充。

1A <- matrix(c(1, 2, 3, 4, 5, 6), 
2            ncol = 3, 
3            nrow = 2, 
4            byrow = TRUE
5            )
6
7class(A)
8dim(A)
9nrow(A)
10ncol(A)

说明

在上述代码中：

• 第 1–5 行：我们使用 matrix() 函数创建一个维度为 $(2 \times 3)$ 的矩阵。
• 第 7 行：使用 class() 函数返回 R 对象的类属性值。
• 第 8 行：使用 dim() 函数打印矩阵的维度。
• 第 9–10 行：nrow() 和 ncol() 函数分别输出矩阵的行数和列数。

返回的对象是一个维度为 $(2 \times 3)$ 的 matrix 类型对象，即 2 行 3 列。

---

矩阵 A $$\begin{bmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \\ \end{bmatrix}$$

矩阵可视化

根据数据科学专业人士的观点，研究数据的第一步是可视化。矩阵代数也是如此。了解数据所采用的结构非常重要。如果我们使用 plot.matrix 包，可以将矩阵以热图的形式进行可视化。

让我们执行以下代码，显示一个随机分布的热图。

1library(plot.matrix)
2# 定义颜色调色板
3pal <- c("gray100", "gray90", "gray80", "gray70", "gray60",
4         "gray50", "gray40", "gray30", "gray20", "gray10")
5# 生成随机数数据集
6x <- runif(n = 81, min = 0, max = 100)
7X <- matrix(x, ncol = 9)
8# 创建绘图区域
9png(file = "output/graph.png")
10par(mar = c(1, 1, 1, 1) + 3.5)
11# 绘制图形
12plot(X, breaks = 10, col = pal)
13# dev.off()

在上述示例中，矩阵没有结构，因为其元素是从均匀分布中随机生成的。

下面的对角矩阵则采用了对角结构：

1library(plot.matrix)
2pal <- c("gray100", "gray90", "gray80", "gray70", "gray60",
3         "gray50", "gray40", "gray30", "gray20", "gray10")
4
5X <- diag(runif(n = 9, min = 1, max = 10))
6
7png(file = "output/graph.png")
8par(mar = c(1, 1, 1, 1) + 3.5)
9
10plot(X, breaks = 10, col = pal)
11# dev.off()

在上述两个示例代码中：

• 第 2 行：我们定义了热图的颜色范围 pal。
• 第 5 行：runif() 方法用于生成均匀分布的随机数。
• 第 8 行：png() 方法将输出图形保存为 graph.png，plot() 函数用于显示结果。
• 第 9 行：使用 par() 函数设置边距。

在 Rcpp 中进行矩阵赋值

在 Rcpp 中，针对不同数据类型提供了多种矩阵类：

• LogicalMatrix
• IntegerMatrix
• NumericMatrix
• ComplexMatrix
• CharacterMatrix

注意：本课程中我们使用的是 NumericMatrix 类。

---

下一个代码示例与上面的 matrix 函数类似。我们需要做一些额外工作来实现 byrow = TRUE，因为在 Rcpp 中没有原生的按行赋值方法。这种缺失可能是因为 Rcpp 是 R 与 C++ 之间的接口，而 C++ 中的矩阵是按列存储的。（这指的是计算机中矩阵数据的存储方式。）

注意：main.cpp 文件包含赋值操作的算法，main.r 文件实现了该算法。

1Rcpp::sourceCpp('/usercode/main.cpp')
2
3A_cpp_col <- mat_assign(v = c(1, 2, 3, 4, 5, 6), 
4                        rows = 2, 
5                        cols = 3, 
6                        byrow = 0)
7
8A_cpp_row <- mat_assign(v = c(1, 2, 3, 4, 5, 6), 
9                        rows = 2, 
10                        cols = 3, 
11                        byrow = 1)
12
13A_cpp_col
14cat("\n")
15A_cpp_row

现在，可以在相同的 R 使用场景中测试 Rcpp 函数，通过插入矩阵元素。首先，我们按列插入元素。

1A_cpp_col <- mat_assign(v = c(1, 2, 3, 4, 5, 6), rows = 2, cols = 3, byrow = 0)

按列赋值 $$\begin{bmatrix} 1 & 3 & 5 \\ 2 & 4 & 6 \\ \end{bmatrix}$$

---

然后我们按行插入元素。

1A_cpp_row <- mat_assign(v = c(1, 2, 3, 4, 5, 6), rows = 2, cols = 3, byrow = 1)

按行赋值 $$\begin{bmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \\ \end{bmatrix}$$

---

在 Armadillo 中进行矩阵赋值

Armadillo 提供了一个矩阵模板类，但为了简化，本课程中我们将使用 typedef mat = Mat<double>。

Armadillo 提供了一个矩阵构造函数，接受向量和期望的矩阵维度作为输入。此外，Armadillo 没有内置的按行赋值方法，因此需要额外处理。

1Rcpp::sourceCpp('/usercode/main.cpp')
2
3A_arma_col <- mat_assign_A(c(1, 2, 3, 4, 5, 6), 2, 3, 0)
4
5A_arma_row <- mat_assign_A(c(1, 2, 3, 4, 5, 6), 2, 3, 1)
6
7A_arma_col
8cat("\n")
9A_arma_row

现在，可以在相同的 R 使用场景中测试 mat_assign_A() 函数，通过插入矩阵元素。首先测试按列赋值：

1A_arma_col <- mat_assign_A(c(1, 2, 3, 4, 5, 6), 2, 3, 0)

注意：书写矩阵时：

1c(1, 2, 3, 4, 5, 6), 2, 3, 0)

或者：

1c(1, 2, 3, 4, 5, 6), nrow = 2, ncol = 3, bycol = 0)

都是个人选择。两种写法都表示创建一个 $(2 \times 3)$ 的矩阵。

使用 Armadillo 按列赋值 $$\begin{bmatrix} 1 & 3 & 5 \\ 2 & 4 & 6 \\ \end{bmatrix}$$

---

接下来测试按行赋值：

1A_arma_row <- mat_assign_A(c(1, 2, 3, 4, 5, 6), 2, 3, 1)

使用 Armadillo 按行赋值 $$\begin{bmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \\ \end{bmatrix}$$

---

在 Eigen 中进行矩阵赋值

Eigen 提供了一个矩阵模板类，但为了简化，本课程中我们使用 double 数据类型，并将维度设为 Dynamic：

typedef Matrix<double, Dynamic, Dynamic> MatrixXd

Eigen 使用逗号初始化器（<<）按列将向量赋值给矩阵。在下面的赋值函数实现中，使用逐元素赋值来处理按行的情况。

1Rcpp::sourceCpp('/usercode/main.cpp')
2
3A_Eigen_col <- mat_assign_E(c(1, 2, 3, 4, 5, 6), 2, 3, 0)
4
5A_Eigen_row <- mat_assign_E(c(1, 2, 3, 4, 5, 6), 2, 3, 1)
6
7A_Eigen_col
8cat("\n")
9A_Eigen_row

现在，可以在相同的 R 使用场景中测试 Eigen 函数 mat_assign_E()，通过插入矩阵元素。首先按列查看：

1A_Eigen_col <- mat_assign_E(c(1, 2, 3, 4, 5, 6), 2, 3, 0)

使用 Eigen 按列赋值 $$\begin{bmatrix} 1 & 3 & 5 \\ 2 & 4 & 6 \\ \end{bmatrix}$$

---

然后按行再次查看：

1A_Eigen_row <- mat_assign_E(c(1, 2, 3, 4, 5, 6), 2, 3, 1)

使用 Eigen 按行赋值 $$\begin{bmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \\ \end{bmatrix}$$