利用转录组分析技术进行差异表达分析的具体步骤是什么？

一、实验设计与样本制备

1. 实验设计原则

• 生物学重复：每组至少设置 3 个生物学重复（如 3 个处理样本、3 个对照样本），降低个体差异导致的假阳性。
• 样本分组：明确分组标签（如 Group A/Group B），记录样本来源（如组织类型、处理时间）。
• 质量控制：提取 RNA 后检测纯度（OD260/280≈2.0，OD260/230≥2.0）和完整性（RNA 完整性指数 RIN≥7）。

2. 测序策略

• 测序深度：哺乳动物样本通常建议 10-30 million reads / 样本，低表达基因研究需更高深度。
• 测序类型：
  • 链特异性测序：保留转录方向信息，避免反义链干扰。
  • 核糖体 RNA 去除：针对 mRNA 分析时，需去除 rRNA 以富集转录本（如使用 Ribo-Zero 试剂盒）。

二、原始数据质控与预处理

1. 数据质控（QC）

• 工具：使用FastQC或MultiQC检测原始数据质量，重点关注：
  • 碱基质量分布（Phred 值≥20，即错误率≤1%）。
  • 接头污染（Adapter Content）：若超过 5% 需进行接头修剪。
  • GC 含量异常：排除基因组污染或样本降解。
• 处理：对低质量 reads 进行修剪（如Trimmomatic），去除接头和低质量碱基（通常设定滑窗质量≤20 时切除末端）。

2. 数据过滤

• 去除长度过短的 reads（如 < 20 bp）。
• 过滤核糖体 RNA（若未在建库时去除）：使用Bowtie2将 clean reads 比对到 rRNA 数据库（如 SILVA），排除比对上的 reads。

三、序列比对与转录本定量

1. 比对到参考基因组 / 转录组

• 工具选择：
  • 基因组比对（适用于已知注释物种）：
    • STAR：速度快，支持剪接位点比对，适合链特异性数据。
    • HISAT2：对可变剪接敏感，适用于复杂转录组。
  • 转录组比对（适用于无参考基因组或新转录本发现）：
    • Kallisto/Salmon：基于 k-mer 的准映射（pseudo-alignment），速度极快，直接输出表达量估计值。
• 输出结果：比对结果保存为SAM/BAM文件，需通过samtools sort排序后用于后续分析。

2. 基因表达定量

• 方法 1：基于比对的计数（适用于差异分析）
  • 工具：HTSeq-count或featureCounts，统计每个基因的 reads 计数（Count）。
  • 原理：将比对到基因外显子区域的 reads 分配至基因，支持多重比对 reads 的过滤或分配策略（如按基因长度均摊）。
• 方法 2：转录本组装与定量（适用于新转录本分析）
  • 工具：StringTie或Salmon，组装转录本并计算表达量（FPKM/TPM，需注意 TPM 在样本间可直接比较，FPKM 需标准化）。

四、差异表达分析核心流程

1. 数据标准化

• 目的：校正测序深度、基因长度和样本间差异，使计数数据可比较。
• 方法：
  • DESeq2：使用几何均值标准化（size factor normalization）。
  • edgeR：采用 Trimmed Mean of M-values（TMM）标准化。
• 输入数据：标准化后得到对数转换的表达矩阵（如 log2 (TPM+1)）。

2. 统计模型构建

• 工具选择：
  • DESeq2：适用于正态分布或负二项分布数据，支持单因素或多因素设计（如包含批次效应的模型）。
  • edgeR：基于广义线性模型（GLM），适合复杂实验设计（如时间序列、配对样本）。
  • limma-voom：适用于 RNA-Seq 或微阵列数据，对低表达基因敏感性较高。
• 模型公式示例（以 DESeq2 为例）：
  r

  ddsc <- DESeqDataSetFromMatrix(countData = counts, 
                               colData = sample_metadata, 
                               design = ~ condition) # condition为分组变量（如处理/对照）
  

   

3. 差异基因筛选

• 步骤 1：估计离散度
  • 计算基因表达的离散度（dispersion），平衡技术误差与生物学变异。
