alterlab-pydeseq2

PyDESeq2

概要

PyDESeq2 は、バルク RNA-seq データを用いた差次的発現解析のための DESeq2 の Python 実装です。データの読み込みから結果の解釈まで、単一因子および多因子のデザイン、多重検定補正を伴う Wald 検定、任意の apeGLM による縮小、pandas および AnnData との統合を含む、完全なワークフローをサポートします。

このスキルを使用する場面

次の場合にこのスキルを使用します：

• バルク RNA-seq のカウントデータを差次的発現解析する場合
• 実験条件間（例：処理群と対照群）で遺伝子発現を比較する場合
• バッチ効果や共変量を考慮した多因子デザインを実行する場合
• R ベースの DESeq2 ワークフローを Python に移行する場合
• 差次的発現解析を Python ベースのパイプラインに統合する場合
• ユーザーが「DESeq2」、「差次的発現」、「RNA-seq 解析」、「PyDESeq2」に言及した場合

対象外

シナリオ	代わりに使用
シングルセルの差次的発現、またはクラスターのマーカー遺伝子	alterlab-scanpy（マーカー）または alterlab-scvi-tools（モデルベースの差次的発現）
FASTQ からカウント行列への変換（salmon/kallisto/STAR による定量、tximport）	alterlab-rnaseq-quant
カウントデータ以外に対する一般的な回帰／GLM／混合モデル	alterlab-statsmodels
シーケンスリードからの体細胞または生殖細胞系列のバリアントコール	alterlab-nf-core-sarek
マイクロバイオームの特徴量テーブルに対する差次的存在量解析	alterlab-qiime2-amplicon

インストールと要件

uv pip install "pydeseq2>=0.5,<0.6"

システム要件（pydeseq2 0.5.x、現行バージョンは 0.5.4）： Python ≥3.11、numpy ≥2.0、pandas ≥2.2、 scipy ≥1.12、scikit-learn ≥1.4、anndata ≥0.11、formulaic ≥1.0.2、formulaic-contrasts ≥0.2 （~ デザイン式の解析とコントラストベクトルの構築に使用）、matplotlib ≥3.9。これらは依存関係として 自動的にインストールされます。

API に関する注意（0.4 以降）： 並列処理は n_cpus= 単独のキーワード引数ではなく、inference オブジェクトで設定します：

from pydeseq2.default_inference import DefaultInference
inference = DefaultInference(n_cpus=8)
dds = DeseqDataSet(counts=counts_df, metadata=metadata, design="~condition", inference=inference)
ds = DeseqStats(dds, contrast=["condition", "treated", "control"], inference=inference)

基本的なワークフロー

1. データを準備する — カウントを サンプル × 遺伝子 の形式で読み込みます（遺伝子 × サンプルで読み込んだ場合は .T で転置します）。低カウント遺伝子（例：総リード数 < 10）をフィルタリングし、メタデータが欠損しているサンプルを除外します。
2. デザインを指定する — Wilkinson 形式の式（"~condition"、"~batch + condition"）を使用します。調整変数は対象の変数より前に指定します。
3. フィッティングする — DeseqDataSet(...).deseq2() でパイプライン全体（サイズファクター → 分散 → LFC → Cook の外れ値）を実行します。
4. 検定する — DeseqStats(dds, contrast=[var, test, ref]).summary() を実行し、results_df を確認します。
5. （任意）縮小する — ds.lfc_shrink() は可視化／順位付け専用です。p 値は縮小されません。
6. 解釈／エクスポートする — padj < 0.05 でフィルタリングし、ボルケーノプロット／MA プロットを作成して、CSV／pickle に保存します。

最小限のひな形：

from pydeseq2.dds import DeseqDataSet
from pydeseq2.ds import DeseqStats

dds = DeseqDataSet(counts=counts_df, metadata=metadata, design="~condition")
dds.deseq2()
ds = DeseqStats(dds, contrast=["condition", "treated", "control"])
ds.summary()
significant = ds.results_df[ds.results_df.padj < 0.05]

コマンドラインスクリプト

このスキルには、標準的な解析向けの完全なスタンドアロンスクリプトが含まれています：

python scripts/run_deseq2_analysis.py \
  --counts counts.csv \
  --metadata metadata.csv \
  --design "~batch + condition" \
  --contrast condition treated control \
  --output results/ \
  --min-counts 10 --alpha 0.05 --n-cpus 4 --plots

