# 第 02 章 生物信息学知识补充

> 本章为零基础读者补齐“生物信息学”这层知识：生信是什么、RNA-seq 数据长什么样、新手最容易卡住的统计学概念，以及 R 语言的第一行代码。学完本章，你应该能看懂后续章节里“limma 差异分析”“FDR 校正”“KM 生存曲线”"ROC"这些词，并亲手跑出第一个 R 脚本。

## 本章目标

1. 能用“图书馆”类比，向别人解释生物信息学是什么、能做什么；
2. 能说出 RNA-seq 的完整流程，并分清 counts、FPKM、TPM 三个量的含义与用途；
3. 能读懂表达矩阵、临床信息表等常见表格文件的格式；
4. 掌握 p 值、FDR/BH 校正、log2FC、limma、富集分析、生存分析、ROC/AUC 这些统计概念的直觉与严谨定义；
5. 能安装 R 与 RStudio，运行第一个 tidyverse 脚本（读取表达矩阵、画散点图）。

## 1. 生物信息学是什么、能做什么

把细胞想象成一座图书馆：基因组（genome）是全部藏书——人类基因组约有 30 亿个碱基对（base pair, bp），相当于一套 23 本“巨著”（23 对染色体）；基因是书里的“章节”；mRNA 是章节的“复印本”；蛋白质是读者拿走复印本后做出来的“成品”。传统生物学家像“手抄员”，一页一页读；高通量测序（high-throughput sequencing）则像一次性把整座图书馆扫描成电子版——一次实验产生 TB 级（1 TB ≈ 10¹² 字节）的数据，人类不可能逐行阅读。

生物信息学（bioinformatics）就是“海量图书的智能检索与统计”：用计算机与统计学方法，从海量生物学数据中提取规律。它能做的事包括：序列比对与基因组注释、基因表达定量、差异表达分析（找出两组样本中显著变化的基因）、功能富集分析（这些基因属于哪些通路）、生存分析（哪些基因与患者生存相关）、预后模型构建（用多个基因预测患者风险），以及公共数据库挖掘（把别人已发表的数据拿来回答自己的问题）。

本教程走的是最后一种模式——"dry lab"（纯计算实验），不需要自己做湿实验。这正是公共数据库时代给零基础新手打开的大门：你不需要一台测序仪，也能完成一个完整的课题。

## 2. 测序与表达矩阵：从 reads 到数字

RNA-seq 的标准流程：样本（肿瘤/癌旁组织）→ 提取 RNA → 建库（逆转录成 cDNA、加接头）→ 高通量测序产生读段（reads）→ 质控（去除低质量读段与接头）→ 比对（alignment）到参考基因组/转录组 → 定量（统计每个基因被多少读段覆盖）→ 得到表达矩阵。整体见图 02（示意图）：左侧是中心法则，右侧是测序把 mRNA 变成数字的过程。

![图02 中心法则与转录组测序（示意图）](../figures/fig02_central_dogma.png)

> 图 02 为教学示意图（模拟绘制），仅演示“中心法则 + 测序定量”的概念流程，不代表任何真实实验数据。

测序得到的原始计数叫 counts（读段计数）：每个基因被多少条 read 命中。counts 是整数，但它同时受两个因素干扰：基因越长越容易被命中（长度偏差）；测序越深、所有基因的计数一起变高（文库大小偏差）。因此需要归一化（normalization）才能比较：

- **FPKM**（Fragments Per Kilobase per Million，每千碱基每百万片段）：先按“每百万读段”（文库大小）归一，再按“每千碱基”（基因长度）归一。双端测序（paired-end）时统计的是片段（fragment）而非单端读段，所以叫 FPKM；
- **TPM**（Transcripts Per Million，每百万转录本）：先按基因长度归一，再按“每百万”归一。FPKM 与 TPM 在数学上可以互相转换，但 TPM 在不同样本间可比性更好，现在更推荐；
- 一句实用口诀：**counts 用于差异分析的统计模型**（limma-voom、DESeq2 等需要 counts 作为输入），**FPKM/TPM 用于展示表达水平和跨样本比较**。

表达矩阵（expression matrix）长什么样？一张表：行是基因（人类约 2 万个蛋白编码基因），列是样本（TCGA-LIHC 约 371 个肿瘤样本 + 50 个癌旁样本），单元格是表达量数字。示意如下：

