R 统计与机器学习研讨会（三）
原文：
annas-archive.org/md5/e61f21209b8b5ad343780cefdc55a102译者：飞龙
第六章：使用 R Markdown 进行有效报告
上一章介绍了不同的绘图技术，所有这些都是静态的。在本章中，我们将更进一步，讨论如何使用 R Markdown 一致地生成图表和表格。
到本章结束时，您将学习到 R Markdown 报告的基础知识，包括如何添加、微调和自定义图表和表格以制作交互式和有效的报告。您还将了解如何生成有效的 R Markdown 报告，这可以为您的演示增添色彩。
本章将涵盖以下主题：
-
R Markdown 基础
-
生成财务分析报告
-
自定义 R Markdown 报告
技术要求
要完成本章的练习，您需要拥有以下软件包的最新版本：
-
rmarkdown, 版本 2.17 -
quantmod, 版本 0.4.20 -
lubridate, 版本 1.8.0
请注意，前面提到的软件包版本是在我编写本书时的最新版本。本章的所有代码和数据均可在以下网址找到：github.com/PacktPublishing/The-Statistics-and-Machine-Learning-with-R-Workshop/tree/main/Chapter_6。
R Markdown 基础
R Markdown 是一种格式化语言，可以帮助您有效地动态地从数据中揭示洞察力，并以 PDF、HTML 文件或网络应用程序的形式生成报告。它允许您通过本书前面介绍的各种图表和表格形式整理您的分析，并以一致、整洁、透明的方式展示，便于其他分析师轻松复制。无论是在学术界还是工业界，证明您分析的可重复性是您工作的一个基本品质。当其他人可以轻松复制并理解您在分析中所做的工作时，这会使沟通更加容易，并使您的工作更加值得信赖。由于所有输出都是基于代码的，因此您在展示初步工作并返回进行进一步修改时，可以轻松地微调分析，这是现实数据分析中常见的迭代过程。
使用 R Markdown，您可以将代码及其输出（包括图表和表格）一起展示，并添加周围文本作为上下文。它与使用 Python 的 Jupyter Notebook 类似，但它在 tidyverse 生态系统支持下具有优势。
R Markdown 基于 Markdown 语法，这是一种易于遵循的标记语言，允许用户从纯文本文件创建类似网页的文件。让我们先下载 R Markdown 软件包，并创建一个简单的起始文件。
开始使用 R Markdown
R Markdown 允许我们创建高效的报告来总结我们的分析，并将结果传达给最终用户。要在 RStudio 中启动 R Markdown，我们首先需要下载 rmarkdown 包并将其加载到控制台，这可以通过以下命令完成：
>>> install.packages("rmarkdown")
>>> library(rmarkdown)
R Markdown 有一种专门的文件类型，以 .Rmd 结尾。要创建 R Markdown 文件，我们可以在 RStudio 中选择 文件 | 新建文件 | R Markdown；这将显示 图 6*.1* 中的窗口。左侧面板包含我们可以选择的不同格式，其中 文档 是一组常见的文件类型，如 HTML、PDF 和 Word，演示 以类似 PowerPoint 的演示模式呈现 R Markdown 文件，Shiny 在 R Markdown 文件中添加交互式 Shiny 组件（交互式小部件），从模板 提供了一系列启动模板以加速报告生成：
图 6.1 – 创建 R Markdown 文件
让我们从 my first rmarkdown 开始，并点击 .Rmd 文件，将创建一个包含基本指令的文件。并非所有这些信息都会被使用，因此熟悉常见组件后，您可以自由删除脚本中的不必要代码。
R Markdown 文档由三个组件组成：文件的元数据、报告的文本和分析的代码。我们将在以下各节中查看这些组件。
了解 YAML 标题
如 图 6*.2* 所示，R Markdown 脚本的顶部是一组由两组三个连字符 --- 包裹的元数据标题信息，并包含在 YAML 标题中。YAML，一种人类可读的数据序列化语言，是用于配置文件中分层数据结构的语法。在这种情况下，默认信息包括标题、输出格式和日期，以键值对的形式表示。标题中的信息会影响整个文档。例如，要生成 PDF 文件，我们只需在输出配置中将 html_document 切换为 pdf_document。这是标题中所需的最小信息集，尽管鼓励您添加作者信息（通过 图 6*.2* 中的相同初始窗口）以显示您的工作版权：
图 6.2 – 默认 R Markdown 脚本的 YAML 标题
在设置好标题信息并假设所有额外的代码都已删除后，我们可以通过点击 test.Rmd 来编译 R Markdown 文件为 HTML 文件：
图 6.3 – 使用 Knit 按钮将 R Markdown 文件转换为 HTML 文件
编译 R Markdown 文件将生成一个在单独的预览窗口中打开的 HTML 文件。它还会在同一文件夹中保存一个名为test.html的 HTML 文件。
接下来，我们将学习更多关于 R Markdown 文件主体结构和语法的知识，包括文本格式化和处理代码块。
格式化文本信息
文本信息的重要性与您为分析和建模编写的代码相当，甚至更高。好的代码通常有很好的文档，当您的最终用户是非技术性的时，这一点尤为重要。在适当的位置放置背景信息、假设、上下文和决策过程是您技术分析的重要伴侣，除了分析的透明性和一致性之外。在本节中，我们将回顾我们可以用来格式化文本的常用命令。
练习 6.1 – 在 R Markdown 中格式化文本
在这个练习中，我们将使用 R Markdown 生成图 6.4中显示的文本：
图 6.4 – 使用 R Markdown 生成的 HTML 文件示例文本
文本包括标题、一些斜体或粗体的单词、一个数学表达式和四个无序列表项。让我们看看如何生成这个文本：
-
使用
#符号编写一级标题：# Introduction to statistical model注意，我们使用的井号越多，标题就会越小。请记住，在井号和文本之间添加一个空格。
-
通过将文本包裹在
* *中以实现斜体和$$以实现数学表达式，来编写中间句子：A *statistical model* takes the form $y=f(x)+\epsilon$, where -
通过在每个项目前使用
*并使用** **将文本包裹起来以实现粗体，来生成无序列表：* $x$ is the **input** * $f$ is the **model** * $\epsilon$ is the **random noise** * $y$ is the **output**注意，我们可以轻松地将输出文件从 HTML 切换到 PDF，只需将
output: html_document更改为output: pdf_document。结果输出显示在图 6.5中：
图 6.5 – 使用 R Markdown 生成的 PDF 文件示例文本
将 R Markdown 文件编译成 PDF 文档可能需要您安装额外的包，例如 LaTeX。当出现错误提示说该包不可用时，只需在控制台中安装此包，然后再进行编译。我们还可以使用编译按钮的下拉菜单来选择所需的输出格式。
此外，YAML 标题中日期键的值是一个字符串。如果您想自动显示当前日期，可以将字符串替换为````pyr Sys.Date()"``.
These are some of the common commands that we can use in a .Rmd file to format the texts in the resulting HTML file. Next, we will look at how to write R code in R Markdown.
Writing R code
In R Markdown, the R code is contained inside code chunks enclosed by three backticks, ```` py ，这在 R Markdown 文件中用于将代码与文本分开。代码块还伴随着对应于所使用语言和其他配置的规则和规范，这些规则和规范位于花括号{}内。代码块允许我们渲染基于代码的输出或在报告中显示代码。
下面的代码片段展示了示例代码块，其中我们指定语言类型为 R 并执行赋值操作：
```{r}
a = 1
```py
除了输入代码块的命令外，我们还可以在工具栏中点击代码图标（以字母c开头）并选择 R 语言的选项，如图图 6.6所示。请注意，您还可以使用其他语言，如 Python，从而使 R Markdown 成为一个多功能的工具，允许我们在一个工作文件中使用不同的编程语言：
图 6.6 – 插入 R 代码块
每个代码块都可以通过点击每个代码块右侧的绿色箭头来执行，结果将显示在代码块下方。例如，图 6.7显示了执行赋值和打印变量后的输出：
图 6.7 – 执行代码块
我们还可以在代码块的大括号中指定其他选项。例如，我们可能不希望在生成的 HTML 文件输出中包含特定的代码块。为了隐藏代码本身并只显示代码的输出，我们可以在代码块的相关配置中添加echo=FALSE，如下面的代码块所示：
```{r echo=FALSE}
a = 1
a
```py
图 6.8显示了生成的 HTML 文件中的两种不同类型的输出：
图 6.8 – 在 HTML 文件中显示和隐藏源代码
此外，当我们加载当前会话中的包时，我们可能在控制台得到一个警告消息。在 R Markdown 中，这样的警告消息也会出现在生成的 HTML 中。要隐藏警告消息，我们可以在配置中添加warning=FALSE。例如，在下面的代码片段中，我们在加载dplyr包时隐藏了警告消息：
```{r warning=FALSE}
library(dplyr)
```py
图 6.9比较了加载包时显示或隐藏警告消息的两个场景：
图 6.9 – 加载包时隐藏警告消息
在这些构建块介绍完毕后，我们将在下一节进行案例研究，使用谷歌股票价格数据生成财务分析报告。
生成财务分析报告
在本节中，我们将分析来自 Yahoo! Finance 的谷歌股票数据。为了方便数据下载和分析，我们将使用quantmod包，该包旨在帮助量化交易者开发、测试和部署基于统计的贸易模型。让我们安装这个包并将其加载到控制台：
>>> install.packages("quantmod")
>>> library(quantmod)
接下来，我们将使用 R Markdown 生成 HTML 报告，并介绍数据查询和分析的基础知识。
获取和显示数据
让我们通过一个练习来生成一个初始报告，该报告会自动从 Yahoo! Finance 查询股票数据，并显示数据集中的基本信息。
练习 6.2 – 生成基本报告
在这个练习中，我们将设置一个 R Markdown 文件，下载谷歌的股价数据，并显示数据集的一般信息：
-
创建一个名为
Financial analysis的空 R Markdown 文件，并在 YAML 文件中设置相应的output、date和author：--- title: "Financial analysis" output: html_document date: "2022-10-12" author: "Liu Peng" --- -
创建一个代码块来加载
quantmod包并使用getSymbols()函数查询谷歌的股价数据。将结果数据存储在df中。同时隐藏结果 HTML 文件中的所有消息，并添加必要的文本说明：# Analyzing Google's stock data since 2007 Getting Google's stock data ```{r warning=FALSE, message=FALSE} library(quantmod) df = getSymbols("GOOG", auto.assign=FALSE) ```py在这里，我们指定
warning=FALSE以隐藏加载包时的警告消息，message=FALSE以隐藏调用getSymbols()函数时生成的消息。我们还指定auto.assign=FALSE将结果 DataFrame 分配给df变量。另外，请注意，我们可以在代码块内添加文本作为注释，这些注释将被视为以井号#开头的典型注释。 -
通过三个单独的代码块计算总行数并显示 DataFrame 的前两行和最后两行。为代码添加相应的文本作为文档：
Total number of observations of `df` ```{r} nrow(df) ```py Displaying the first two rows of `df` ```{r} head(df, 2) ```py Displaying the last two rows of `df` ```{r} tail(df, 2) ```py注意，我们使用
` `来表示文本中的内联代码。到目前为止，我们可以编织 R Markdown 文件以观察生成的 HTML 文件，如图图 6.10所示。经常检查输出是一个好习惯，这样就可以及时纠正任何潜在的不期望的错误：
图 6.10 – 显示 HTML 输出
-
使用
chart_Series()函数绘制每日收盘价的时间序列图：Plotting the stock price data ```{r} chart_Series(df$GOOG.Close,name="Google Stock Price") ```py将此代码块添加到 R Markdown 文档中并编织它，将生成与图 6.11中所示相同的输出文件，并增加一个额外的图形。
chart_Series()函数是quantmod提供的用于绘图的实用函数。我们也可以根据前一章讨论的ggplot包来绘制它：
图 6.11 – 自 2017 年以来谷歌的每日股价
除了从代码生成图形外，我们还可以将链接和图片包含在输出中。此图片可以从本地驱动器或从网络加载。在下面的代码片段中，我们添加了一行带有超链接的文本，指向一个示例图片，并在下一行直接从 GitHub 读取显示该图片：
The following image can be accessed [here](https://github.com/PacktPublishing/The-Statistics-and-Machine-Learning-with-R-Workshop/blob/main/Chapter_6/Image.png).
注意，我们通过将单词 here 放在方括号内并跟在括号中的超链接后面来添加一个超链接。要添加图片，我们可以在方括号前添加一个感叹号。我们还可以通过在图片链接后添加 {width=250px} 来指定图片的大小。
在 R Markdown 中编译前面的代码生成 图 6*.12*：
图 6.12 – 从 GitHub 可视化图像
接下来，我们将执行数据分析并将结果以文本形式显示。
执行数据分析
加载数据集后，我们可以在生成的输出文档中执行数据分析并展示洞察力，所有这些都是自动且一致的。例如，我们可以展示特定时期内股票价格的高级统计数据，如平均、最大和最小价格。这些统计数据可以嵌入到文本中，使演示风格更加自然和自包含。
练习 6.3 – 执行简单的数据分析
在这个练习中，我们将提取 Google 的年度最高、平均和最低股价。为此，我们首先将数据集从其原始的 xts 格式转换为 tibble 对象，然后使用 dplyr 总结这些统计量。最后，我们将在 HTML 文档的文本中显示这些信息：
-
加载
dplyr和tibble包，并将df从xts格式转换为tibble格式。将生成的tibble对象存储在df_tbl中：library(dplyr) library(tibble) df_tbl = df %>% as_tibble() %>% add_column(date = index(df), .before = 1)在这里，我们将使用
as_tibble()函数将xts对象转换为tibble格式，然后使用add_column()函数在 DataFrame 的开头插入一个日期列。日期信息作为索引存在于原始的xts对象中。 -
存储自 2022 年以来的年度最高、平均和最低收盘价，分别存储在
max_ytd、avg_ytd和min_ytd中：max_ytd = df_tbl %>% filter(date >= as.Date("2022-01-01")) %>% summarise(price = max(GOOG.Close)) %>% .$price avg_ytd = df_tbl %>% filter(date >= as.Date("2022-01-01")) %>% summarise(price = mean(GOOG.Close)) %>% .$price min_ytd = df_tbl %>% filter(date >= as.Date("2022-01-01")) %>% summarise(price = min(GOOG.Close)) %>% .$price对于每个统计量，我们首先按日期过滤，然后根据
GOOG.Close列提取相关统计量。最后，我们将结果作为单个标量值返回，而不是 DataFrame。 -
以文本形式显示这些统计量：
Google's **highest** year-to-date stock price is `r max_ytd`. Google's **average** year-to-date stock price is `r avg_ytd`. Google's **lowest** year-to-date stock price is `r min_ytd`.如 图 6*.13* 所示，编译文档将统计量输出到 HTML 文件中，这使得我们可以在 HTML 报告中引用代码结果：
图 6.13 – 提取简单统计量并在 HTML 格式中显示
在下一节中，我们将探讨如何在 HTML 报告中添加图表。
在报告中添加图表
在 HTML 报告中添加图表的方式与在 RStudio 控制台中相同。我们只需在代码块中编写绘图代码，然后编译 R Markdown 文件后，图表就会出现在生成的报告中。让我们通过一个练习来可视化使用上一章中介绍的 ggplot2 包的股票价格。
练习 6.4 – 使用 ggplot2 添加图表
在这个练习中，我们将通过线形图可视化过去三年的平均月收盘价。我们还将探索报告中图表的不同配置选项：
-
创建一个包含 2019 年至 2021 年间月度平均收盘价的数据库：
library(ggplot2) library(lubridate) df_tbl = df_tbl %>% mutate(Month = factor(month(date), levels = as.character(1:12)), Year = as.character(year(date))) tmp_df = df_tbl %>% filter(Year %in% c(2019, 2020, 2021)) %>% group_by(Year, Month) %>% summarise(avg_close_price = mean(GOOG.Close)) %>% ungroup()在这里，我们首先创建两个额外的列，分别称为
Month和Year，这些列是基于日期列通过lubridate包中的month()和year()函数派生出来的。我们还把Month列转换为因子类型的列，其级别在 1 到 12 之间，这样当我们在后面绘制月度价格图时，这个列可以遵循特定的顺序。同样，我们将Year列设置为字符类型的列，以确保它不会被ggplot2解释为数值变量。接下来，我们通过
Year对df_tbl变量进行筛选，按Year和Month分组，并计算GOOG.Close的平均值，然后使用ungroup()函数从保存在tmp_df中的结果 DataFrame 中移除分组结构。 -
在线形图中将每年的月度平均收盘价作为单独的线条绘制。更改相应的图表标签和文本大小：
p = ggplot(tmp_df, aes(x = Month, y = avg_close_price, group = Year, color = Year)) + geom_line() + theme(axis.text=element_text(size=16), axis.title=element_text(size=16,face="bold"), legend.text=element_text(size=20)) + labs(titel = "Monthly average closing price between 2019 and 2021", x = "Month of the year", y = "Average closing price") p在代码块中运行前面的命令将生成图 6**.14中显示的输出。请注意，我们还添加了标题和一些文本，以指出代码的目的和上下文。在编织 R Markdown 文件后，代码和输出会自动显示，这使得 R Markdown 成为生成透明、吸引人和可重复的技术报告的绝佳选择：
图 6.14 – 添加图表以显示过去三年的月度平均收盘价
我们还可以配置图表的大小和位置。
-
通过在代码块的配置部分设置
fig.width=5和fig.height=3来缩小图表的大小，并显示输出图形：Control the figure size via the `fig.width` and `fig.height` parameters. ```{r fig.width=5, fig.height=3} p ```py使用这些添加的命令编织文档会产生图 6.15：
图 6.15 – 改变图表的大小
-
将图表的位置对齐，使其位于文档的中心：
Align the figure using the `fig.align` parameter. ```{r fig.width=5, fig.height=3, fig.align='center'} p ```py使用这些添加的命令编织文档会产生图 6.16：
图 6.16 – 改变图表的位置
-
为图表添加标题：
Add figure caption via the `fig.cap` parameter. ```{r fig.width=5, fig.height=3, fig.align='center', fig.cap='图 1.1 2019 年至 2021 年间的月度平均收盘价'} p ```py使用这些添加的命令编织文档会产生图 6.17：
图 6.17 – 为图表添加标题
除了图形外，表格也是报告中常用的一种用于呈现和总结信息的方式。我们将在下一节中探讨如何生成表格。
向报告中添加表格
当报告用户对深入了解细节或进一步分析感兴趣时，以表格形式呈现信息是图形对应物的良好补充。对于最终用户来说，能够访问和使用报告中的数据起着关键作用，因为这给了他们更多控制权，可以控制报告中已经预处理好的信息。换句话说，基于 R Markdown 的 HTML 报告不仅以图形形式总结信息以便于消化，还提供了关于特定数据源的详细信息作为表格，以促进即席分析。
我们可以使用 knitr 包中的 kable() 函数添加表格，该函数是支持在每个代码块中执行代码的核心引擎，然后在对 R Markdown 文档进行编织时进行动态报告生成。请注意，在通过 kable() 将数据作为表格展示之前进行预处理和清理数据是一个好的实践；这项任务应该只涉及展示一个干净且有序的表格。
让我们通过一个练习来看看如何向报告中添加干净的表格。
练习 6.5 – 使用 kable() 添加表格
在这个练习中，我们将以表格形式展示 tmp_df 变量的前五行，然后演示不同的表格显示配置选项：
-
使用
knitr包中的kable()函数显示tmp_df的前五行：# Adding tables Printing `tmp_df` as a static summary table via the `kable()` function. ```{r} library(knitr) kable(tmp_df[1:5,]) ```py使用这些添加的命令编织文档会产生 图 6*.18*：
图 6.18 – 向报告中添加表格
-
使用
col.names参数更改表格的列名：Changing column names via the `col.names` parameter. ```{r} kable(tmp_df[1:5,], col.names=c("Year", "Month", "Average closing price")) ```py使用这些添加的命令编织文档会产生 图 6*.19*：
图 6.19 – 更改表格中的列名
我们还可以使用 align 参数修改表格内的列对齐方式。默认情况下，数值列的列对齐在右侧，其他所有类型的列对齐在左侧。如 图 6*.19* 所示，Year（字符类型）和 Month（因子类型）列左对齐，而 Average closing price（数值）列右对齐。对齐方式按列指定，使用单个字母表示，其中 "l" 表示左对齐，"c" 表示居中对齐，"r" 表示右对齐。
-
使用
align参数将所有列对齐到中心：Align the table via the `align` argument. ```{r} kable(tmp_df[1:5,], col.names=c("Year", "Month", "Average closing price"), align="ccc") ```py在这里，我们指定
align="ccc"以将所有列对齐到中心。使用这些添加的命令编织文档会产生 图 6*.20*：
图 6.20 – 使表格的所有列居中
最后，我们还可以为表格添加一个标题。
- ```py Add table caption via the `caption` parameter. ```{r} date: "2022-10-12" ```py Knitting the document with these added commands produces *Figure 6**.21*:
Figure 6.21 – Adding a caption to the table
In the next section, we will discuss some common options we can use to modify the code chunk outputs after knitting the R Markdown document.
Configuring code chunks
We have seen several options from previous exercises that we can use to control the output style of a code chunk. For example, by setting warning=FALSE and message=FALSE, we could hide potential warnings and messages in the resulting output document.
There are other commonly used options. For example, we can use the include option to decide whether the code and results appear in the output report or not. In other words, setting include=FALSE will hide the code and results of the specific code chunk in the report, although the code will still be executed upon knitting the R Markdown document. By default, we have include=TRUE and all the code and execution results will appear in the report.
Another related option is echo, where setting echo=FALSE hides the code and only shows the execution outputs in the report. We can consider this option when we’re generating plots in the report since most users are more interested in the graphical analysis compared to the process that generates the graph. Again, by default, we have echo=TRUE, which displays the code in the report before the plots.
Besides this, we may only be interested in showing some code instead of executing all of it. In this case, we can set eval=FALSE to make sure that the code in the code chunk does not impact the overall execution and result of the report. This is in contrast to setting include=FALSE, which hides the code but still executes it in the backend, thus bearing an effect on the subsequent code. By default, we have eval=FALSE, which evaluates all the code in the code chunk. Figure 6*.22* summarizes these three options:
| | Code execution | Code appearance | Result appearance |
| include=FALSE | Yes | No | No |
| echo=FALSE | Yes | No | Yes |
| eval=FALSE | No | Yes | No |
Figure 6.22 – Common options for configuring code chunks
Next, we will go over an exercise to practice different options.
Exercise 6.6 – configuring code chunks
In this exercise, we will go through a few ways to configure the code chunks we covered previously:
-
Display the maximum closing price for the past five years in a table. Show both the code and the result in the report:
# author: "刘鹏" html_document: ```py{r} tmp_df = df_tbl %>% mutate(Year = as.integer(Year)) %>% filter(Year >= max(Year)-5, Year < max(Year)) %>% group_by(Year) %>% summarise(max_closing = max(GOOG.Close)) kable(tmp_df)Here, we first convert `Year` into an integer-typed variable, then subset the DataFrame to keep only the last five years of data, followed by extracting the maximum closing price for each year. The result is then shown via the `kable()` function. Knitting the document with these added commands produces *Figure 6**.23*. The result shows that Google has been making new highs over the years:
Figure 6.23 – Displaying the maximum closing price for the past five years
-
Obtain the highest closing price in a code chunk with the code and result hidden in the report by setting
include=FALSE. Display the result in a new code chunk:执行代码块，但通过设置`include=FALSE`在输出中隐藏代码和结果。 ```py{r include=FALSE} total_max_price = max(df_tbl$GOOG.Close)author: "刘鹏"
total_max_priceKnitting the document with these commands produces *Figure 6**.24*:
Figure 6.24 – Hiding the code and result in one code chunk and displaying the result separately
-
For the running table, hide the code chunk and only display the result in the report by setting
echo=FALSE:--- ```py{r echo=FALSE} kable(tmp_df)Knitting the document with these commands produces *Figure 6**.25*:
Figure 6.25 – Hiding the code chunk and displaying the result in the report
-
Only display the code on table generation in the code chunk and do not execute it in the report by setting
eval=FALSE:不执行代码块，并通过设置`eval=FALSE`在输出中仅显示代码。 ```py{r eval=FALSE} kable(tmp_df)Knitting the document with these commands produces *Figure 6**.26*:
Figure 6.26 – Displaying the code chunk without executing it in the report
-
Print a test message and a warning message in separate blocks. Then, put the same contents in a single block by setting
collapse=TRUE:使用这些命令编织文档会产生*图 6.37**。37*： ```py{r} print("This is a test message") warning("This is a test message")在本章中，我们介绍了 R Markdown，这是一个灵活、透明且一致的报告生成工具。我们首先回顾了 R Markdown 的基础知识，包括 YAML 标题和代码块等基本构建块，然后介绍了文本格式化技巧。
print("This is a test message") warning("This is a test message")Knitting the document with these commands produces *Figure 6**.27*, which shows that both the printed and warning messages are shown in a single block together with the code:
Figure 6.27 – Displaying the code and results in one block
In addition, we can hide the warning by configuring the warning attribute in the code chunk.
-
Hide the warning by setting
warning=FALSE:body { ```py{r collapse=TRUE, warning=FALSE} print("This is a test message") warning("This is a test message")Knitting the document with these commands produces *Figure 6**.28*, where the warning has now been removed from the report:
Figure 6.28 – Hiding the warning in the report
Setting the display parameters for each code chunk becomes troublesome when we need to repeat the same operation for many chunks. Instead, we can make a global configuration that applies to all chunks in the R Markdown document by using the knitr::opts_chunk() function at the beginning of the document. For example, the following code snippet hides the warnings for all following code chunks:
显示代码和结果。
```py{r include=FALSE}
knitr::opts_chunk$set(warning=FALSE)
In the next section, we will look at how to customize R Markdown reports, such as by adding a table of contents and changing the report style.
# Customizing R Markdown reports
In this section, we will look at adding metadata such as a table of contents to the report, followed by introducing more options for changing the report style.
## Adding a table of contents
When reading a report for the first time, a table of contents provides an overview of the report and thus helps readers quickly navigate the different sections of the report.
To add a table of contents, we can append a colon to the `html_document` field in the YAML header and set `toc: true` as a separate line with one more indentation than the `html_document` field. This is shown in the following code snippet:
}
代码块选项
font-size: 20px;
执行代码块，并通过设置echo=FALSE在输出中仅显示结果。
date: "2022-10-12"
设置全局选项。
Knitting the document with these commands produces *Figure 6**.29*, where a table of contents is now displayed at the top of the report. Note that when you click on a header in the table of contents, the report will directly jump to that section, which is a nice and user-friendly feature:
Figure 6.29 – Adding a table of contents to the report
We can also set `toc_float=true` to make the table of contents float. With this property specified, the table of contents will remain visible as the user scrolls through the report. The following code snippet includes this property in the YAML header:
body {
}
background-color: #F5F5F5;
params:
通过设置warning=FALSE隐藏警告。
toc_float: true
date: "2022-10-12"
output:
pre {
Knitting the document with these commands produces *Figure 6**.30*, where the table of contents appears on the left-hand side and remains visible as the user navigates different sections:
Figure 6.30 – Setting up a floating table of contents in the report
Next, we will look at creating a report with parameters in the YAML header.
## Creating a report with parameters
Recall that our running dataset contains the daily stock prices of Google since 2007\. Imagine that we need to create a separate annual report for each year; we may need to manually edit the `year` parameter for each report, which would be a repetitive process. Instead, we can set an input parameter in the YAML header as a global variable that’s accessible to all code chunks. When generating other similar reports, we could simply change this parameter and rerun the same R Markdown file.
We can set a parameter input by adding the `params` field, followed by a colon in the YAML header. Then, we must add another line, indent it, and add the key and value of the parameter setting, which are separated by a colon. Note that the value of the parameter is not wrapped in quotations.
Let’s go through an exercise to illustrate this.
### Exercise 6.7 – generating reports using parameters
In this exercise, we will configure parameters to generate reports that are similar and only differ in the parameter setting:
1. Add a `year` parameter to the YAML header and set its value to `2020`:
```
---
html_document:
background-color: #F5F5F5;
html_document:
toc: true
toc_float: true
date: "2022-10-12"
}
params:
year: 2020
---
```py
Here, we use the `params` field to initiate the parameter setting and add `year: 2020` as a key-value pair.
2. Extract the summary of the closing price using the `summary()` function for 2020:
```
# 使用参数生成报告
```
```py{r}
df_tbl %>%
filter(Year == params$year) %>%
select(GOOG.Close) %>%
summary()
```
```py
Knitting the document with these commands produces *Figure 6**.31*, where we use `Year` `== params$year` as a filtering condition in the `filter()` function:
Figure 6.31 – Generating summary statistics of the closing price for 2020 using parameters
1. Change the parameter setting and generate the same statistics for 2021:
```
---
---
output:
color: blue;
toc: true
默认情况下，所有结果都在单独的块中。
title: "财务分析"
author: "刘鹏"
author: "刘鹏"
year: 2021
}
```py
Knitting the document with these commands produces *Figure 6**.32*. With a simple change of value in the parameters, we can generate a report for a different year without editing the contents following the YAML header:
Figure 6.32 – Generating summary statistics of the closing price for 2021 using parameters
We can also create a report based on multiple parameters, which can be appended as key-value pairs in the YAML header.
1. Generate the same statistics for the closing price for Q1 2021:
```
年份`r params$year`和季度`r params$quarter`的摘要统计
```py{r}
df_tbl %>%
mutate(Qter = quarters(date)) %>%
filter(Year == params$year,
Qter == params$quarter) %>%
select(GOOG.Close) %>%
summary()
```
```py
Here, we create a new column to represent the quarter using the `quarters()` function based on the date, followed by filtering using the `year` and `the` `quarter` parameters set in the YAML header. Knitting the document with these commands produces *Figure 6**.33*:
Figure 6.33 – Generating summary statistics of the closing price for 2021 Q1 using multiple parameters
In the following section, we will look at the style of the report using **Cascading Style Sheets** (**CSS**), a commonly used web programming language to adjust the style of web pages.
## Customizing the report style
The report style includes details such as the color and font of text in the report. Like any web programming framework, R Markdown offers controls that allow attentive users to make granular adjustments to the report’s details. The adjustable components include most HTML elements in the report, such as the title, body text, code, and more. Let’s go through an exercise to learn about different types of style control.
### Exercise 6.8 – customizing the report style
In this exercise, we will customize the report style by adding relevant configurations within the `<style>` and `</style>` flags. The specification starts by choosing the element(s) to be configured, such as the main body (using the `body` identifier) or code chunk (using the `pre` identifier). Each property should start with a new line, have the same level of indentation, and have one more level of indentation than the preceding HTML element.
In addition, all contents to be specified are key-value pairs that end with a semicolon and are wrapped within curly braces, `{}`. The style configuration can also exist anywhere after the YAML header. Let’s look at a few examples of specifying the report style:
1. Change the color of the text in the main body to red and the background color to `#F5F5F5`, the hex code that corresponds to gray:
```
# 自定义报告样式
<style>
body {
---
color: orange;
图 6.37 – 改变报告标题的颜色、字体大小和透明度
</style>
```py
Here, we directly use the word `blue` to set the color attribute of the text in the body and the hex code to set its background color; these two approaches are equivalent. Knitting the document with these commands produces *Figure 6**.34*:
Figure 6.34 – Changing the color of the text in the body of the report
1. Change the color of the code in the code chunks to red by specifying `color: red` in the `pre` attribute:
```
# 自定义报告样式
<style>
默认显示代码和结果。
color: blue;
title: "财务分析"
color: green;
使用`caption`参数为表格添加标题：
kable(tmp_df[1:5,], col.names=c("Year", "Month", "Average closing price"), align="ccc", caption="表 1.1 平均收盘价")
color: red;
</style>
```py
Knitting the document with these commands produces *Figure 6**.35*:
Figure 6.35 – Changing the color of the code in the report
1. For the table of contents, change the color of the text and border to `green`, and set the font size to `16px`:
```
# 自定义报告样式
<style>
body {
color: blue;
}
}
border-color: green;
color: red;
background-color: #F5F5F5;
#TOC {
pre {
font-size: 16px;
border-color: green;
年份`r params$year`的摘要统计
output:
```py
Note that the style for the table of contents is specified using `#TOC` without any space in between. Knitting the document with these commands produces *Figure 6**.36*:
Figure 6.36 – Changing the color and font size of the table of contents in the report
1. For the header, change the color to `orange`, the opacity to `0.9`, and the font size to `20px`:
```
# 自定义报告样式
title: "财务分析"
</style>
color: blue;
background-color: #F5F5F5;
}
output:
color: red;
}
html_document:
color: green;
font-size: 16px;
toc: true
toc: true
#header {
toc_float: true
opacity: 0.8;
通过设置`collapse=TRUE`将所有结果合并到一个块中。
}
</style>
}
注意，标题的样式是通过使用`#header`来指定的，其中没有任何空格。
#TOC {
这项练习到此结束。现在，让我们总结本章内容。
# 摘要
title: "财务分析"
接下来，我们通过使用谷歌的股票数据进行了案例研究。在从网络下载股票数据后，我们生成了一份报告来总结每日收盘价的统计数据，向报告中添加了图表和表格，进行了数据处理，并使用不同的样式选项展示了结果。我们还探索了几种配置代码块的不同方法。
最后，我们讨论了如何自定义 R Markdown 报告。我们涵盖的主题包括在报告中添加目录，使用 YAML 标头中的参数创建重复的报告，以及通过编辑不同组件的视觉属性来更改报告的视觉风格，使用 CSS 进行编辑。
在下一章中，我们将开始本书的**第二部分**，并使用 R 语言介绍线性代数和微积分的基础知识。
# 第二部分：R 中线性代数和微积分的基础
随着我们进入本书的第二部分，我们的目标是加深你对支撑数据科学领域的数学基础的理解。虽然第一部分通过统计学基础和 R 编程基础奠定了基础，但这一部分通过深入探讨线性代数和微积分的更复杂方面，并通过 R 实现来展示，从而提升你的技能集。
到本部分结束时，你将全面理解对高级数据科学研究至关重要的数学基础。本书的这一部分将为你提供参与统计学和机器学习更高级主题所需的先进数学和计算工具。
本部分包含以下章节：
+ *第七章*, *R 中的线性代数*
+ *第八章*, *R 中的中级线性代数*
+ *第九章*, *R 中的微积分*
# 第七章：R 中的线性代数
上一章介绍了使用 R Markdown 的高效有效的报告方法。本书的*第一部分*主要涵盖了使用 R 完成实际工作的各个方面。本书的*第二部分*回归到基础，涵盖了数学的两个基本支柱：线性代数和微积分。理解这些基础知识将更好地为我们欣赏和操作常见的数学运算做好准备，使这些运算对我们来说变得自然。*第二部分*的目标是通过本章对线性代数的根本性回顾，开始培养这种水平的素养。
到本章结束时，你将了解线性代数的基本概念，包括向量、矩阵和方程组系统。你还将能够解释线性代数中的基本符号，并使用 R 处理常见的矩阵。
在本章中，我们将涵盖以下主题：
+ 介绍线性代数
+ 常见矩阵运算和属性
+ 解线性方程组
# 技术要求
本章的所有代码和数据均可在[`github.com/PacktPublishing/The-Statistics-and-Machine-Learning-with-R-Workshop/blob/main/Chapter_7/working.R`](https://github.com/PacktPublishing/The-Statistics-and-Machine-Learning-with-R-Workshop/blob/main/Chapter_7/working.R)找到。
# 介绍线性代数
本章深入探讨了数学最重要的分支之一：**线性代数**。线性代数处理数学对象的线性运算，包括向量、矩阵和张量（高维矩阵），这些是最常见的数据形式。例如，我们通常在 Excel 中用来存储数据的典型表格由一系列列组成。每一列被称为一个向量，它存储一定数量的元素，默认情况下以列的形式存在而不是行。这些列向量的集合形成一个矩阵，一个二维的 Excel 表格或 DataFrame，正如我们在前几章中所称呼的。我们也可以将同一表格视为行向量的集合，其中每个向量以行的形式存在。
让我们将这些内容置于上下文中。以下代码片段加载了睡眠数据集，并打印出前六行和三列。我们在以下阐述中使用 A 表示这个 6x3 矩阵：
```py
>>> data(sleep)
>>> head(sleep)
extra group ID
1 0.7 1 1
2 -1.6 1 2
3 -0.2 1 3
4 -1.2 1 4
5 -0.1 1 5
6 3.4 1 6
让我们在这里引入一些符号。我们使用粗体小写字母 a 来表示一个列向量，例如extra变量，它由六个元素组成。换句话说，a 是一个六维向量。列向量 a 可以转置成一个行向量，表示为 a T。由于 a T 是一个行向量，我们可以写成 a T = [0.7, − 1.6, − 0.2, − 1.2, − 0.1, 3.4]，这也是许多书中书写行向量（或转置的列向量）的典型方式。同样，一个行向量可以转置回（原始的）列向量，给我们(a T) T = a。
一个图形说明将有助于理解。图 7*.1*展示了将列向量 a 转换为行向量 a T，然后再将其转换回原始列向量(a T) T = a 的过程。因此，我们可以将矩阵 A 视为水平连接的列向量集合，或者垂直连接的行向量集合：
图 7.1 – 将列向量转换为行向量，然后再将其转换回原始列向量
注意，行也被称为观测值，而列被称为特征或属性。
在下一节中，我们将从处理向量开始。
处理向量
R 中的向量将相同类型的数据元素连接在一起。多个向量携手，通常是并排，形成一个矩阵。因此，它是线性代数的基石，我们必须从这里开始。
创建向量的方法有多种。在接下来的练习中，我们将探索其中的一些选项。
练习 7.1 – 创建向量
这个练习将介绍创建向量的三种常见方法。让我们开始吧：
-
使用
c()函数从 1 到 6 创建一个整数向量，其中c代表连接：x = c(1, 2, 3, 4, 5, 6) >>> x [1] 1 2 3 4 5 6 x = c(1:6)我们还可以使用
seq()函数来创建一个整数序列，该函数根据by参数生成等间距的数字序列。 -
使用
seq()创建相同的向量：x = seq(1, 6, by=1)这三种实现方式创建的是相同的向量。
此外，我们还可以使用
rep()函数来创建重复数字的列表。以下代码片段生成一个包含所有元素都是 1 的 6 元素向量：y = rep(1, 6) >>> y [1] 1 1 1 1 1 1
R 中的向量是可变的，这意味着我们可以更改向量中特定位置的数据。让我们看看如何实现这一点。
练习 7.2 – 修改向量
修改向量的内容是一个常见的操作，尤其是在使用列过滤条件处理 DataFrame 的上下文中。向量的代数主要遵循与标量相同的原理。为了看到这一点，我们必须首先在向量中索引所需的元素，然后使用赋值操作来覆盖特定元素的值。让我们看看：
-
将向量
x中第二个元素的值改为20：x[2] = 20 >>> x [1] 1 20 3 4 5 6因此，我们可以通过将位置包裹在方括号中来访问向量中的元素。我们还可以执行对向量中所有元素产生相同效果的批量操作。
-
将向量
x中的所有元素都乘以二：>>> x*2 [1] 2 40 6 8 10 12在这里，乘数
2被广播到向量的每个元素上，以执行相应的乘法。相同的广播机制也适用于加法操作等。 -
将
x中的所有元素增加一：>>> x + 1 [1] 2 21 4 5 6 7当然，我们可以通过将
x加到y上，一个长度相同的全一向量，来获得相同的结果：>>> x + y [1] 2 21 4 5 6 7
现在我们已经回顾了向量操作的基础，我们可以进入矩阵的领域。
矩阵操作
矩阵是一系列向量叠加而成的。一个 m x n 矩阵可以被认为是水平堆叠的 n 个 m 维列向量，或者垂直堆叠的 m 个 n 维行向量。对于我们所处理的每个 DataFrame，每一行是一个观察值，每一列代表一个特征。
让我们通过完成一个练习来看看如何创建矩阵。
练习 7.3 – 创建矩阵
我们可以使用matrix()函数来创建矩阵。这个函数接受三个参数：通过data传递到矩阵中的数据，通过nrow传递的行数，以及通过ncol传递的列数：
-
创建一个填充数字
2的 3x2 矩阵X：X = matrix(data=2, nrow=3, ncol=2) >>> X [,1] [,2] [1,] 2 2 [2,] 2 2 [3,] 2 2注意这里使用了广播机制。因为我们只传递了一个数字，它被复制到矩阵的所有单元格中。同时，观察结果中的索引模式。逗号左侧的位置索引每个行，逗号右侧的位置索引每个列。我们可以通过索引第二行第一列的元素来验证这个观察：
>>> X[1,2] 2让我们验证这个对象的类别：
>>> class(X) [1] "matrix" "array"结果显示
X既是矩阵也是（多维）数组。我们也可以基于向量创建矩阵。
-
使用变量
x创建一个 3x2 矩阵。按行填充矩阵：>>> matrix(x, nrow=3, ncol=2, byrow=TRUE) [,1] [,2] [1,] 1 20 [2,] 3 4 [3,] 5 6在这里，按行填充意味着按顺序将原始向量
x中的每个元素按行放置，从左到右的第一行开始，当达到当前行的末尾时跳到下一行。我们也可以设计按列填充。 -
改为按列填充矩阵：
>>> matrix(x, nrow=3, ncol=2, byrow=FALSE) [,1] [,2] [1,] 1 4 [2,] 20 5 [3,] 3 6与向量一样，我们可以通过索引定位矩阵中的特定元素，并分配新值来更改矩阵中的元素，如下所示：
X[1,1] = 10 >>> X [,1] [,2] [1,] 10 2 [2,] 2 2 [3,] 2 2
线性代数中的一种常见操作是将矩阵与向量相乘。这种操作产生了许多关于矩阵操作的有趣和重要的解释。我们将在下一节中查看这个操作。
矩阵向量乘法
当将矩阵与向量相乘时，一个总的规则是内维度的相等。换句话说，矩阵的列维度（即列数）乘以右侧的向量时，必须等于乘法向量的行维度。例如，给定一个 m x n 矩阵，它只能乘以一个大小为 n x 1 的列向量，这样的乘法结果将是一个大小为 m x 1 的向量。
在这里，让我们看一个具体的例子。回想一下，X 是一个 3x2 的矩阵。将其乘以一个 2x1 的向量 y，应该产生一个 3x1 的向量 z。从原始向量 y 到新向量 z 的转变在这里有额外的意义：y 改变了空间，现在生活在一个三维世界中而不是二维！我们也可以说，旧的向量 y 由于投影矩阵 X 而被投影和拉伸到新的向量 z。这种投影或变换构成了现代神经网络中大多数操作。通过在不同层之间投影矩阵（也称为表示），神经网络可以学习不同层次的概念抽象。
假设 y 是一个全为 1 的向量。矩阵-向量乘法的规则表明，结果向量 z 中的每个元素都是对应向量的内积。两个向量之间的内积，也称为点积，是每个向量中对应元素乘积的总和。例如，z 中的第一个元素是 12。它的位置索引[1,1]表示它需要 X 矩阵的第一行[10, 2]和 y 矩阵的第一列[1, 1]来进行内积操作。将每个向量中对应元素相乘并求和，我们得到 101+21=12。同样，第二个元素的计算为 21+21=4。图 7.2总结了这种矩阵-向量乘法：
图 7.2 – 矩阵-向量乘法过程示意图
现在，让我们看看如何在 R 中执行矩阵-向量乘法。
练习 7.4 – 应用矩阵-向量乘法
在 R 中，矩阵-向量乘法操作是通过%*%符号实现的。%*%符号与单独的*符号非常不同，后者代表逐元素乘法。以下练习说明了这种差异：
-
将之前的矩阵 X 乘以一个 2x1 的全为 1 的向量：
>>> X %*% c(1,1) [,1] [1,] 12 [2,] 4 [3,] 4结果显示生成的向量是 3x1。请注意，如果维度不匹配，乘法将无法进行：
>>> X %*% c(1,1,1) Error in X %*% c(1, 1, 1) : non-conformable arguments错误信息表明向量
c(1,1,1)的维度与X矩阵中的维度不匹配。我们可以验证结果向量中每个单元格的计算。
-
分别计算结果向量中的第一个元素：
>>> X[1,] %*% c(1,1) [,1] [1,] 12结果显示，第一个元素是第一行
X[1,]和第一列（也是唯一的一列）c(1,1)的点积。我们也可以将点积重新表达为显式的逐元素乘法和求和。 -
将之前的点积重新表达为显式的逐元素乘法和求和：
>>> sum(X[1,] * c(1,1)) [1] 12注意，两个结果返回相同的值，但假设不同的数据结构：点积操作返回一个矩阵，而重新表达的操作返回一个向量。
-
让我们使用结果向量的第二个元素来验证相同的计算：
>>> X[2,] %*% c(1,1) [,1] [1,] 4 >>> sum(X[2,] * c(1,1)) [1] 4我们还可以检查矩阵与标量相乘的行为。
-
将矩阵 X 中的每个元素加倍：
>>> X * 2 [,1] [,2] [1,] 20 4 [2,] 4 4 [3,] 4 4注意，当切换到逐元素乘法符号
*时，结果矩阵的维度保持不变。此外，这里还涉及广播机制，其中 2 这个标量乘以矩阵 X 中的每个元素。效果就像我们用另一个相同维度的矩阵乘以 X 一样：>>> X * matrix(2, 3, 2) [,1] [,2] [1,] 20 4 [2,] 4 4 [3,] 4 4
现在，让我们来讨论矩阵-矩阵乘法。我们将简称为矩阵乘法。
矩阵乘法
矩阵乘法在许多领域中被广泛使用。以神经网络模型为例。图 7.3展示了简单网络架构，称为全连接神经网络，其中所有神经元都是完全连接的。表示为 10 行 3 列的输入数据将进入一系列矩阵乘法（加上这里忽略的非线性变换）来学习有用的表示（Z1）和，因此，准确的预测（Z2）。这里有两个矩阵乘法。第一个矩阵乘法发生在 10x3 的输入数据 X 和权重矩阵 W1 之间，在隐藏层中产生一个 10x4 的隐藏表示 Z1。第二个矩阵乘法使用另一个 4x2 的权重矩阵 W2 将 Z1 转换为最终的 10x2 输出 Z2。我们也可以将一系列矩阵乘法解释为将输入数据 X 转换/投影到隐藏表示 Z1，然后到输出 Z2：
图 7.3 – 简单两层全连接神经网络的示意图
如图 7.3所示，输入数据，由 10 行 3 列的输入层表示，将进入一系列矩阵乘法（加上这里忽略的非线性变换）来学习有用的表示（Z1）和，因此，准确的预测（Z2）。第一个矩阵乘法发生在 10x3 的输入数据 X 和权重矩阵 W1 之间，产生一个 10x4 的隐藏表示 Z1。第二个矩阵乘法使用另一个 4x2 的权重矩阵 W2 将 Z1 转换为最终的 10x2 输出 Z2。
这些矩阵乘法不仅仅是这样。假设输入数据具有许多特征。通过应用这些矩阵乘法与自动学习的权重矩阵，这些特征可以相应地加权，以在输出层产生准确的预测。自动特征学习，包括隐藏层中的那些节点，是现代神经网络与传统手动特征学习相比的一个显著特点。
让我们更详细地看看单个矩阵乘法。图 7.4说明了将 2x2 矩阵[1 3 2 4]与另一个 2x2 矩阵[2 2 2 2]相乘的计算过程。我们还可以将结果矩阵中的每一列视为前面提到的矩阵-向量乘法的结果；然后这些列被连接起来形成新的输出矩阵：
图 7.4 – 矩阵乘法过程的分解
让我们通过一个练习来将实际部分的内容置于上下文中。
练习 7.5 – 矩阵乘法操作
在 R 中，矩阵乘法仍然依赖于%*%符号。在这个练习中，我们将重新生成图 7.4中的示例，在乘法顺序上略有扩展：
-
重新生成之前的矩阵乘法示例：
>>> matrix(1:4, 2, 2) %*% matrix(2, 2, 2) [,1] [,2] [1,] 8 8 [2,] 12 12在这里，我们使用
matrix()函数创建了两个矩阵，第一个是通过转换向量生成的，第二个是通过广播重复标量值生成的。在矩阵代数中，乘法的顺序至关重要；改变顺序在大多数情况下会产生不同的结果。让我交换这两个矩阵。
-
交换这两个矩阵的乘法顺序：
>>> matrix(2, 2, 2) %*% matrix(1:4, 2, 2) [,1] [,2] [1,] 6 14 [2,] 6 14请验证矩阵乘法的结果，并欣赏乘法顺序在两个矩阵乘法中的重要性：矩阵乘法不是交换的。
此外，元素级联乘法是通过
*运算符进行的。 -
将左侧前一个矩阵中的每个元素都乘以 2：
>>> matrix(1:4, 2, 2) * matrix(2, 2, 2) [,1] [,2] [1,] 2 6 [2,] 4 8
几种特殊的矩阵特别引人关注。我们将在下一节中回顾它们。
单位矩阵
矩阵有很多特殊类型。第一种特殊矩阵是单位矩阵，其对角线位置上的值为 1，其他位置为 0。单位矩阵的最大特点是保持原始乘法矩阵的恒等性——也就是说，任何能够成功与单位矩阵相乘的矩阵都会得到与自身相同的结果。这听起来像是乘以 1，而当我们处理 1x1 的单位矩阵时，我们确实是这样做的。
练习 7.6 – 单位矩阵操作
让我们通过一个练习来看看单位矩阵是如何工作的：
-
使用
diag()函数创建一个 2x2 的单位矩阵：>>> diag(2) [,1] [,2] [1,] 1 0 [2,] 0 1要创建一个单位矩阵，我们只需要传递矩阵对角线上的 1 的数量。对角矩阵也可以从一个向量中创建：
>>> diag(c(1,2,3), nrow=3, ncol=3) [,1] [,2] [,3] [1,] 1 0 0 [2,] 0 2 0 [3,] 0 0 3在这里，向量中的元素用于填充对角矩阵的对角线位置，而将所有其他单元格保留为 0。
让我们将 2x2 的单位矩阵与之前的运行矩阵相乘。
-
将这个单位矩阵与之前编号为 1 到 4 的 2x2 矩阵相乘：
>>> matrix(1:4, 2, 2) %*% diag(2) [,1] [,2] [1,] 1 3 [2,] 2 4我们可以看到结果矩阵没有发生变化。让我们再次验证这一点：
>>> matrix(1:4, 2, 2) [,1] [,2] [1,] 1 3 [2,] 2 4最后，让我们验证在交换乘法顺序后结果是否仍然相同。
-
将乘法单位矩阵移到左边并执行矩阵乘法：
>>> diag(2) %*% matrix(1:4, 2, 2) [,1] [,2] [1,] 1 3 [2,] 2 4结果显示原始矩阵没有变化。如果我们用 X 表示原始矩阵，用 I 表示单位矩阵，我们将得到 XI = IX = X。
其他操作可以根据原始矩阵派生出新的矩阵，例如转置或求逆。
转置一个矩阵
我们已经有一些转置向量的经验，这是转置矩阵的特殊情况。转置矩阵意味着翻转原始矩阵 X，生成一个新的矩阵 Xᵀ，其列和行现在是交换的。转置一个转置矩阵 Xᵀ 会返回原始矩阵，(Xᵀ)ᵀ = X。
一种特殊的矩阵称为对称矩阵，即 Xᵀ = X。对称矩阵也是一个方阵；否则，转置的维度将不会匹配。
让我们看看实际操作。
练习 7.7 – 转置矩阵
在这个练习中，我们将使用 t() 函数转置一个矩阵，然后再转置一次，看看它是否与原始矩阵匹配：
-
创建一个元素从一到四按列填充的方阵
X：X = matrix(1:4, 2, 2) >>> X [,1] [,2] [1,] 1 3 [2,] 2 4 -
使用
t()转置矩阵：>>> t(X) [,1] [,2] [1,] 1 2 [2,] 3 4在这里，我们可以看到转置矩阵的行和列与原始矩阵相反。尽管如此，对角元素保持不变。
-
再次转置转置的矩阵并验证它是否等于原始矩阵：
>>> t(t(X)) [,1] [,2] [1,] 1 3 [2,] 2 4视觉检查显示它确实是同一个原始矩阵。然而，我们也可以使用内置的
all.equal()函数进行系统检查：>>> all.equal(X, t(t(X))) [1] TRUE结果显示这两个矩阵彼此相等，因此是同一个矩阵。
让我们回顾最后一种操作：逆一个矩阵。
求逆矩阵
逆一个标量数是直观的。对于任何数 x，其逆（或倒数）是 1/x，如果 x ≠ 0。这个条件确保 x 是可逆的。同样，并非所有矩阵都是可逆的。
形式上，我们说一个矩阵 X 是可逆的，如果它可以乘以其逆矩阵 X⁻¹，从而产生一个单位矩阵 I。换句话说，XX⁻¹ = X⁻¹X = I。
可逆矩阵也有一些有趣的性质。例如，逆一个逆矩阵会给我们原始矩阵：(X⁻¹)⁻¹ = X。另外，由于单位矩阵与自身的乘积是一个单位矩阵，因此单位矩阵的逆就是单位矩阵本身，这是一个与其他情况不同的特殊情况。
矩阵的逆可以使用 R 中的 solve() 函数获得，如果矩阵不可逆，它将给出错误。让我们通过以下练习来练习这个。
练习 7.8 – 求逆矩阵
在这个练习中，我们将使用 solve() 函数逆一个单位矩阵和一个标准方阵：
-
逆一个二维单位矩阵：
>>> solve(diag(2)) [,1] [,2] [1,] 1 0 [2,] 0 1结果显示单位矩阵的逆就是它自己。
-
从上一个例子中逆矩阵 X：
Xinv = solve(X) >>> Xinv [,1] [,2] [1,] -2 1.5 [2,] 1 -0.5我们可以通过将这个逆矩阵乘以原始矩阵并检查它是否给出单位矩阵来验证这个逆矩阵的有效性。
-
根据矩阵逆的定义，通过乘以原始矩阵来验证它：
>>> Xinv %*% X [,1] [,2] [1,] 1 0 [2,] 0 1 >>> X %*% Xinv [,1] [,2] [1,] 1 0 [2,] 0 1结果表明
Xinv确实是X的逆矩阵。
在下一节中，我们将从解线性方程组的视角来探讨矩阵-向量乘法，这在许多机器学习算法中是一个基本任务。
解线性方程组
矩阵-向量乘法操作产生了一个方程组。在典型的机器学习算法中，数据以矩阵 X 的形式出现，目标结果是向量 y。当使用的模型是一个简单的线性模型时，我们假设输入-输出关系为 Xw = y，其中 w 代表特征/系数向量。一个 n x p 的输入数据矩阵乘以一个 p x 1 的特征向量 w，将产生一个预期的 n x 1 输出向量 y。因此，线性回归的本质就是求解 w 的确切值，使得 Xw = y 中的线性方程组得到满足。
矩阵-向量乘法与线性方程组之间的等价性可能需要一段时间才能变得明显。让我们暂停一下，看看这种等价性。
线性方程组
我们已经熟悉了计算矩阵-向量乘法操作的过程。一个 2x2 的矩阵 X，当乘以一个 2x1 的向量 w 时，将得到一个 2x1 的向量 y。y 中的第一个元素，位于 (1, 1) 位置，来自 X 的第一行与 y 的第一列之间的点积（加权求和）。同样，y 中的第二个元素，位于 (2, 1) 位置，来自 X 的第二行与 y 的第一列之间的点积。y 中每个元素的位置索引决定了相应点积操作中使用的成分。
然而，这种解释是原始和低级的。一个更高级的解释在于 X 的列空间以及 X 中列向量的相关线性组合，这些列向量在相同的矩阵-向量乘法中由 w 中的权重加权。一个具体的例子将有助于说明这一点。
假设我们有一个简单的矩阵，X = [1 3 2 4]，以及一个向量，w T = [1,1]。根据通常的计算方法，输出向量 y = [1 * 1 + 3 * 1 2 * 1 + 4 * 1] = [4 6]。从列的角度看，给出了另一种计算过程：y = 1 [1 2] + 1 [3 4] = [1 2] + [3 4] = [4 6]。在这里，这两个列，[1 2] 和 [3 4]，分别由 w 中的每个元素加权。换句话说，这两个列通过线性组合产生输出列向量 y。这构成了线性方程组的基础。图 7.5 总结了本例中矩阵-向量乘法的这两种不同视角：
图 7.5 – 阐述矩阵-向量乘法的两种不同计算方式
现在，假设权重向量 w 是未知的，有两个未知元素，w1 和 w2——换句话说，wT = [w1, w2]。同时假设输入矩阵 X 和输出向量 y 是已知的。这是一个常见的情况，我们被给出输入-输出数据对，并要求估计线性模型的系数，使得 Xw = y。现在出现了线性方程组。将之前的列向量加权组合插入其中，我们得到
w1 * [1 2] + w2 * [3 4] = [w1 2 w1] + [3 w2 4 w2] = [w1 + 3 w2 2 w1 + 4 w2] = [4 6]。我们有一个方程组，如下所示：
{w1 + 3w2 = 4 2w1 + 4w2 = 6
解这个线性方程组得到 w1 = w2 = 1。矩阵-向量乘法和线性方程组的等价性在这里得到了体现。
让我们通过一个练习来欣赏代码中的等价性。
练习 7.9 – 理解线性方程组
在这个练习中，我们将首先使用常规的矩阵-向量乘法过程来计算结果，然后提供一个详细的、手动的列视图实现相同的操作。我们将看到两者给出相同的结果：
-
创建一个由四个数字 1 到 4 组成的 2x2 矩阵
X，按列填充：X = matrix(c(1:4), nrow=2, ncol=2) >>> X [,1] [,2] [1,] 1 3 [2,] 2 4 -
创建一个长度为 2 的向量
w，两个元素都填充值为 1：w = c(1,1) >>> w [1] 1 1注意，默认情况下向量是列向量。当 2x2 矩阵
X与 2x1 向量w相乘，如下面的代码所示时，向量w在乘法操作之前被表示为列向量。 -
将
X和w相乘并将结果保存在y中：y = X %*% w >>> y [,1] [1,] 4 [2,] 6现在，
y向量以显式的列向量形式显示。最后，让我们使用列视图验证计算。 -
使用列视图执行相同的矩阵-向量乘法，即列向量的加权求和：
>>> w[1] * X[,1] + w[2] * X[,2] [1] 4 6在这里，我们使用
X[,1]和X[,2]分别访问第一列和第二列。结果与之前的一致，尽管现在以行向量的格式显示，格式不同。
矩阵-向量乘法产生一个线性方程组。然而，这个线性方程组可能有解也可能没有解。即使有解，也可能不是唯一的。在下一节中，我们将检查矩阵-向量方程的潜在解。
矩阵-向量方程的解
让我们继续一个例子，其中 2x2 矩阵 X 与 2x1 权重向量 w 相乘，生成一个 2x1 输出向量 y。在这里我们将讨论三种情况：无解、有唯一解和有无穷多解。
我们将从第二种情况开始。回想一下我们之前在 X、w 和 y 中的值。在[1 3 2 4] x [1 1] = [4 6]中，我们将前两个（X 和 w）相乘以产生第三个（y），并展示了如何通过在[1 3 2 4] x [w 1 w 2] = [4 6]中解线性方程组来得到 w 的值。换句话说，解是唯一的，我们可以通过写出显式的线性方程组并求解未知数来解析地解决它。当我们有更多的未知数要解决时，这个过程会更复杂，但过程是相同的。在一个典型的线性回归设置中，当未知数的数量（权重向量 w 的长度）与方程的数量（X 的行数）相同时，并且输入矩阵 X 是一个好的矩阵（例如，如果其列是不相关的），我们可以找到线性方程组的唯一解。
现在，让我们调整输入矩阵 X，使得无法找到任何解。一种方法是将 X 的第二行向量设为零，得到 X = [1 3 0 0]。因此，我们有以下情况：
[1 3 0 0]x[w 1 w 2] = [4 6]
{w 1 + 3 w 2 = 4 0 = 6
显然，第二个方程失败了，因此这个线性方程组没有解。这被称为不一致系统，因为 0 不能等于 6。请注意，在第一个方程中，有无数对( w 1, w 2)满足该方程。这也导致了我们的第三种情况。
现在，假设我们稍微改变输出向量，使 y 的第二个元素为零，得到 y = [4 0]。因此，我们有以下情况：
[1 3 0 0]x[w 1 w 2] = [4 0]
{w 1 + 3 w 2 = 4 0 = 0
现在，这个线性方程组至少有一个解，并且是一致的，但有无穷多个解。从解这个方程组中我们无法得到太多信息；由于有无穷多个解，我们无法评估哪个是最好的报告。现实世界的优化问题通常关注于找到唯一的最优化解或最接近经验上无法达到的最优解的最佳次优解。
图 7*.6*总结了这三个情况：
图 7.6 – 解线性方程组的三个不同情况
事实上，我们也可以从几何的角度来看解线性方程组的过程。让我们深入探讨。
解线性方程组的几何解释
再次，让我们从具有两个未知数（w 1 和 w 2）和两个方程的线性方程组的唯一解的情况开始：
[1 3 2 4] x [w 1 w 2] = [4 6]
{ w 1 + 3 w 2 = 4 2 w 1 + 4 w 2 = 6
如果我们引入一个二维坐标系，并将 w₁和 w₂分别放在x轴和y轴上，我们将在坐标系上看到两条线，每条线代表系统中的一个线性方程。现在，我们可以将之前的方程组重新表达如下：
{w₂ = −w₁ * 3 + 4 * 3 w₂ = −w₁ * 2 + 3 * 2
因此，解线性方程组相当于找到这两条线的交点，因为只有在这个交点处，两个方程都得到了满足。以下代码帮助我们绘制这两条线：
plot(x=1, y=1, xlab="w1", ylab="w2", xlim=c(0, 5), ylim=c(0, 5))
abline(a=4/3, b=-1/3)
abline(a=3/2, b=-1/2)
在这里，我们使用plot()函数绘制一个二维坐标系，其中用圆表示点(1,1)。我们还使用abline()函数添加了两条线，该函数接受线的截距和斜率作为输入参数。
运行此代码生成图 7.7。我们可以看到这两条线恰好相交于点(1,1)，这并非巧合！
图 7.7 – 将线性方程组绘制为相交线
现在，我们将探讨方程组无解的情况。前面的例子显示了一个不一致的方程组，其中第二个方程根本不成立。在坐标系的几何解释中，我们可以画出一条与第一条线平行的线，例如通过将线向上移动一个单位，来产生一个无解的情况。具体来说，我们可以有以下方程组：
{w₁ + 3w₂ = 4 w₁ + 3(w₂ − 1) = 4
在这里，我们从 w₂中减去 1，将线向上移动一个单位。同样，我们可以将这些方程表达为 w₁的函数：
{w² = −w₁ * 3 + 4 * 3 w² = −w₁ * 3 + 7 * 3
让我们通过以下代码在坐标系上绘制这两条线。唯一的变化是第二条线的截距比第一条线大：
plot(x=1, y=1, xlab="w1", ylab="w2", xlim=c(0, 5), ylim=c(0, 5))
abline(a=4/3, b=-1/3)
abline(a=7/3, b=-1/3)
运行此代码生成图 7.8。由于这两条线是平行的，它们之间不会有交点，因此线性方程组将没有解：
图 7.8 – 绘制两条平行线
现在，让我们继续讨论有无穷多个解的情况。正如你可能猜到的，我们只需要将第二个方程与第一个方程相同，从而在坐标系中创建两条重叠的线。这两条重叠线上的任何点都是有效的解，并且在线上存在无限多个这样的点。
与线性方程组一起工作并不保证有解。即使未知变量的数量与系统中的行数相同，这也并不能保证我们得到一个解。在下一节中，我们将探讨在什么条件下我们保证得到一个唯一解。
获取线性方程组的唯一解
回顾我们的分析框架：Xw = y，其中我们给出了输入输出对(X, y)，并希望求解未知向量 w。如果事情很简单，并且所有这些变量都是标量，我们就会简单地两边除以 X，得到解 w = X⁻¹y，前提是 X 是可逆的。注意，我们会将两边都乘以 X⁻¹的左边来进行除法。这使我们来到了第一个条件：矩阵 X 需要有一个相应的逆矩阵；也就是说，它是可逆的。
- 图 7.9说明了常规代数中的简单标量计算与矩阵代数中的矩阵操作之间的对应关系。注意，根据矩阵逆的定义，我们有 X⁻¹X = I。我们在推导中使用了这个结果。单位矩阵是特殊的；任何矩阵乘以单位矩阵都将保持不变。因此，我们有 Iw = w。再次，我们假设矩阵 X 是可逆的：
图 7.9 – 比较解线性方程组时标量和矩阵操作
第二个条件是矩阵 X 的行列式不能为零。行列式是矩阵大小的总结性度量；我们将在下一章中更详细地讨论这个问题。第三个条件是我们之前提到的：矩阵 X 的行和列可以形成对应行和列空间的基。矩阵中不应存在相关的行或列。向量空间的基是一组向量，可以通过线性组合唯一地生成同一向量空间中的所有其他向量。
此外，我们还需要一个相对简单的条件：当 Xw = 0 时，我们有 w = 0。
在 R 中计算矩阵的逆和行列式很简单。在下面的代码中，我们使用solve()函数计算矩阵 X 的逆，使用det()函数计算其行列式：
>>> solve(X)
[,1] [,2]
[1,] -2 1.5
[2,] 1 -0.5
>>> det(X)
[1] -2
下一章将更详细地讨论矩阵行列式、范数、迹和其他特殊性质。
现在我们已经能够计算矩阵的逆，让我们学习如何使用输入矩阵 X 和输出向量 y 的逆来获得线性方程组的解。
练习 7.10 – 获得线性方程组的唯一解
在这个练习中，我们将利用矩阵逆来获得线性方程组的解：
-
根据前面的示例构建输入矩阵
X和输出向量y：X = matrix(c(1:4), nrow=2, ncol=2) >>> X [,1] [,2] [1,] 1 3 [2,] 2 4 y = c(4, 6) >>> y [1] 4 6注意，我们这里选择的数据是基于之前的运行示例，其中未知变量的解是 w = [1,1]。
-
计算线性方程组的解：
w_hat = solve(X) %*% y >>> w_hat [,1] [1,] 1 [2,] 1如我们所见，解，可以称为系数向量，与真实结果完全匹配。我们还可以将输入矩阵乘以这个解，以获得一个估计的输出，并检查它是否与给定的输出向量相同。
-
通过输入矩阵和估计的系数向量之间的矩阵-向量乘法来计算输出：
>>> X %*% w_hat [,1] [1,] 4 [2,] 6结果与给定的输出向量匹配。
使用平方输入矩阵求解线性方程组是一个简单的情况。当处理实际数据时，输入矩阵很可能是非方阵，拥有更多的行或列。当行数多于列数时，我们面临的是一个需要满足更多方程的有限未知变量，这种情况称为超定系统。另一方面，当需要满足的未知变量多于方程时，我们有一个欠定系统。我们将在下一节讨论这两种情况。
线性方程组的超定和欠定系统
我们之前讨论的例子，由两个方程和两个未知变量组成，给出了一个通过坐标系统上两条线的解。这被称为插值，其中解完美地满足两个方程，没有任何误差。在机器学习的语言中，插值意味着训练好的模型可以在训练数据集上获得 100%的分数。这并不一定是一件好事，因为模型可能会遇到过拟合的风险——也就是说，它在训练集上表现得很好，但在测试集上表现不佳。
复杂模型，如神经网络和其他非线性模型，往往具有非常低甚至为零的训练误差——这意味着它们很可能会插值训练数据。当模型变得足够复杂以至于可以插值训练数据时，就会发生零训练成本。当模型使用的系数数量，p，等于或大于训练输入中的观测数量，n，时，训练数据中的完美预测很可能会发生。在线性回归的情况下，当 p = n 时，我们解决一个线性方程组，其中自由变量的数量，即我们正在求解其最优值的数量，等于系统中的线性方程数量。
一个简单的例子是二维训练集，在坐标系统上有两个观测点，正如我们之前看到的。当只有一个参数可用时，模型在坐标系统中表现为一条水平的直线。这种单变量模型太简单，无法拟合这两个点，除非它们恰好位于拟合线上。由此产生的方程组被称为超定和欠参数化，因为我们有比系统中的未知变量更多的线性方程。
当我们有两个参数来拟合模型时，问题可以用精确解来解决：两个线性方程和两个未知数，这是一个直接的练习。我们可以通过调整其截距和斜率在任何坐标系中拟合一条直线。因此，问题变成了拟合通过两个点的直线。我们可以训练一个通过这两个点并产生零训练误差的模型，从而在这两点之间进行插值。
当存在超过两个未知变量时——比如说模型中有三个权重——问题变为欠定和过参数化。当可解时，解将是无限的。超过两个参数对应于用曲线拟合两个点。任何通过这两个点的曲线都会产生零训练误差，但曲线可以任意扭曲和复杂，这取决于使用的参数数量。图 7.10 总结了这三种情况：
图 7.10 – 模型复杂性的三种情况
图 7.10 显示了使用不同的模型复杂度拟合两个二维观测值。左图包含一个 0 次多项式模型，它只包含一个参数，因此是一条水平线。中间的图具有与观测数量相等的参数数量，使得解是唯一和精确的。右图包含具有超过两个参数的模型，由于是欠定系统，解是非唯一的。
当 p = n 时，即特征数量等于观测数量时，发生完美插值。当方程组变为欠定和过参数化时，即 p > n，完美插值继续进行，导致无限多个解，这些解对应于任意形状的模型。所有这些模型都通过观测数据点，因此给出零训练误差。
摘要
在本章中，我们介绍了线性代数的基础知识，包括处理向量和矩阵以及执行矩阵-向量乘法。我们强调了几个特殊的矩阵，例如单位矩阵，以及一些常见的操作，例如矩阵的转置和求逆。
接下来，我们使用矩阵-向量乘法在不同的设置下解决线性方程组。我们介绍了与线性方程组相对应的几何解释，以及如何使用矩阵逆和乘法运算来获得解。
最后，我们简要介绍了机器学习环境中输入矩阵的常见设置，包括欠定和过定系统。当我们深入到本书第三部分的统计建模和机器学习时，这种理解将至关重要。
在下一章中，我们将讨论矩阵代数中的一些更高级的概念以及在 R 中的实现。
第八章：R 中的中级线性代数
上一章介绍了线性代数的基础及其在 R 中的计算。本章将进一步扩展到中级线性代数，并涵盖行列式、秩和迹、特征值和特征向量以及主成分分析（PCA）等主题。除了提供对这些抽象但重要的数学概念的直观理解外，我们还将涵盖在 R 中计算这些量的实际应用。
到本章结束时，你将掌握重要的矩阵性质，如行列式和秩，并在计算这些量方面获得实践经验。
在本章中，我们将涵盖以下主题：
-
介绍矩阵行列式
-
介绍矩阵迹
-
理解矩阵范数
-
了解特征值和特征向量
-
介绍主成分分析
技术要求
要运行本章中的代码，你需要以下要求：
-
Matrix包的最新版本，写作时为 1.5.1 -
factoextra包的最新版本，写作时为 1.0.7
本章的所有代码和数据均可在github.com/PacktPublishing/The-Statistics-and-Machine-Learning-with-R-Workshop/blob/main/Chapter_8/working.R找到。
介绍矩阵行列式
矩阵的行列式是一个可以从矩阵中计算出的特殊标量值。在这里，矩阵必须是方阵，即行数和列数相等。对于一个 2x2 的方阵，行列式简单地计算为对角线元素乘积与非对角线元素乘积的差。
从数学上讲，假设我们的 2x2 矩阵是 A = [a b c d ]。因此，其行列式|A|的计算如下：
det(A) = |A| = ad − bc
请不要将这些垂直线与绝对值操作符混淆。它们代表矩阵上下文中的行列式，矩阵的行列式也可以是负数。
假设我们的 2x2 矩阵是 A = [2 6 1 8]。我们可以这样找到它的行列式：
|A| = 2 * 8 − 6 * 1 = 10
计算矩阵的行列式是容易的部分，但理解其用途同样重要。在我们介绍其性质之前，首先，我们将回顾在 R 中的计算方法，以获得对标量输出值的直观理解。
在下面的代码片段中，我们从一个向量创建矩阵 A，并使用适当的配置（两行，按行填充）。像往常一样，我们通过在控制台打印矩阵内容来验证其内容。然后我们调用det()函数来计算其行列式：
A <- matrix(c(2, 6, 1, 8), nrow=2, byrow=TRUE)
>>> A
[,1] [,2]
[1,] 2 6
[2,] 1 8
>>> det(A)
[1] 10
对于 3x3 矩阵或更高维度的矩阵，也有相应的公式来计算行列式。在这里，我们不会讨论这些情况，因为在这个阶段，理解行列式的性质更为重要。
行列式的解释
回想一下，任何矩阵都可以被视为一种变换或投影，它将输入从一种空间变换到另一种空间。对于这种变化，有两点需要注意：数量和方向。数量衡量的是矩阵原始尺寸大小的百分比变化，而方向表示变换的符号，可以是正数或负数。在这里，矩阵大小可以被认为是 2x2 矩阵的面积或 3x3 矩阵的体积。
矩阵中的列代表了一组线性变换，这些变换要么拉伸要么压缩原始输入空间，从而改变矩阵的大小。因此，行列式衡量了这组线性变换拉伸或压缩输入的程度。它给出了一个因子，表示区域面积或体积的增加或减少。此外，由于方向性很重要，变化可能会导致输入翻转，如负行列式所示。
让我们来看一个例子。假设我们有一个 2x2 的输入矩阵[1 0 0 1]，我们希望通过另一个 2x2 的矩阵[3 0 0 2]对其进行变换。直接相乘得到输出[3 0 0 2]，这可以通过执行矩阵乘法规则来验证：
[1 0 0 1][3 0 0 2] = [3 0 0 2]
输出没有变化，因为输入矩阵本质上是一个单位矩阵，我们知道任何矩阵乘以单位矩阵都不会改变。这并不奇怪。然而，当我们把[1 0 0 1]看作输入矩阵，将[3 0 0 2]作为左边的变换矩阵，将[3 0 0 2]作为右边的输出矩阵时，我们可以看到变换矩阵将输入矩阵的面积增加了六倍。
要看到这一点，想象一下输入矩阵[1 0 0 1]在二维坐标系上。输入矩阵的面积是 1 * 1 = 1，而输出矩阵的面积是 3 * 2 = 6，这恰好是变换矩阵的行列式。
这不是巧合。变换矩阵的净效应是将输入矩阵的面积放大 6 倍，保持相同的方向。当将变换矩阵改为[− 3 0 0 2]或[3 0 0 − 2]时，也不难获得相同面积的增加，但方向相反：
图 8.1 – 阐述矩阵行列式在确定输入矩阵面积变化中的作用
图 8*.1* 总结了这个重要的性质。
与矩阵秩的联系
矩阵 A 的秩是矩阵中线性无关列的最大数量。这个数字与矩阵的行列式有关。具体来说，秩是 A 中最大的平方子矩阵的行数（或列数），其行列式不为零。
让我们看看一个例子。假设 A 是一个 2x3 矩阵，[1 2 3 3 2 4]。首先，我们找到最大的平方子矩阵，即[1 2 3 2]或[2 3 2 4]。这两个矩阵的行列式都不为零。因此，A 的秩是。这意味着我们可以使用这种技术来找到矩阵的秩。图 8*.2*总结了这种方法：
图 8.2 – 使用行列式推导矩阵的秩
让我们看看如何计算矩阵的秩：
-
为了做到这一点，我们需要在 R 中加载
Matrix包并调用rankMatrix()函数。 -
如以下代码片段所示，首先，我们创建 3x2 矩阵 A 并打印出来。在设计这个矩阵时，我们只是简单地填充一个由 A 的行拼接而成的向量：
library(Matrix) A = matrix(c(1,2,3,3,2,4), nrow=2, byrow=TRUE) >>> A [,1] [,2] [,3] [1,] 1 2 3 [2,] 3 2 4 -
接下来，我们调用
rankMatrix()函数来获取它的秩：>>> rankMatrix(A) [1] 2 attr(,"method") [1] "tolNorm2" attr(,"useGrad") [1] FALSE attr(,"tol") [1] 6.661338e-16 -
返回多个属性。我们可以如下访问第一个属性：
>>> rankMatrix(A)[1] [1] 2
在下一节中，我们将探讨矩阵的另一个重要属性：迹。
矩阵迹的介绍
迹只适用于方阵，例如在机器学习中经常遇到的协方差矩阵。它用 tr(A)表示一个方阵 A，并计算为方阵中所有对角元素的和。让我们看看：
-
在下面的代码片段中，我们创建了一个 3x3 矩阵 A，并使用
diag()函数提取对角元素并将它们相加以获得矩阵的迹。请注意，我们首先创建一个包含三个列，每列有三个元素的 DataFrame，然后将它转换为矩阵格式存储在A中：A = as.matrix(data.frame("c1"=c(1,2,3),"c2"=c(2,5,2),"c3"=c(-1,8,3))) >>> A c1 c2 c3 [1,] 1 2 -1 [2,] 2 5 8 [3,] 3 2 3 >>> diag(A) [1] 1 5 3 >>> sum(diag(A)) [1] 9 -
由于没有内置函数可以一次性计算迹，我们可以构建一个自定义函数来执行此任务。如下面的代码片段所示，自定义函数
trace()本质上遍历输入方阵的所有对角元素并将它们相加作为返回值：trace <- function(A) { # get matrix dimension n = dim(A)[1] # track trace value tr = 0 # add diagonal elements to trace for(k in 1:n) { l = A[k,k] tr = tr + l } return(tr[[1]]) } -
测试这个函数给出了与之前相同的迹：
>>> trace(A) 9
关于矩阵的迹有一些有趣的属性，我们将在下一节中看到。
矩阵迹的特殊属性
为了说明这些属性，我们首先创建另一个矩阵 B：
B = as.matrix(data.frame("c1"=c(1,0,1),"c2"=c(1,1,2),"c3"=c(-1,2,0)))
>>> B
c1 c2 c3
[1,] 1 1 -1
[2,] 0 1 2
[3,] 1 2 0
我们将介绍五个在统计建模中常用到的属性。所有这些属性都将通过我们用矩阵 A 和 B 的例子进行验证：
- 属性 1：两个平方矩阵之和的迹是这两个矩阵迹的和：
tr(A + B) = tr(A) + tr(B)
>>> trace(A + B) == trace(A) + trace(B)
[1] TRUE
这里，我们使用等号来检查左边是否等于右边。返回 TRUE 表示属性已被验证。
- 属性 2：矩阵的迹等于矩阵转置的迹：
tr(A) = tr(A T)
>>> trace(A) == trace(t(A))
[1] TRUE
注意我们在使用 t() 函数来获取矩阵的转置时：
>>> t(A)
[,1] [,2] [,3]
c1 1 2 3
c2 2 5 2
c3 -1 8 3
- 属性 3：对于一个乘以标量值的矩阵，其迹与原始迹乘以相同的标量相同：
tr(αA) = αtr(A)
>>> trace(2*A) == 2*trace(A)
[1] TRUE
这里，我们选择了一个标量系数 2。您可以随意更改此值并验证属性是否仍然成立。
- 属性 4：迹是循环的：
tr(AB) = tr(BA)
>>> trace(A %*% B) == trace(B %*% A)
[1] TRUE
这个属性表明，当我们乘以 A 和 B 时，结果矩阵的迹与乘以 B 和 A 的迹相同。
注意在执行矩阵乘法时使用 %*% 符号。
- 属性 5：迹是不变的：
tr(A) = tr(BA B −1)
>>> trace(A) == trace(crossprod(crossprod(B,A),solve(B)))
[1] TRUE
这个属性表明，如果我们先乘以 B 和 A，然后乘以逆矩阵 B −1，结果矩阵的迹与矩阵 A 的迹相同。图 8*.3* 总结了这五个属性：
图 8.3 – 矩阵迹的五个属性
在下一节中，我们将介绍矩阵的另一个重要汇总度量：矩阵范数。
理解矩阵范数
矩阵的范数是一个标量值，用于衡量矩阵的大小。因此，范数是衡量向量或矩阵的大小或长度的方法。例如，深度神经网络的权重存储在矩阵中，我们通常会约束权重的范数以保持较小，以防止过拟合。这使我们能够量化大小，这在比较由多个元素组成的不同向量或矩阵时很有用。由于它从向量范数推广而来，我们将首先介绍向量范数的基本知识。
理解向量范数
假设我们有一个向量，a = [1,0, − 1]，另一个向量，b = [1,2,0]。为了评估这两个向量之间的相似性，我们可以认为它们在第一个元素上是相同的，而在剩余的两个元素上不同。为了全面比较这两个向量，我们需要一个单一的度量标准——一个可以总结整个向量的度量标准。范数是前进的一种方式。
对于长度为 n 的任意向量，有不同的范数。所有形式都来自以下 L p-范数的广义形式：
‖x‖ p = (∑ i=1 n |x i| p) 1/p
这里，双竖线表示范数。这被称为 L p-范数，因为 p 用作占位符来表示特定的范数类型。p 的常见值包括 1、2 和 ∞，尽管理论上它可以取任何正整数值。例如，为了计算向量的 L 1-范数，我们可以简单地将在公式中 p = 1。
我们将在接下来的几节中介绍这些常见的范数。
计算向量的 L 1-范数
将前一个方程中的 p 替换为 1，我们得到 L1-范数：
(||x||1 = \sum{i=1}^{n} |x_i|)
这可以被认为是向量 x 的总长度，它是向量中所有项的绝对值之和。当我们有一个二维向量，(x = [x_1, x_2])，L1-范数将是(||x||_1 = |x_1| + |x_2|)。
首先，让我们创建一个 3x1 矩阵来表示一个列向量：
a = as.matrix(c(1,2,3))
>>> a
[,1]
[1,] 1
[2,] 2
[3,] 3
我们可以在 R 中使用norm()函数来计算 L1-范数，如下所示：
>>> norm(a)
[1] 6
注意，norm()函数默认计算 L1-范数。要明确设置范数类型，我们可以传递type="1"参数，如下所示：
>>> norm(a, type="1")
[1] 6
接下来，我们将继续讨论 L2-范数。
计算向量的 L2-范数
L2-范数是我们通常使用的最常见范数。也称为欧几里得范数，L2-范数衡量两点之间的通常距离。将 p 替换为 2 到前一个公式中，我们得到以下 L2-范数的定义：
(||x||2 = \sqrt{\sum{i=1}^{n} |x_i|²} = \sqrt{\sum_{i=1}^{n} x_i²})
计算过程涉及对每个元素进行平方，将它们相加，然后取平方根。同样，当我们有一个二维向量，(x = [x_1, x_2])，L2-范数将是(||x||_2 = \sqrt{x_1² + x_2²})。
我们可以通过指定type="2"来计算向量的 L2-范数，如下所示：
>>> norm(a, type="2")
[1] 3.741657
接下来是最大范数。
计算向量的 L∞-范数
L∞-范数，或称为最大范数，找出向量中所有元素的最大绝对值。这种范数常用于最坏情况场景——例如，表示信号中注入的最大噪声。其定义如下：
(||x||{\infty} = \max{1 \leq i \leq n} |x_i|)
因此，计算过程涉及成对比较绝对值以寻找最大值。
要计算 L∞-范数，我们可以指定type="2"，如下所示：
>>> norm(a, type="I")
3
现在，我们将继续讨论矩阵范数。
理解矩阵范数
norm()函数来完成这项工作。我们将使用 X 表示一个 m×n 的一般矩阵，其中(X_{ij})表示位于第 i 行第 j 列的元素。所有矩阵范数形式都来源于以下 Lp-范数的推广形式：
(||X||p = \left(\sum{i=1}^{m} \sum_{j=1}^{n} |X_{ij}|^p\right)^{1/p})
首先，让我们创建一个 3x3 矩阵：
X = as.matrix(data.frame("c1"=c(1,2,3),"c2"=c(2,5,2),"c3"=c(-1,8,3)))
>>> X
c1 c2 c3
[1,] 1 2 -1
[2,] 2 5 8
[3,] 3 2 3
接下来，我们将查看 X 矩阵的 L1-范数。
计算矩阵的 L1-范数
矩阵的 L1-范数与其向量形式类似，但略有不同。如图所示，为了计算矩阵的 L1-范数，我们必须首先对每一列的绝对值进行求和，然后将最大的求和值作为 L1-范数：
(||X||1 = \max{1 \leq j \leq n} \sum_{i=1}^{m} |X_{ij}|)
由于求和是按列进行的，因此矩阵的 L1-范数也称为列和范数。
我们可以使用之前使用的相同命令来计算 L1-范数：
>>> norm(X, type="1")
[1] 12
通过视觉检查，我们可以看到第三列给出了绝对值求和的最大值 12。我们也可以快速检查矩阵的列绝对值求和，如下所示：
>>> colSums(abs(X))
c1 c2 c3
6 9 12
计算矩阵的 Frobenius 范数
矩阵的 L2 范数在这个阶段更为复杂，所以我们将关注一个在实践中广泛使用的类似概念：Frobenius 范数。Frobenius 范数是通过计算矩阵中所有元素的平方和然后取平方根得到的：
‖X‖F = √∑ i=1 m ∑ j=1 n |Xij|²
我们可以通过设置type="f"来计算 Frobenius 范数，如下所示：
>>> norm(X, type="f")
[1] 11
让我们通过手动平方所有项、求和并取平方根的过程来验证计算：
>>> sqrt(sum(X²))
[1] 11
现在，我们将研究矩阵的无限范数。
计算矩阵的无限范数
矩阵的无限范数与矩阵的 L1 范数类似，尽管序列的顺序不同。特别是，我们将对每一行的绝对值进行求和，然后返回最大的求和值：
‖X‖∞ = max 1≤j≤m ∑ i=1 n |Xij|
因此，无限范数也被称为norm()函数中的type="I"：
>>> norm(X, type="I")
[1] 15
再次强调，手动执行计算过程以验证结果是一个好习惯：
>>> max(rowSums(abs(X)))
[1] 15
图 8.4 总结了这些三种范数，包括向量和矩阵：
图 8.4 – 向量和矩阵的常见范数
在了解了这些基础知识之后，让我们继续下一个重要主题：特征值和特征向量。
了解特征值和特征向量
特征值，通常用标量λ表示，和特征向量，通常用 v 表示，是方阵 A 的基本属性。理解特征值和特征向量的目的需要两个核心思想。第一个思想是矩阵 A 是一个将一个输入向量映射到另一个输出向量的变换，这可能会改变方向。第二个思想是特征向量是一个特殊的向量，在经过 A 诱导的变换后不会改变方向。相反，特征向量沿着原始方向被相应的标量特征值的倍数缩放。以下方程总结了这一点：
Av = λv
这两点概括了特征分解的精髓，它用特征值和特征向量来表示原始矩阵 A，从而在许多情况下简化矩阵运算。让我们从一个简单的案例开始理解：标量-向量乘法。
理解标量-向量乘法
矩阵-向量乘法可以导致多种形式的变换，如旋转、反射、膨胀、收缩、投影以及这些操作的组合。有了特征值和特征向量，我们可以将这些操作分解成一系列更简单的操作。对于标量-向量乘法的案例，当标量值在 0 和 1 之间时，这将使向量的元素变小，从而收缩向量。
假设我们想要将一个向量 v 乘以一个标量λ，得到λv。由于 v 包含一个或多个元素，乘法本质上应用于向量的每个元素。以下代码片段显示了标量乘以向量的结果，其中每个元素都加倍：
v = c(1,2,3)
lambda = 2
>>> lambda * v
[1] 2 4 6
现在介绍一个关键技术：引入单位矩阵 I。由于向量有三个元素，我们可以将一个 3x3 的单位矩阵引入方程中。在上一章中，我们了解到乘以单位矩阵不会改变结果，因此我们可以继续这样做：
λv = λIv
在这里，我们首先使用diag()函数创建一个 3x3 的单位矩阵 I，如下所示：
I = diag(3)
>>> I
[,1] [,2] [,3]
[1,] 1 0 0
[2,] 0 1 0
[3,] 0 0 1
乘以标量λ，将λI 的对角线元素都变为 2：
>>> lambda * I
[,1] [,2] [,3]
[1,] 2 0 0
[2,] 0 2 0
[3,] 0 0 2
由于λI 是一个 3x3 矩阵，之前的标量-向量乘法现在变成了矩阵-向量乘法。这意味着我们需要切换到%*%符号来执行内积运算：
>>> (lambda * I) %*% v
[,1]
[1,] 2
[2,] 4
[3,] 6
结果与之前的标量-向量积相同，尽管现在它被表示为一个列向量而不是行向量。从矩阵-向量乘法的角度来看，矩阵通过将向量中的每个元素加倍来变换向量。
接下来，我们将正式定义特征值和特征向量的概念。
定义特征值和特征向量
通过引入单位矩阵，我们成功地将标量-向量乘法λv 转换成了矩阵-向量乘法λIv。这使理解关键方程 Av = λv 变得至关重要，其中左侧是矩阵-向量乘法，右侧是标量-向量乘法。通过写成 Av = λIv，这个方程突然在两侧具有相同类型的操作，因此更有意义。
现在，我们必须定义λ和 v，给它们合适的名字。对于一个方阵 A，我们说标量λ是 A 的特征值，与一个相关的特征向量 v ≠ 0 相关联，如果 Av = λv 成立。这个方程说明矩阵-向量乘法在 Av 中产生的向量与标量-向量乘法相同。λ和 v 一起被称为特征对。
我们可以快速验证 Av = λv 中提到的等式。假设 A = [2 3 0 1]，λ = 2，v T = [1,0]。左侧的计算给出以下结果：
Av = [2 3 0 1][1 0] = [2 0]
右侧的计算给出以下结果：
λv = 2[1 0] = [2 0]
因此，等式成立。
几何上，特征向量是在我们对矩阵 A 应用矩阵变换时保持其原始方向的向量。它在同一直线上，并且在乘法后保持不变。对于方阵，通常有一组这样的特征向量（以及相应的特征值）。
以下代码片段验证了相同的结果：
A = matrix(c(2,3,0,1), byrow=TRUE, nrow=2)
lambda = 2
v = c(1,0)
>>> A%*%v
[,1]
[1,] 2
[2,] 0
>>> lambda*v
[1] 2 0
注意，特征向量完全关注于变换不变向量的方向，而不是其大小。为了看到这一点，我们可以将特征向量加倍，并发现等式仍然成立：
Av = [2 3 0 1][2 0] = [4 0]
λv = 2[2 0] = [4 0]
我们可以通过取它们的差来验证等式：
>>> A%*%v - lambda*v
[,1]
[1,] 0
[2,] 0
图 8*.5* 总结了我们到目前为止的理解：
图 8.5 – 总结我们对特征值和特征向量的理解
接下来，我们将探讨如何计算方阵的特征值和特征向量。
计算特征值和特征向量
之前的例子假设我们能够访问特征值和特征向量。在实践中，这些需要从原始方阵中计算出来。本节重点介绍如何获得特征值和特征向量的解。
让我们从我们上次停止的地方开始。通过引入单位矩阵，我们设法获得了 Av = λIv。重新排列这些项，我们得到 Av − λIv = 0。根据定义，我们知道 v ≠ 0。根据矩阵的可逆性，如果一个矩阵的零空间（所有最终变为零的向量的集合）中存在非零向量，那么这个矩阵是不可逆的。因此，A − λI 是不可逆的。
存在一个方便的性质将矩阵的行列式与可逆性联系起来——那就是，当一个矩阵不可逆时，它的行列式必须为零。同样，如果一个矩阵是可逆的，它的行列式不能为零。因此，我们有以下结果：
det(A − λI) = 0
这给我们一个线性方程组，我们可以使用它来求解λ的值，正如前一章所示。让我们看一个具体的例子。
假设我们有一个 2x2 的方阵，A = [ 0 1 − 2 − 3]。将这个值代入前面的方程，我们得到以下结果：
det([ 0 1 − 2 − 3] − λ[1 0 0 1]) = 0
通过将标量λ乘以单位矩阵，我们得到以下结果：
det([ 0 1 − 2 − 3] − [ λ 0 0 λ ]) = 0
将两个矩阵相加得到以下结果：
det([ − λ 1 − 2 − 3 − λ]) = 0
通过应用矩阵行列式的定义，我们得到以下结果：
− λ(− 3 − λ) − 1 *(− 2) = 0
λ 2 + 3λ + 2 = 0
(λ + 1)(λ + 2) = 0
λ = − 1, − 2
因此，我们有两个解：λ = − 1 和 λ = − 2。下一步是找到这两个对应的特征向量；在接下来的阐述中，我们将关注 λ = − 1。
由于方阵 A 是 2x2 的，我们知道特征向量 v 需要是二维的。通过表示 v T = [ v 1, v 2] 并将 λ = − 1 和 A = [ 0 1 − 2 − 3] 代入 Av = λIv，我们得到以下结果：
[ 0 1 − 2 − 3][v 1 v 2] = [− v 1 − v 2]
这给我们以下方程组：
{ v 2 = − v 1 − 2 v 1 − 3 v 2 = − v 2
解这个方程组给出 v2 = − v1，这对应于无限多个解。这是有意义的，因为特征向量更关注方向而不是绝对大小。在这种情况下，方向由二维坐标系中的一条线，y = − x，表示。图 8*.6*总结了寻找方阵特征值和特征向量的过程：
图 8.6 – 推导方阵的特征值和特征向量的过程
假设我们取 v1 = 1。得到的特征向量变为 vT = [1, − 1]。当 v1 = 2 时，我们得到 vT = [2, − 2]。要使用特征值分解计算特征值和特征向量，我们可以在 R 中简单地调用eigen()函数，如下面的代码片段所示：
>>> eigen(A)
eigen() decomposition
$values
[1] 2 1
$vectors
[,1] [,2]
[1,] 1 -0.9486833
[2,] 0 0.3162278
values属性中有两个条目，vectors属性中有两个列向量，表明矩阵中总共有两个特征对。我们可以如下访问第一个特征值：
>>> eigen(A)$values[1]
[1] 2
同样，特征向量作为一组成列向量返回。因此，我们可以如下访问第一个特征向量：
>>> eigen(A)$vectors[,1]
[1] 1 0
我们可以从特征值分解中得出一些有用的性质。以一个 n × n 的方阵 A 为例。不同的特征值的数量最多是 n。
我们也可以通过使用之前推导出的条件来验证所得特征值和特征向量的正确性：det(A − λI) = 0。对于第一个特征值和特征向量，以下代码片段进行了验证：
>>> det(eigen_rst$values[1] * diag(2) - A)
[1] 0
The second eigenpair can also be verified:
>>> det(eigen_rst$values[2] * diag(2) - A)
[1] 0
我们还可以根据原始方程验证特征值分解 – 即，Av = λv：
>>> A%*%eigen_rst$vector[,1] - eigen_rst$values[1]*eigen_rst$vector[,1]
[,1]
[1,] 0
[2,] 0
>>> A%*%eigen_rst$vector[,2] - eigen_rst$values[2]*eigen_rst$vector[,2]
[,1]
[1,] -1.110223e-16
[2,] 0.000000e+00
在这里，第二个命令返回一个非常小的数字，这是由于数值近似造成的，可以认为是零。再次注意，标量-向量乘法符号*和矩阵-向量乘法符号%*%的使用。
现在我们对特征值和特征向量有了更好的理解，让我们看看一个流行的应用：主成分分析（PCA）。
介绍主成分分析
当构建机器学习模型时，用于训练模型的训练集可能包含预测变量中的冗余信息。数据集中预测变量/列的冗余来自特征的相关性，在使用某些类模型时需要加以注意。在这种情况下，主成分分析（PCA）是一种流行的技术，可以解决这类挑战，因为它减少了数据集的特征维度，从而减少了冗余。共线性问题，即模型中两个或多个预测变量线性相关，可以通过使用 PCA 进行降维来缓解。
在构建机器学习模型时，预测变量之间的共线性通常被认为是一个大问题。使用皮尔逊相关系数，它是一个介于 -1 和 1 之间的数字，其中接近 0 的系数表示两个变量线性独立，而接近 -1 或 1 的系数表示两个变量线性相关。
当两个独立变量线性相关时，例如 x₂ = 2x₁，x₂ 不提供额外信息。x₁ 和 x₂ 之间的完美相关性使得 x₂ 在解释结果变量方面变得无用。一个自然的选择是从独立变量集中删除 x₂。然而，当相关性不是完美的时候，由于删除，我们将失去一些信息。
PCA 为我们提供了一种对抗预测变量之间相关性的另一种方法。它允许我们从原始数据集中提取有意义的、不相关的信息。具体来说，它揭示了数据集背后的隐藏的低维特征。这些低维隐藏特征使得可视化解释变得方便。
为了帮助我们理解这一技术，我们将首先介绍方差-协方差矩阵的概念。
理解方差-协方差矩阵
所有机器学习模型都建立在训练数据集之上。在监督学习的背景下，数据集由输入-输出对组成。输入也称为设计矩阵，包含 n 行观测值和 p 列特征。这个 n × p 的设计矩阵 X 是我们在主成分分析（PCA）中的研究对象。
假设我们想了解每对特征之间的相关性。相关系数的范围在 -1 到 1 之间，它是基于两个变量的协方差计算的。协方差是一个标量值，用于衡量两个变量共同运动的强度。因此，设计矩阵的协方差矩阵衡量了每对独特特征之间共同运动的强度。它是一个 p × p 的方阵，cov(X)，其中位于第 i 行、第 j 列的项表示 x_i 和 x_j 特征之间的协方差值。
为了获得这个协方差矩阵，我们必须对所有特征进行均值化处理——也就是说，从给定列的每个元素中减去列均值。这消除了中心趋势，并表明了相对于均值的相对偏差量。这导致了一个均值化后的 n × p 设计矩阵 (X − μ_X)。我们必须执行以下操作：
-
将均值化后的设计矩阵转置，以获得一个 p × n 的矩阵，(X − μ_X)ᵀ。
-
将均值化后的设计矩阵的转置与原始设计矩阵相乘，两者都进行了均值化处理，以获得一个 p × p 的方阵，(X − μ_X)ᵀ(X − μ_X)。
-
将结果除以 n - 1 以对矩阵中的项进行归一化（而不是在总体协方差中使用的 n）。
完成这些操作后，将生成一个 p × p 的方差-协方差矩阵 (X − μ_X)ᵀ(X − μ_X) / (n − 1)，其中第 i 个对角元素是原始设计矩阵中第 i 列的方差。
让我们简要说明方差-协方差矩阵 (X − _ X) T(X − _ X) _ n − 1 的内容。我们将首先展示该矩阵中的项，如下所示：
(X − _ X) T(X − _ X) ____________ n − 1 = ⎡ ⎢ ⎣ ∑ (x 1 − _ x 1) 2 / (n − 1) ⋯ ∑ ( x 1 − _ x 1)( x p − _ x p) / (n − 1) ⋮ ⋱ ⋮ ∑ ( x p − _ x p)( x 1 − _ x 1) / (n − 1) ⋯ ∑ (x p − _ x p) 2 / (n − 1) ⎤ ⎥ ⎦
第一项，即 ∑ (x 1 − _ x 1) 2 / (n − 1)，是 x 1 特征的样本方差的定义。这来自于两个长度为 n 的变量的点积，后来除以 n − 1。因此，求和作用于两列向量中的所有 n 个元素。
移动到第一行的最右侧，我们有 ∑ ( x 1 − _ x 1)( x p − _ x p) / (n − 1) 作为 x 1 和 x p 变量之间的协方差。这正是我们计算这两个变量样本协方差的方式，其中求和应用于两列向量中的所有 n 个元素。图 8.7 展示了计算过程：
图 8.7 – 基于训练数据集中设计矩阵计算方差-协方差矩阵的计算过程总结
注意，虽然 n × p 设计矩阵 X 不一定是方阵，但我们设法通过将转置去均值 p × n 矩阵 (X − _ X) T 与去均值 n × p 矩阵 (X − _ X) 相乘，得到了一个 p × p 的方阵。
让我们计算给定矩阵的方差-协方差：
-
在以下代码片段中，我们首先生成一个虚拟矩阵，其中第二列是第一列的两倍：
X = matrix(c(1:5,2*(1:5)), byrow=FALSE, nrow=5) >>> X [,1] [,2] [1,] 1 2 [2,] 2 4 [3,] 3 6 [4,] 4 8 [5,] 5 10 -
接下来，我们对两列进行去均值处理：
X[,1] = X[,1] - mean(X[,1]) X[,2] = X[,2] - mean(X[,2]) >>> X [,1] [,2] [1,] -2 -4 [2,] -1 -2 [3,] 0 0 [4,] 1 2 [5,] 2 4 -
最后，我们可以这样计算方差-协方差矩阵：
>>> t(X)%*%X / (nrow(X)-1) [,1] [,2] [1,] 2.5 5 [2,] 5.0 10
结果是一个 2x2 矩阵，这与我们之前的讨论一致。
注意，我们还可以手动计算特定的方差或协方差项。例如，以下命令计算第二个变量的方差：
>>> var(X[,2])
[1] 10
我们还可以计算两个变量之间的协方差：
>>> cov(X[,1], X[,2])
[1] 5
在下一节中，我们将把方差-协方差矩阵与主成分分析（PCA）联系起来。
与 PCA 的联系
从上一节中得到的方差-协方差矩阵 (X − _ X) T(X − _ X) _ n − 1 可以用于特征值分解，这将生成一系列特征值以及相关的特征向量。注意，这些特征值将是实数标量，相关的特征向量将彼此正交，每个都指向一个不同的方向。
让我们看看与这些特征值和特征向量相关的 PCA 的一些性质。首先，数据集的总方差是这些特征值的和。因此，我们可以按降序排列这些特征值，只保留前几个，并将相关的特征向量作为下游建模的隐藏特征集。通过这样做，我们可以实现降维。
此外，这些特征向量被称为主成分，并指向特定的方向。因此，对于特定的特征向量 v i，它指向的方向可以解释数据集总方差的λ i。因此，我们可以通过累加相应的特征值，并决定保留隐藏变量的数量和解释的总方差百分比之间的权衡来解释累积方差（作为百分比）。
现在，让我们对之前的方差-协方差矩阵进行特征分解：
>>> eigen(t(X)%*%X / (nrow(X)-1))
eigen() decomposition
$values
[1] 1.250000e+01 1.110223e-16
$vectors
[,1] [,2]
[1,] 0.4472136 -0.8944272
[2,] 0.8944272 0.4472136
结果显示，第一个特征值在输出值方面显著主导第二个特征值，这表明第一个特征向量足以解释原始设计矩阵的总方差。这是有意义的，因为第二个变量只是第一个变量的两倍，因此没有提供额外的信息。尽管数据集有两个列，但它只有一个列的信息量，第二个列是完全冗余的。第二个特征向量也位于第一个特征向量所在的同一线上，因此在特征空间中没有变异性。
在下一节中，我们将介绍一个用于执行 PCA 的函数。
执行 PCA
我们可以使用prcomp()函数执行 PCA。在我们的练习中，我们将使用 Iris 数据集的前四列来执行 PCA：
-
首先，让我们加载数据集：
X = iris[,c(1:4)] >>> head(X) Sepal.Length Sepal.Width Petal.Length Petal.Width 1 5.1 3.5 1.4 0.2 2 4.9 3.0 1.4 0.2 3 4.7 3.2 1.3 0.2 4 4.6 3.1 1.5 0.2 5 5.0 3.6 1.4 0.2 6 5.4 3.9 1.7 0.4 -
prcomp()函数将自动执行所有上述涉及特征分解的步骤：X_pca = prcomp(X) >>> X_pca Standard deviations (1, .., p=4): [1] 2.0562689 0.4926162 0.2796596 0.1543862 Rotation (n x k) = (4 x 4): PC1 PC2 PC3 PC4 Sepal.Length 0.36138659 -0.65658877 0.58202985 0.3154872 Sepal.Width -0.08452251 -0.73016143 -0.59791083 -0.3197231 Petal.Length 0.85667061 0.17337266 -0.07623608 -0.4798390 Petal.Width 0.35828920 0.07548102 -0.54583143 0.7536574注意，所解释的方差（换句话说，特征值）是以标准差而不是方差的形式表示的。
现在我们已经可以访问 PCA 结果，我们可以将这些结果可视化，使其更直观。为此，我们将使用
factoextra包。请记住，首先安装并加载此包。 -
我们将要使用的第一个函数是
fviz_eig()函数，它在散点图中显示了每个主成分解释的方差百分比：>>> fviz_eig(X_pca)执行此命令生成图 8.8：
图 8.8 – 在散点图中可视化 PCA 结果
-
我们还可以使用
fviz_pca_ind()函数显示单个观测值的图表，其中具有相似轮廓的个体被分组：>>> fviz_pca_ind(X_pca, col.ind = "cos2", # Color by the quality of representation gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"), repel = TRUE # Avoid text overlapping )执行此命令生成图 8.9。请注意，质量较低（解释的方差较少）的观测值被分配到调色板的下半部分：
图 8.9 – 可视化解释总方差的个体贡献
-
最后，我们还可以可视化变量的方向。正相关变量将指向图表的同一侧，而负相关变量将指向图表的相对两侧：
>>> fviz_pca_var(X_pca, col.var = "contrib", # Color by contributions to the PC gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"), repel = TRUE # Avoid text overlapping )执行此命令生成图 8.10：
图 8.10 – 可视化变量方向
图表显示Petal.Width和Petal.Length几乎指向重叠的方向。
摘要
在本章中，我们介绍了中级线性代数及其在 R 中的实现。我们首先介绍了矩阵行列式，这是数值分析中广泛使用的属性。我们强调了矩阵行列式的直觉及其与矩阵秩的联系。
我们还介绍了额外的属性，包括矩阵迹和范数。特别是，我们介绍了三种流行的范数：L 1-范数、L 2-范数和 L ∞-范数。我们详细介绍了它们的数学结构和计算过程。
接下来，我们介绍了特征分解，它导致一个平方矩阵的特征值和特征向量集合。我们提供了核心方程的逐步推导和分析，以及计算它们的途径。
最后，我们介绍了 PCA，这是一种用于降维的流行技术。具体来说，我们强调了它在去除数据集中的共线性中的作用，并提供了计算和可视化 PCA 结果的一些方法。
在下一章中，我们将转换方向，介绍数学的另一个关键分支：微积分。
