RNA-seq 详细教程： `DESeq2` 差异表达分析（7）

学习目标

1. 了解如何设计公式
2. 了解如何使用 DESeq2 执行差异表达分析

1. DE 分析

差异表达分析工作流程的最后一步是将原始计数拟合到 NB 模型并对差异表达基因进行统计检验。在这一步中，我们本质上是想确定不同样本组的平均表达水平是否存在显著差异。

DESeq2 论文发表于 2014 年，但该软件包不断更新并通过 Bioconductor 在 R 中使用。它建立在分散估计和 DSS 和 edgeR 中的广义线性模型之上。

使用 DESeq2 进行差异表达分析涉及多个步骤，如下面流程图中蓝色部分所示。简而言之，DESeq2 将对原始计数进行建模，使用归一化因子（大小因子）来解决文库深度的差异。然后，它将估计基因方面的分散并缩小这些估计以生成更准确的分散估计来模拟计数。最后，DESeq2 将拟合负二项式模型并使用 Wald 检验或似然比检验进行假设检验。

2. 设计公式

在执行差异表达分析之前，最好通过 QC 期间的探索或先验知识了解数据中存在哪些变异来源。一旦了解了主要的变异来源，就可以在分析之前将其移除，或者通过将它们包含在设计公式中来在统计模型中对其进行控制。

设计公式告诉统计软件要控制的已知变异来源，以及差异表达测试期间要测试的感兴趣因素。例如，如果您知道 sex 是数据变异的重要来源，那么您的模型中就应该包含 sex。设计公式应该包含元数据中的所有因素，这些因素可以解释数据中主要的变化来源。公式中输入的最后一个因素应该是感兴趣的条件。

例如，假设您有以下元数据：

如果您想检查处理之间的表达差异，并且您知道变异的主要来源包括 sex 和 age，那么您的设计公式将是：

design = ~ sex + age + treatment

波浪号 (~) 应始终位于您的因子之前，并告诉 DESeq2 使用以下公式对计数进行建模。请注意，设计公式中包含的因素需要与元数据中的列名称相匹配。

2.1. 复杂情况

DESeq2 还允许分析更加复杂的情况。您可以通过在设计公式中指定来探索相互作用或“差异中的差异”。例如，如果你想探索 sex 对治疗效果的影响，你可以在设计公式中指定如下：

design = ~ sex + age + treatment + sex:treatment

由于交互项 sex:treatment 在公式的最后，因此 DESeq2 输出的结果将输出该项的结果。

3. MOV10 DE 分析

现在我们知道如何指定 DESeq2 使用的模型，可以在原始计数上运行差异表达管道。

要从我们的原始计数数据中得到我们的差异表达结果，只需要运行 2 行代码！

首先，我们像在“计数归一化”课程中所做的那样创建一个 DESeqDataSet，并指定包含我们的原始计数的 txi 对象、元数据变量，并提供我们的设计公式：

# 创建 DESeq2Dataset 对象
dds <- DESeqDataSetFromTximport(txi, colData = meta, design = ~ sampletype)

然后，为了运行差异表达分析，我们使用 DESeq() 函数。

# 运行
dds <- DESeq(dds)

通过将函数的结果重新分配回相同的变量名 (dds)，我们可以填充 DESeqDataSet 对象。

从归一化到线性建模，一切都是通过使用上面这个函数进行的！此函数将针对它执行的各个步骤打印下面的消息：

estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing

我们将在接下来的几节课中讨论每个步骤中发生的事情，但是执行这些步骤的代码在上面的两行中。