| 基因 | TCGA-XX-0001（肿瘤） | TCGA-XX-0002（肿瘤） | … | TCGA-XX-00NN（癌旁） |
| --- | --- | --- | --- | --- |
| GSDMD | 8.21 | 7.94 | … | 5.12 |
| NLRP3 | 4.05 | 4.66 | … | 2.30 |

第一列通常是基因名（gene symbol，如 GSDMD）或 Ensembl ID，后续每列对应一个样本。之后所有分析（差异、富集、生存）都在这张表上展开。补充一点：TCGA 的 RNA-seq 数据提供 FPKM/FPKM-UQ 等归一化值，而 GEO 中有些队列（如 GSE14520）用的是表达芯片（microarray），单位是荧光信号强度而非读段计数——第 06 章会具体处理。

## 3. 常用数据文件格式

生信里 90% 的数据都是“表格”，最常见的两种纯文本格式：

- **CSV**（Comma-Separated Values，逗号分隔值）：列与列之间用英文逗号分隔；
- **TSV**（Tab-Separated Values，制表符分隔）：列与列之间用 Tab 分隔（基因名里偶尔会含有逗号，TSV 更保险）。

两者都可以用 Excel/WPS 打开，但文件很大时强烈建议用 R 或专业文本编辑器处理。三个新手常踩的坑：① 编码用 UTF-8，否则中文注释会乱码；② 缺失值统一记为 NA；③ 列名不要用空格和特殊符号（用下划线代替）。

教程会用到两类核心表格：

1. **表达矩阵（expression matrix）**：如上节所示，基因 × 样本；
2. **临床信息表（clinical information table）**：一行一个患者，列包含样本 ID、年龄、性别、TNM 分期、生存时间（overall survival time, OS time，单位月或天）、生存状态（OS status，0 = 存活或删失，1 = 死亡）等。生存分析全靠这张表的最后两列。

GEO 的原始数据包（SOFT/MINiML 格式）和 TCGA 的下载清单会在第 06 章详细介绍；本章先掌握“表格式文件”就够用了。

## 4. 统计分析必备概念（新手最容易卡住的地方）

这一节是本教程的“统计急救包”。每个概念都先给直觉，再给严谨定义，并说明它在示例课题里出现在哪一步。

#### 4.1 p 值

严谨定义：在原假设（null hypothesis，例如“肿瘤组与癌旁组表达无差异”）为真的前提下，观测到“当前数据或更极端数据”的概率。直觉：如果两组其实没有差异，纯靠随机抽样，得到现在这么大差异的概率有多大？这个概率就是 p 值。

新手最容易犯的两个错误：① p 值不是“两组有差异的概率”，它衡量的是“证据的强度”，而不是“效应的大小”；② p > 0.05 不代表“没有差异”，可能只是样本量太小、检验功效（power）不足。惯例阈值是 p < 0.05，意思是“如果原假设为真，出现这种结果的概率不到 5%”——这只是约定俗成的门槛，不是数学铁律。

#### 4.2 多重检验与 FDR/BH 校正

差异表达分析一次要检验约 2 万个基因。设想所有基因其实都没有差异：按 p < 0.05 的标准，纯随机也会有约 5% 的基因“碰巧”显著——也就是约 1000 个假阳性（false positive）。检验做得越多，假阳性越泛滥，这就是多重检验问题（multiple testing problem）。

解决办法是校正（adjustment）。我们控制的是 FDR（False Discovery Rate，错误发现率）：“在被判定为显著的基因里，预期有多少比例其实是假阳性”。最常用的校正是 BH 方法（Benjamini–Hochberg 法），校正后的值叫 adjusted p value（记为 adj.P 或 padj）。本教程的差异基因筛选标准统一为：`|log2FC| > 1 且 adj.P < 0.05（BH 校正）`。

#### 4.3 log2FC

FC 是倍数变化（fold change）：处理组均值 ÷ 对照组均值。FC 有一个毛病：上调 2 倍（FC = 2）和下调一半（FC = 0.5）在数值上不对称。取以 2 为底的对数后：log2(2) = 1，log2(0.5) = −1，对称且直观。所以 log2FC（log2 fold change）> 1 表示上调超过 2 倍，< −1 表示下调超过 2 倍。

#### 4.4 t 检验与 limma 的思想

t 检验（Student's t-test）比较两组均值是否显著不同，核心思想是：差异（两组均值之差）除以这个差异的“不确定性”（标准误）。差异越大、数据越稳定（标准误越小），t 值越大，p 值越小。

