RNA-seq 2次解析 完全ガイド：発現変動遺伝子解析から機能解析まで
はじめに
本記事は、RNA-seq（トランスクリプトーム解析）における2次解析（Secondary Analysis）を体系的に解説するガイドです。1次解析（リードのQC・マッピング・定量）で得られたカウントデータを出発点として、発現変動遺伝子（DEG）の同定、機能エンリッチメント解析、可視化までの一連の流れを、初学者にもわかりやすく、中級者にも実践的な知見を提供する形でまとめました。
2次解析はRNA-seq実験の「答え」を導き出す中核的なステップです。1次解析がデータの前処理であるのに対し、2次解析では統計モデルを使って生物学的に意味のある遺伝子やパスウェイを抽出します。
本記事の対象読者
1次解析（FASTQからカウントデータの取得まで）を終えた方、またはこれから2次解析に取り組む予定の方。Rの基本操作（パッケージのインストール、data.frameの操作）ができることを前提としています。
目次
1. 2次解析の全体像と位置づけ
2. 1次解析から2次解析への橋渡し：tximport / tximeta
3. カウントデータの前処理と品質確認
4. 正規化の基礎と手法
5. 発現変動遺伝子（DEG）解析の三大ツール
6. DESeq2 による DEG 解析
7. edgeR による DEG 解析
8. limma-voom による DEG 解析
9. DEG 解析ツールの比較と選択基準
10. 多重検定補正とフィルタリング
11. DEG 結果の可視化
12. 遺伝子IDの変換：biomaRt
13. Gene Ontology（GO）エンリッチメント解析
14. KEGG パスウェイ解析
15. Gene Set Enrichment Analysis（GSEA）
16. 機能解析結果の可視化
17. 実践パイプライン：1次解析結果からの完全ワークフロー
18. バッチ効果の検出と補正
19. 複雑な実験デザインへの対応
20. よくあるトラブルと対処法
21. まとめと次のステップ
22. 参考文献・リソース
1. 2次解析の全体像と位置づけ
RNA-seq解析は大きく3段階に分けられます。1次解析（Primary Analysis）はリードのQC・トリミング・マッピング・定量を行う段階、2次解析（Secondary Analysis）はカウントデータから生物学的知見を抽出する段階、そして3次解析（Tertiary Analysis）はネットワーク解析やマルチオミクス統合など、より発展的な解析を行う段階です。
2次解析の主要ステップ
2次解析は以下のステップで構成されます。
カウントデータ（1次解析の出力）
│
▼
tximport / tximeta（定量結果の読み込み）
│
▼
品質確認・フィルタリング
│
▼
正規化
│
▼
発現変動遺伝子（DEG）解析
├── DESeq2
├── edgeR
└── limma-voom
│
▼
DEG結果の可視化
├── Volcano Plot
├── MA Plot
└── Heatmap
│
▼
機能エンリッチメント解析
├── GO解析
├── KEGGパスウェイ解析
└── GSEA
│
▼
生物学的解釈
初学者向けポイント
2次解析のゴールは「どの遺伝子が条件間で発現量に差があるか」「その遺伝子群はどのような生物学的機能に関わるか」を明らかにすることです。統計的な裏付けのある結果を得ることで、実験結果の信頼性が格段に高まります。
2. 1次解析から2次解析への橋渡し：tximport / tximeta
1次解析でSalmon、kallisto、RSEM等のツールを使った場合、トランスクリプトレベルの定量結果が得られます。2次解析では通常、遺伝子レベルのカウントデータを使用するため、トランスクリプトレベルから遺伝子レベルへの変換が必要です。ここで活躍するのがtximportとtximetaです。
tximport
バージョン：v1.38.2（Bioconductor 3.22）
論文：Soneson C, Love MI, Robinson MD (2015). F1000Research.
機能：Salmon、kallisto、RSEM、Sailfish等の出力を読み込み、遺伝子レベルに集約する
特徴：トランスクリプトレベルの推定値を遺伝子レベルに変換する際、トランスクリプト長のバイアスを適切に補正する（lengthScaledTPM等）
# tximportの基本的な使い方
library(tximport)
# Salmonの出力ファイルパスを指定
files <- file.path("salmon_output", samples$name, "quant.sf")
names(files) <- samples$name
# tx2gene: トランスクリプトIDと遺伝子IDの対応表
tx2gene <- read.csv("tx2gene.csv")
# インポート（遺伝子レベルに集約）
txi <- tximport(files, type = "salmon", tx2gene = tx2gene)
# 結果の確認
names(txi)
# [1] "abundance" "counts" "length"
# [4] "countsFromAbundance"
tximeta
バージョン：v1.28.3（Bioconductor 3.22）
論文：Love MI et al. (2020). PLOS Computational Biology.
機能：tximportの機能に加え、使用したリファレンスのメタデータ（ゲノムバージョン、アノテーション情報）を自動的に付与する
利点：再現性の向上。リファレンスのバージョンが自動記録されるため、解析の追跡が容易になる
# tximetaの使い方
library(tximeta)
# coldata: サンプル情報のデータフレーム
coldata <- data.frame(
files = file.path("salmon_output", samples$name, "quant.sf"),
names = samples$name,
condition = samples$condition
)
# Salmonの出力を自動でメタデータ付きで読み込み
se <- tximeta(coldata)
# 遺伝子レベルに集約
gse <- summarizeToGene(se)
# SummarizedExperimentオブジェクトとして利用可能
class(gse)
# [1] "RangedSummarizedExperiment"
中級者向けポイント
tximportでcountsFromAbundanceに"lengthScaledTPM"を指定すると、DESeq2やedgeRに渡す際にトランスクリプト長バイアスが補正されたカウントが得られます。Salmonの出力をそのままHTSeqのカウントのように扱うのではなく、このオフセット補正を行うことが精度向上の鍵です。
3. カウントデータの前処理と品質確認
DEG解析に入る前に、カウントデータの品質を確認し、低発現遺伝子のフィルタリングを行うことが重要です。
サンプル間の品質チェック
まず、サンプル間の全体的な関係を把握するために、主成分分析（PCA）やサンプル間相関のヒートマップを確認します。
# DESeq2でPCAプロットを作成
library(DESeq2)
# DESeqDataSetの作成
dds <- DESeqDataSetFromTximport(txi, colData = coldata, design = ~ condition)
# 分散安定化変換（VST）
vsd <- vst(dds, blind = TRUE)
# PCAプロット
plotPCA(vsd, intgroup = "condition")
PCAで確認すべきこと
（1）同じ条件のサンプルがクラスタリングしているか、（2）異なる条件のサンプルが分離しているか、（3）外れ値のサンプルがないか。もしバッチ効果が見られる場合は、後述のバッチ補正を検討します。
低発現遺伝子のフィルタリング
ほとんどのサンプルでカウントが0または極めて低い遺伝子は、統計的検定力が低く、多重検定の負担を増やすだけです。一般的に、最低でも数サンプル以上で一定のカウント数（例：10以上）がある遺伝子のみを残します。
# フィルタリングの例：少なくとも3サンプル以上でカウント10以上の遺伝子を残す
keep <- rowSums(counts(dds) >= 10) >= 3
dds <- dds[keep, ]
# フィルタリング前後の遺伝子数を確認
cat("フィルタリング前:", nrow(dds_original), "遺伝子\n")
cat("フィルタリング後:", nrow(dds), "遺伝子\n")
初学者向けポイント
フィルタリングの閾値は厳密なルールがあるわけではありません。実験のリプリケート数やシーケンス深度に応じて調整します。DESeq2では内部的にindependent filteringが行われますが、事前に明らかな低発現遺伝子を除くことで計算効率も改善します。
4. 正規化の基礎と手法
RNA-seqのカウントデータは、サンプル間のライブラリサイズ（総リード数）の違いや、遺伝子長の違いにより、そのままでは比較できません。正規化はこれらの技術的バイアスを除去し、生物学的な差異を正しく検出するために不可欠なステップです。
主要な正規化手法
DESeq2のMedian of Ratios：各サンプルのカウントを全サンプルの幾何平均に対する比率の中央値で補正する。外れ値遺伝子の影響を受けにくい堅牢な手法
edgeRのTMM（Trimmed Mean of M-values）：上位・下位の極端な発現変動を除いた後の対数比の加重平均でスケーリングファクターを算出する。Robinson & Oshlack (2010)により提案
TPM / FPKM / RPKM：遺伝子長とライブラリサイズを補正する手法。サンプル間比較にはTPMが推奨されるが、DEG解析にはカウントベースの正規化（上記2手法）を使うべき
重要な注意点
DEG解析にはTPMやFPKMではなく、生のカウントデータを使用してください。DESeq2やedgeRは内部で適切な正規化を行うため、事前にTPM変換されたデータを入力すると統計モデルが正しく機能しません。TPMは結果の報告や可視化の際に使用します。
5. 発現変動遺伝子（DEG）解析の三大ツール
RNA-seqにおけるDEG解析は、主に3つのBioconductorパッケージが広く使用されています。いずれもカウントデータの統計的分布を適切にモデリングし、条件間の発現差を検定します。
DESeq2：負の二項分布モデル。shrinkage estimatorによるfold change推定が特徴。初学者にも使いやすいインターフェース
edgeR：負の二項分布モデル。2025年にv4が発表され、大規模データセットへの対応が強化。柔軟なGLMフレームワーク
limma-voom：線形モデル＋voom重み付け。マイクロアレイ解析で培われた成熟したフレームワーク。複雑な実験デザインに強い
これら3ツールは、適切に使用すればほぼ同等の結果を返すことが多くの比較研究で示されています。それぞれの詳細を以下で解説します。
6. DESeq2 による DEG 解析
バージョン：v1.50.2（Bioconductor 3.22）
論文：Love MI, Huber W, Anders S (2014). Genome Biology, 15, 550.
ライセンス：LGPL (>= 3)
統計モデル：負の二項分布、Wald検定またはLRT（尤度比検定）
DESeq2の特徴
DESeq2は、カウントデータを負の二項分布でモデリングし、分散の推定にempirical Bayes shrinkageを適用します。これにより、リプリケート数が少ない場合でも安定した推定が可能です。
また、log2 fold change（LFC）のshrinkage推定機能があり、カウントが少ない遺伝子で生じがちな極端なfold changeを抑制できます。これはapeglm、ashr等のshrinkageメソッドを通じて利用可能です。
DESeq2の基本ワークフロー
library(DESeq2)
# 1. DESeqDataSetの作成
# tximportから
dds <- DESeqDataSetFromTximport(txi, colData = coldata, design = ~ condition)
# またはカウントマトリックスから
# dds <- DESeqDataSetFromMatrix(countData = count_matrix,
# colData = coldata,
# design = ~ condition)
# 2. リファレンスレベルの設定（対照群を指定）
dds$condition <- relevel(dds$condition, ref = "control")
# 3. 低発現遺伝子のフィルタリング
keep <- rowSums(counts(dds) >= 10) >= 3
dds <- dds[keep, ]
# 4. DESeq2の実行（正規化・分散推定・検定を一括で実行）
dds <- DESeq(dds)
# 5. 結果の取得
res <- results(dds, contrast = c("condition", "treated", "control"))
# 6. LFC shrinkage（推奨）
res_shrunk <- lfcShrink(dds, coef = "condition_treated_vs_control", type = "apeglm")
# 7. 結果のサマリー
summary(res_shrunk)
# 8. 有意な遺伝子の抽出（padj < 0.05, |LFC| > 1）
sig_genes <- subset(res_shrunk, padj < 0.05 & abs(log2FoldChange) > 1)
cat("有意な遺伝子数:", nrow(sig_genes), "\n")
中級者向けポイント：apeglm shrinkage
lfcShrinkのtypeにはnormal、apeglm、ashrの3つがあります。apeglm（Zhu et al., 2019）はベイズ推定に基づき、通常のshrinkageよりも正確で高速です。特に理由がなければapeglmの使用を推奨します。apeglmではcontrastではなくcoef引数を使う必要がある点に注意してください。
7. edgeR による DEG 解析
バージョン：v4.8.2（Bioconductor 3.22）
論文：Chen Y et al. (2025). Nucleic Acids Research, 53(2), gkaf018.（edgeR v4論文）
ライセンス：GPL (>= 2)
統計モデル：負の二項分布、exact test、準尤度（quasi-likelihood）GLM
edgeR v4の新機能
2025年に発表されたedgeR v4では、大規模データセットへの対応が大幅に強化されました。従来のexact testに加え、準尤度F検定（quasi-likelihood F-test）がデフォルトの推奨手法となり、小さいカウントや大規模データセットでの性能が向上しています。
edgeRの基本ワークフロー
library(edgeR)
# 1. DGEListの作成
y <- DGEList(counts = count_matrix, group = coldata$condition)
# tximportの結果を使う場合
# y <- DGEList(counts = txi$counts, group = coldata$condition)
# 2. 低発現遺伝子のフィルタリング
keep <- filterByExpr(y)
y <- y[keep, , keep.lib.sizes = FALSE]
# 3. TMM正規化
y <- calcNormFactors(y)
# 4. MDS（多次元尺度法）プロットで品質確認
plotMDS(y, col = as.numeric(y$samples$group))
# 5. 分散推定
design <- model.matrix(~ condition, data = coldata)
y <- estimateDisp(y, design)
# 6. 準尤度（quasi-likelihood）GLM fitting
fit <- glmQLFit(y, design)
# 7. 検定
qlf <- glmQLFTest(fit, coef = 2)
# 8. 結果の取得
topTags(qlf, n = 20)
# 全結果をデータフレームとして取得
res_edger <- topTags(qlf, n = Inf)$table
sig_edger <- subset(res_edger, FDR < 0.05 & abs(logFC) > 1)
cat("有意な遺伝子数:", nrow(sig_edger), "\n")
filterByExprについて
edgeRのfilterByExpr関数は、実験デザインを考慮した上で適切なフィルタリング閾値を自動的に決定します。最小のグループサイズに基づいてCPM閾値を算出するため、手動で閾値を設定するよりも合理的です。DESeq2で解析する場合にもこの関数は利用可能です。
8. limma-voom による DEG 解析
バージョン：v3.66.0（Bioconductor 3.22）
論文：Ritchie ME et al. (2015). Nucleic Acids Research, 43(7), e47.
ライセンス：GPL (>= 2)
統計モデル：線形モデル＋経験ベイズ（voomによる重み付け）
limma-voomの特徴
limmaは元々マイクロアレイ解析のために開発されたパッケージで、20年以上の歴史があります。voom関数により、カウントデータを対数変換し、精度重み（precision weights）を付与することで、limmaの線形モデルフレームワークをRNA-seqに適用可能にしました。
limmaの最大の強みは、複雑な実験デザイン（マルチファクター、交互作用項、ペアードデザイン、タイムコースなど）を柔軟に扱えることです。
limma-voomの基本ワークフロー
library(limma)
library(edgeR)
# 1. DGEListの作成（edgeRを使用）
y <- DGEList(counts = count_matrix, group = coldata$condition)
# 2. フィルタリングと正規化
keep <- filterByExpr(y)
y <- y[keep, , keep.lib.sizes = FALSE]
y <- calcNormFactors(y)
# 3. デザイン行列の作成
design <- model.matrix(~ condition, data = coldata)
# 4. voom変換
v <- voom(y, design, plot = TRUE)
# voomプロットで平均-分散の関係を確認
# 5. 線形モデルのフィット
fit <- lmFit(v, design)
# 6. 経験ベイズによる分散の安定化
fit <- eBayes(fit)
# 7. 結果の取得
res_limma <- topTable(fit, coef = 2, number = Inf)
sig_limma <- subset(res_limma, adj.P.Val < 0.05 & abs(logFC) > 1)
cat("有意な遺伝子数:", nrow(sig_limma), "\n")
voomプロットの読み方
voomプロットは横軸に平均対数発現量、縦軸に平方根標準偏差をプロットしたものです。適切なデータであれば、右下がりの滑らかなカーブが見られます。ばらつきが大きい場合や不規則なパターンが見られる場合は、データの品質やフィルタリングの見直しが必要です。
9. DEG 解析ツールの比較と選択基準
3つのツールはそれぞれ異なる特徴を持っています。以下の観点から適切なツールを選択します。
DESeq2が適している場合
edgeRが適している場合
limma-voomが適している場合
実践的なアドバイス
迷った場合はDESeq2から始めることをお勧めします。ドキュメントとチュートリアルが最も充実しており、Bioconductor Support Siteでの質問にも開発者のMichael Love氏が積極的に回答しています。
10. 多重検定補正とフィルタリング
RNA-seqでは数千〜数万の遺伝子を同時に検定するため、多重検定の問題が生じます。個々のp値が0.05を下回っていても、全体として偽陽性が大量に含まれてしまいます。
Benjamini-Hochberg（BH）法
最も広く使われる多重検定補正法です。False Discovery Rate（FDR）を制御し、偽陽性の割合を指定した水準（通常5%）以下に抑えます。DESeq2ではpadj、edgeRではFDR、limmaではadj.P.Valとして出力されます。
DEGの抽出基準
一般的に以下の基準が用いられます。ただし、これらは絶対的なルールではなく、実験の目的やデータの特性に応じて調整してください。
# DESeq2での有意遺伝子抽出例
# 基本的な閾値
sig <- subset(as.data.frame(res_shrunk), padj < 0.05 & abs(log2FoldChange) > 1)
# アップレギュレーション遺伝子
up_genes <- subset(sig, log2FoldChange > 0)
# ダウンレギュレーション遺伝子
down_genes <- subset(sig, log2FoldChange < 0)
cat("Up:", nrow(up_genes), " Down:", nrow(down_genes), "\n")
Independent Filtering（DESeq2）
DESeq2はindependent filteringを自動的に適用し、平均発現量の低い遺伝子をFDR補正前に除外します。これにより、検出力のない遺伝子が多重検定の負担を増やすことを防ぎ、真のDEGをより多く検出できます。
11. DEG 結果の可視化
RNA-seq 2次解析 完全ガイド：発現変動遺伝子解析から機能解析まで
はじめに
本記事は、RNA-seq（トランスクリプトーム解析）における2次解析（Secondary Analysis）を体系的に解説するガイドです。1次解析（リードのQC・マッピング・定量）で得られたカウントデータを出発点として、発現変動遺伝子（DEG）の同定、機能エンリッチメント解析、可視化までの一連の流れを、初学者にもわかりやすく、中級者にも実践的な知見を提供する形でまとめました。
2次解析はRNA-seq実験の「答え」を導き出す中核的なステップです。1次解析がデータの前処理であるのに対し、2次解析では統計モデルを使って生物学的に意味のある遺伝子やパスウェイを抽出します。
目次
1. 2次解析の全体像と位置づけ
2. 1次解析から2次解析への橋渡し：tximport / tximeta
3. カウントデータの前処理と品質確認
4. 正規化の基礎と手法
5. 発現変動遺伝子（DEG）解析の三大ツール
6. DESeq2 による DEG 解析
7. edgeR による DEG 解析
8. limma-voom による DEG 解析
9. DEG 解析ツールの比較と選択基準
10. 多重検定補正とフィルタリング
11. DEG 結果の可視化
12. 遺伝子IDの変換：biomaRt
13. Gene Ontology（GO）エンリッチメント解析
14. KEGG パスウェイ解析
15. Gene Set Enrichment Analysis（GSEA）
16. 機能解析結果の可視化
17. 実践パイプライン：1次解析結果からの完全ワークフロー
18. バッチ効果の検出と補正
19. 複雑な実験デザインへの対応
20. よくあるトラブルと対処法
21. まとめと次のステップ
22. 参考文献・リソース
1. 2次解析の全体像と位置づけ
RNA-seq解析は大きく3段階に分けられます。1次解析（Primary Analysis）はリードのQC・トリミング・マッピング・定量を行う段階、2次解析（Secondary Analysis）はカウントデータから生物学的知見を抽出する段階、そして3次解析（Tertiary Analysis）はネットワーク解析やマルチオミクス統合など、より発展的な解析を行う段階です。
2次解析の主要ステップ
2次解析は以下のステップで構成されます。
2. 1次解析から2次解析への橋渡し：tximport / tximeta
1次解析でSalmon、kallisto、RSEM等のツールを使った場合、トランスクリプトレベルの定量結果が得られます。2次解析では通常、遺伝子レベルのカウントデータを使用するため、トランスクリプトレベルから遺伝子レベルへの変換が必要です。ここで活躍するのがtximportとtximetaです。
tximport
バージョン：v1.38.2（Bioconductor 3.22）
論文：Soneson C, Love MI, Robinson MD (2015). F1000Research.
機能：Salmon、kallisto、RSEM、Sailfish等の出力を読み込み、遺伝子レベルに集約する
特徴：トランスクリプトレベルの推定値を遺伝子レベルに変換する際、トランスクリプト長のバイアスを適切に補正する（lengthScaledTPM等）
tximeta
バージョン：v1.28.3（Bioconductor 3.22）
論文：Love MI et al. (2020). PLOS Computational Biology.
機能：tximportの機能に加え、使用したリファレンスのメタデータ（ゲノムバージョン、アノテーション情報）を自動的に付与する
利点：再現性の向上。リファレンスのバージョンが自動記録されるため、解析の追跡が容易になる
3. カウントデータの前処理と品質確認
DEG解析に入る前に、カウントデータの品質を確認し、低発現遺伝子のフィルタリングを行うことが重要です。
サンプル間の品質チェック
まず、サンプル間の全体的な関係を把握するために、主成分分析（PCA）やサンプル間相関のヒートマップを確認します。
低発現遺伝子のフィルタリング
ほとんどのサンプルでカウントが0または極めて低い遺伝子は、統計的検定力が低く、多重検定の負担を増やすだけです。一般的に、最低でも数サンプル以上で一定のカウント数（例：10以上）がある遺伝子のみを残します。
4. 正規化の基礎と手法
RNA-seqのカウントデータは、サンプル間のライブラリサイズ（総リード数）の違いや、遺伝子長の違いにより、そのままでは比較できません。正規化はこれらの技術的バイアスを除去し、生物学的な差異を正しく検出するために不可欠なステップです。
主要な正規化手法
DESeq2のMedian of Ratios：各サンプルのカウントを全サンプルの幾何平均に対する比率の中央値で補正する。外れ値遺伝子の影響を受けにくい堅牢な手法
edgeRのTMM（Trimmed Mean of M-values）：上位・下位の極端な発現変動を除いた後の対数比の加重平均でスケーリングファクターを算出する。Robinson & Oshlack (2010)により提案
TPM / FPKM / RPKM：遺伝子長とライブラリサイズを補正する手法。サンプル間比較にはTPMが推奨されるが、DEG解析にはカウントベースの正規化（上記2手法）を使うべき
5. 発現変動遺伝子（DEG）解析の三大ツール
RNA-seqにおけるDEG解析は、主に3つのBioconductorパッケージが広く使用されています。いずれもカウントデータの統計的分布を適切にモデリングし、条件間の発現差を検定します。
DESeq2：負の二項分布モデル。shrinkage estimatorによるfold change推定が特徴。初学者にも使いやすいインターフェース
edgeR：負の二項分布モデル。2025年にv4が発表され、大規模データセットへの対応が強化。柔軟なGLMフレームワーク
limma-voom：線形モデル＋voom重み付け。マイクロアレイ解析で培われた成熟したフレームワーク。複雑な実験デザインに強い
これら3ツールは、適切に使用すればほぼ同等の結果を返すことが多くの比較研究で示されています。それぞれの詳細を以下で解説します。
6. DESeq2 による DEG 解析
バージョン：v1.50.2（Bioconductor 3.22）
論文：Love MI, Huber W, Anders S (2014). Genome Biology, 15, 550.
ライセンス：LGPL (>= 3)
統計モデル：負の二項分布、Wald検定またはLRT（尤度比検定）
DESeq2の特徴
DESeq2は、カウントデータを負の二項分布でモデリングし、分散の推定にempirical Bayes shrinkageを適用します。これにより、リプリケート数が少ない場合でも安定した推定が可能です。
また、log2 fold change（LFC）のshrinkage推定機能があり、カウントが少ない遺伝子で生じがちな極端なfold changeを抑制できます。これはapeglm、ashr等のshrinkageメソッドを通じて利用可能です。
DESeq2の基本ワークフロー
7. edgeR による DEG 解析
バージョン：v4.8.2（Bioconductor 3.22）
論文：Chen Y et al. (2025). Nucleic Acids Research, 53(2), gkaf018.（edgeR v4論文）
ライセンス：GPL (>= 2)
統計モデル：負の二項分布、exact test、準尤度（quasi-likelihood）GLM
edgeR v4の新機能
2025年に発表されたedgeR v4では、大規模データセットへの対応が大幅に強化されました。従来のexact testに加え、準尤度F検定（quasi-likelihood F-test）がデフォルトの推奨手法となり、小さいカウントや大規模データセットでの性能が向上しています。
edgeRの基本ワークフロー
8. limma-voom による DEG 解析
バージョン：v3.66.0（Bioconductor 3.22）
論文：Ritchie ME et al. (2015). Nucleic Acids Research, 43(7), e47.
ライセンス：GPL (>= 2)
統計モデル：線形モデル＋経験ベイズ（voomによる重み付け）
limma-voomの特徴
limmaは元々マイクロアレイ解析のために開発されたパッケージで、20年以上の歴史があります。voom関数により、カウントデータを対数変換し、精度重み（precision weights）を付与することで、limmaの線形モデルフレームワークをRNA-seqに適用可能にしました。
limmaの最大の強みは、複雑な実験デザイン（マルチファクター、交互作用項、ペアードデザイン、タイムコースなど）を柔軟に扱えることです。
limma-voomの基本ワークフロー
9. DEG 解析ツールの比較と選択基準
3つのツールはそれぞれ異なる特徴を持っています。以下の観点から適切なツールを選択します。
DESeq2が適している場合
リプリケート数が少ない（3以下）実験
初めてDEG解析を行う場合（最もドキュメントが豊富）
LFC shrinkageを利用したい場合
tximport / tximetaとの連携をスムーズに行いたい場合
edgeRが適している場合
大規模データセット（v4で性能向上）
GLMベースの柔軟な解析を行いたい場合
filterByExprによる合理的なフィルタリングを使いたい場合
exact testでシンプルな2群比較を行いたい場合
limma-voomが適している場合
複雑な実験デザイン（交互作用項、ネスト構造、タイムコース等）
サンプル数が多い場合（計算が高速）
マイクロアレイデータとの統合解析
contrast行列を柔軟に設計したい場合
10. 多重検定補正とフィルタリング
RNA-seqでは数千〜数万の遺伝子を同時に検定するため、多重検定の問題が生じます。個々のp値が0.05を下回っていても、全体として偽陽性が大量に含まれてしまいます。
Benjamini-Hochberg（BH）法
最も広く使われる多重検定補正法です。False Discovery Rate（FDR）を制御し、偽陽性の割合を指定した水準（通常5%）以下に抑えます。DESeq2ではpadj、edgeRではFDR、limmaではadj.P.Valとして出力されます。
DEGの抽出基準
一般的に以下の基準が用いられます。ただし、これらは絶対的なルールではなく、実験の目的やデータの特性に応じて調整してください。
FDR（adjusted p-value）< 0.05：統計的有意性の基準。より厳しくする場合は0.01を使用
|log2FoldChange| > 1：2倍以上の発現変動。効果量の基準。探索的な解析では0.585（1.5倍）に緩和することも
11. DEG 結果の可視化