このスクリプトは、データの読み込み／検証、遺伝子とサンプルのフィルタリング、DESeq2 パイプライン全体、パラメーターをカスタマイズできる統計検定、結果のエクスポート（CSV、pickle）、および任意のボルケーノプロット／MA プロットに対応しています。複数のデータセットをバッチ処理する場合は、scripts/run_deseq2_analysis.py を参照してください。

参照先の案内

• 標準的な解析（読み込み → フィッティング → 検定 → エクスポート）または特定のデザイン（2 群間比較、複数比較、バッチ、共変量）の実行 → references/pipeline_steps.md。
• 結果の解釈、遺伝子の順位付け、ボルケーノプロット／MA プロットの作成、品質指標の確認 → references/interpretation_and_plots.md。
• エラーが発生した場合（インデックスの不一致、全カウントがゼロ、「not full rank」、有意な遺伝子がない） → references/interpretation_and_plots.md のトラブルシューティング。
• クラス／メソッドの正確なパラメーターまたはオブジェクト属性の確認 → references/api_reference.md。
• 複雑な実験デザインまたは詳細なワークフロー → references/workflow_guide.md。

重要な注意事項

1. データの向きは重要です： 通常、カウントデータは遺伝子 × サンプルの形式で読み込まれますが、サンプル × 遺伝子の形式が必要です。必要に応じて .T で転置してください。
2. サンプルのフィルタリング： 解析前に、メタデータが欠損しているサンプルを除外してください。
3. 遺伝子のフィルタリング： 検出力を高めるため、カウント数の少ない遺伝子（例：総リード数が 10 未満）を除外してください。
4. デザイン式の順序： 関心のある変数より前に調整変数を指定します（"~batch + condition"）。
5. LFC 縮小のタイミング： 検定後に、可視化や順位付けのためだけに縮小します。p 値は縮小されません。
6. 有意性： 生の p 値ではなく、padj < 0.05（Benjamini-Hochberg FDR）を使用してください。
7. コントラストの形式： [variable, test_level, reference_level]。contrast には、生物学的デザイン行列の列に対応する numpy の生コントラストベクトルも指定できます。これは、3 要素のリストでは表現できない比較（例：交互作用項、複数の水準の平均化）に使用します。
8. 中間データを保存する： 計算コストの高いフィッティングを再実行せずに済むよう、DeseqDataSet を pickle 形式で保存してください。
9. ゼロではなく、フォールドチェンジのしきい値を基準に検定する： DeseqStats(..., lfc_null=1.0, alt_hypothesis="greaterAbs") は、モデル内で「|LFC| > 1 か？」を検定します。効果量を求める統計的に適切な方法です。lfc_null=0 の結果を得た後で `abs(log2FoldChange)

   1` でフィルタリングしても、その主張に対する FDR は制御されません。

10. ゼロが多い、または疎なカウントデータ： DeseqDataSet(..., size_factors_fit_type="poscounts") は、すべてのサンプルで検出される遺伝子がない場合に失敗するデフォルトの中央値比推定器ではなく、正のカウントを使う推定器を使用します。control_genes= を指定すると、サイズファクターの推定対象をスパイクイン遺伝子またはハウスキーピング遺伝子に限定できます。

参照インデックス

• references/pipeline_steps.md — クイックスタート、コード付きのパイプライン 6 ステップ（データ準備、デザイン、フィッティング、検定、縮小、エクスポート）、および一般的な実験デザイン 4 種類。
• references/interpretation_and_plots.md — 有意な遺伝子のフィルタリング／順位付け、品質指標、ボルケーノプロットと MA プロット、トラブルシューティングガイド。
• references/api_reference.md — PyDESeq2 のクラス、メソッド、パラメーター、およびデータ構造に関する包括的なドキュメント。
• references/workflow_guide.md — 詳細な完全ワークフロー、データ読み込みパターン、多因子デザイン、ベストプラクティス。

追加リソース

• 公式ドキュメント： https://pydeseq2.readthedocs.io
• GitHub リポジトリ： https://github.com/owkin/PyDESeq2
• 論文： Muzellec et al. (2023) Bioinformatics, DOI: 10.1093/bioinformatics/btad547
• オリジナルの DESeq2（R）： Love et al. (2014) Genome Biology, DOI: 10.1186/s13059-014-0550-8