limma 是生信领域差异分析的“事实标准”包：它把每个基因的检验写成一个线性模型（linear model），再用经验贝叶斯（empirical Bayes）方法“借力”——用全基因组所有基因的信息来稳定单个基因的小样本方差估计，避免“样本只有 3 个、方差算不准”的尴尬。对 RNA-seq 的 counts 数据，用 limma-voom 流程（先做 voom 变换再跑 limma）。你现在不需要会推导，只要记住：**limma 的本质是“更稳的 t 检验”**。

#### 4.5 富集分析的统计学思想

拿到一张差异基因列表后，自然的问题：这些基因是不是“扎堆”出现在某些通路（如焦亡、IL-1β 信号、细胞因子）里？富集分析（enrichment analysis）回答的就是这个问题。

它的统计思想是“抽球问题”：假设全基因组有 N 个基因，其中某通路有 M 个成员；你挑了 n 个差异基因，其中属于该通路的有 k 个。问：如果差异基因是随机挑的，出现“k 个或更多通路基因”的概率是多少？这个概率由超几何分布（hypergeometric distribution）给出；等价的实现方式是 2×2 列联表上的 Fisher 精确检验（Fisher's exact test）或卡方检验（chi-square test）。概率小，说明差异基因在该通路里“富集”了。常用工具如 clusterProfiler 包（支持 GO 与 KEGG 富集），第 09 章会手把手教。

#### 4.6 生存分析基本概念

生存分析（survival analysis）研究“事件发生的时间”——在本课题里事件就是患者死亡。三个核心概念：

- **KM 生存曲线**（Kaplan–Meier curve）：把患者按某个特征（如焦亡基因高/低表达）分组，逐段估计“活到各时间点的比例”，画出阶梯状曲线；用 log-rank 检验比较两条曲线是否有显著差异（给出 p 值）。曲线越高，生存越好；
- **删失**（censoring）：患者随访结束或失访时还没观察到死亡，我们不能说“他活着”，只能说“他至少活到了这个时点”。删失数据仍然有效，KM 法就是为处理删失而设计的；
- **HR 与 Cox 回归**：风险比（hazard ratio, HR）是 Cox 比例风险回归（Cox proportional hazards regression）输出的核心量。直觉：HR = 该组单位时间死亡风险 ÷ 对照组。HR > 1 表示该因素增加死亡风险（预后差），HR < 1 表示降低风险（预后好）；通常报告 95% 置信区间（confidence interval, CI），如 HR = 1.8（95% CI 1.2–2.7）。多因素 Cox 则可以同时校正年龄、分期等混杂因素（confounder），看看焦亡基因是否“独立”于这些因素仍有预后价值。

#### 4.7 ROC 与 AUC

当我们用某个指标（如风险评分）预测“患者 3 年内会不会死亡”时，怎么评价预测准不准？ROC 曲线（receiver operating characteristic curve，受试者工作特征曲线）横轴是“1 − 特异度”（假阳性率），纵轴是“灵敏度”（真阳性率），曲线上每个点对应一个判断阈值。AUC（Area Under the Curve）是曲线下面积：AUC = 0.5 相当于瞎猜（纯随机），AUC 越接近 1 预测越准，通常 AUC > 0.7 认为有较好的区分能力。它回答的问题是：“把高风险和低风险患者分开的能力有多强？”

## 5. R 语言零基础入门

#### 5.1 为什么生信用 R，RStudio 长什么样

R 是免费开源的统计编程语言，生信圈选它的原因很实在：统计生态最全（t 检验、Cox 回归开箱即用）、有专门的生信包仓库 Bioconductor（limma、DESeq2、clusterProfiler 都在里面）、绘图能力强（ggplot2）、且完全可复现——脚本跑一遍，结果一致，这是论文写作的基本要求。

RStudio 是 R 的“驾驶舱”，界面分四块：左上脚本编辑器（source，写代码的地方，保存为 .R 文件）、左下控制台（console，代码运行与结果显示）、右上环境窗格（environment，显示已创建的变量/数据）、右下文件/绘图/包/帮助窗格（plots 里显示画出的图）。`#` 开头是注释，不会被运行；Windows 下按 Ctrl+Enter 运行光标所在行。

#### 5.2 数据结构与 tidyverse 精神