• 步骤 2：统计检验
  • 进行 Wald 检验（DESeq2）或 likelihood ratio test（LRT，适用于复杂模型），计算每个基因的差异显著性（p 值）。
• 步骤 3：多重检验校正
  • 使用 Benjamini-Hochberg 法将 p 值转换为 FDR（错误发现率），设定阈值（如 FDR < 0.05 且 |log2FC| ≥ 1）。

五、结果可视化与功能注释

1. 可视化分析

• 火山图（Volcano Plot）：
  • X 轴：log2 (FC)（差异倍数），Y 轴：-log10 (FDR)（显著性）。
  • 用颜色标记显著基因（如红色为上调，蓝色为下调），工具：ggplot2或EnhancedVolcano。
• 热图（Heatmap）：
  • 展示 Top DEG 的表达模式，聚类样本或基因，工具：pheatmap或ComplexHeatmap。
• MA 图（MA Plot）：
  • 展示基因表达量（M 值，log2 (均值)）与差异倍数（A 值，log2 (FC)）的关系，识别高表达基因的差异趋势。

2. 功能富集分析

• GO/KEGG 富集：
  • 使用clusterProfiler或GOplot对 DEG 进行功能注释，识别富集的生物学过程或通路（如 “细胞凋亡”“PI3K-AKT 信号通路”）。
• 基因集富集分析（GSEA）：
  • 适用于整体趋势一致但单个基因差异不显著的情况，工具：GSEA软件或fgsea包。

六、结果验证与生物学解读

1. 实验验证

• qRT-PCR：随机选取 5-10 个 DEG 进行定量验证，确保 RNA-Seq 结果的可靠性（相关系数 R² 应 > 0.8）。
• 蛋白质水平验证：通过 Western blot 检测关键基因的蛋白表达，排除转录后调控的影响。

2. 生物学意义解读

• 结合文献和数据库（如 GeneCards、OMIM）分析差异基因的已知功能，构建调控网络（如使用 Cytoscape 绘制共表达网络或通路图）。
• 关注枢纽基因（Hub Gene）：在共表达网络中连接度高的基因，可能是潜在的调控核心（如使用 WGCNA 识别模块枢纽基因）。

七、常用工具与代码示例（以 R 语言为例）

1. DESeq2 差异分析流程

# 安装包
if (!require("DESeq2")) install.packages("DESeq2")

# 输入数据：count矩阵（行=基因，列=样本），分组信息（如condition向量）
dds <- DESeqDataSetFromMatrix(countData = counts, 
                             colData = data.frame(condition), 
                             design = ~ condition)

# 标准化与差异分析
dds <- DESeq(dds)
res <- results(dds, contrast = c("condition", "处理组", "对照组"))
res <- res[order(res$padj), ] # 按FDR排序

# 保存结果
write.csv(as.data.frame(res), "differential_genes.csv")

2. 火山图绘制

library(ggplot2)
volcano_data <- as.data.frame(res)
volcano_data$sig <- ifelse(volcano_data$padj < 0.05 & abs(volcano_data$log2FoldChange) >= 1, "Significant", "Non-significant")

ggplot(volcano_data, aes(x = log2FoldChange, y = -log10(padj), color = sig)) +
  geom_point(alpha = 0.6, size = 2) +
  scale_color_manual(values = c("Non-significant" = "gray", "Significant" = "red")) +
  labs(x = "Log2(Fold Change)", y = "-Log10(FDR)", title = "Volcano Plot") +
  theme_bw()

关键注意事项

1. 数据分布假设：DESeq2/edgeR 基于负二项分布，适用于计数数据；若使用 TPM/FPKM 数据，需先转换为对数正态分布（如 log2 (TPM+1)）再用 limma-voom 分析。
2. 批次效应校正：若样本存在技术批次（如不同测序批次），需在模型中加入批次变量（如design = ~ batch + condition）。
3. 低表达基因过滤：分析前过滤掉在多数样本中表达量为 0 的基因（如保留至少 20% 样本中 Count≥10 的基因），减少噪声。