R 里最常用的两种数据结构：**向量**（vector，一维数据，用 `c()` 创建）和**数据框**（data frame，二维表格，行 = 观测、列 = 变量）。tidyverse 是一套“数据科学全家桶”（ggplot2、dplyr、readr、tidyr 等），它的精神是**整洁数据**（tidy data）：每一行是一个观测，每一列是一个变量，每个单元格是一个值。tidyverse 还有一个标志性符号 **`%>%` 管道**（读作“然后”）：把左边的结果传给右边函数的第一个参数，让代码像流水线一样从左往右读。

```r
# 向量：一维数据
x <- c(1, 2, 3, 4)

# 数据框：二维表格
df <- data.frame(
  gene  = c("GSDMD", "NLRP3"),
  logFC = c(1.2, -0.8)
)

# 三个最常用的"查看"函数
str(df)   # 查看结构
head(df)  # 查看前几行
dim(df)   # 查看维度（行数 × 列数）
```

#### 5.3 第一个脚本：读取表达矩阵、画散点图

下面这个脚本是完整的、可运行的（先给“模拟数据”版本，无需任何文件即可跑通；再给“读取真实文件”版本，第 06 章会教你如何得到那个文件）。安装包需要联网，`install.packages("tidyverse")` 执行一次即可。

```r
# ============ 第一个脚本：表达矩阵 + 散点图 ============

# 0) 安装（联网时执行一次即可，tidyverse 是真实存在的 R 包合集）
# install.packages("tidyverse")

# 1) 加载
library(tidyverse)

# 2) 方式 A：用模拟数据先跑通流程（无需任何外部文件）
set.seed(42)  # 固定随机种子，保证每次运行结果一致（可复现）
expr_demo <- tibble(
  sample = 1:50,
  GSDMD  = rnorm(50, mean = 8, sd = 1.5),   # 模拟 50 个样本的 GSDMD 表达
  NLRP3  = rnorm(50, mean = 8, sd = 1.5),   # 模拟 NLRP3 表达
  group  = rep(c("tumor", "normal"), each = 25)  # 前 25 个是肿瘤，后 25 个是癌旁
)

# 3) 画散点图：X = GSDMD，Y = NLRP3，按分组着色
#    %>%，读作"然后"：把 expr_demo 传给 ggplot 的第一个参数
expr_demo %>%
  ggplot(aes(x = GSDMD, y = NLRP3, color = group)) +
  geom_point(size = 2, alpha = 0.7) +
  labs(title = "GSDMD 与 NLRP3 表达散点图（示例数据）",
       x = "GSDMD 表达量", y = "NLRP3 表达量") +
  theme_minimal()

# 4) 方式 B：读取真实的表达矩阵（第 06 章会教你如何得到这个文件）
# expr <- read_csv("data/expr_matrix.csv")  # 文件存在时取消注释即可运行
# head(expr)
# dim(expr)
```

跑完后，右下角 plots 窗格会出现一张按“肿瘤/癌旁”着色的散点图。如果看到图，恭喜——你已经完成了生信入门的第一公里。后续章节会在这个基础上逐步加入 limma、survival、clusterProfiler 等更专业的包，但“读数据 → 处理 → 画图”这个循环，从现在起会一直陪着你。

## 本章小结

- 生物信息学 = 用计算机与统计学处理海量生物学数据；"dry lab"让新手也能用公共数据做完整课题。
- RNA-seq 流程：样本 → 建库 → 测序 → 质控 → 比对 → 定量 → 表达矩阵；counts 做统计，FPKM/TPM 做展示。
- 表达矩阵（基因 × 样本）与临床信息表（一行一患者）是本课题的两张核心表格。
- 统计急救包：p 值衡量证据强度而非效应大小；多重检验需 FDR/BH 校正；log2FC 让上下调对称；limma 是“更稳的 t 检验”；富集分析基于超几何分布/Fisher 精确检验；生存分析看 KM 曲线、HR 与 Cox 回归；ROC/AUC 评价预测区分能力。
- R + RStudio + tidyverse 是零基础最平滑的入口；`%>%` 管道让代码像流水线；第一个脚本能画出第一张图。

## 思考题

1. counts、FPKM、TPM 分别回答什么问题？为什么做差异分析通常用 counts 而不是 TPM？
2. 如果 2 万个基因其实都没有差异，用 p < 0.05 筛选，预期会有大约多少个假阳性？FDR 校正解决的是什么问题？
3. 一句报告写着"HR = 1.8（95% CI 1.2–2.7, p = 0.003）"，该怎么向非专业人士解读？
4. 一个预后模型的 AUC = 0.85，另一个 AUC = 0.55，分别意味着什么？AUC = 0.5 又代表什么？
