alterlab-rnaseq-quant

RNA-seq 定量 — salmon と kallisto による転写産物量の推定

バルク RNA-seq のコマンドライン定量エントリーポイントです。生の FASTQ リードと参照トランスクリプトームを入力し、salmon（選択的アラインメント）または kallisto（kb-python による疑似アラインメント）で転写産物レベルの量の推定値（カウント＋TPM）を生成します。その後、tximport／tximeta で遺伝子レベルに集約し、差次的発現解析のために alterlab-pydeseq2 へ引き渡します。このリポジトリの Python 解析スキルが、すでに実行済みであることを前提としている、生データからカウントマトリックスを作成するパイプラインです。

クイックスタート

Quantify these RNA-seq FASTQs with salmon and a decoy-aware index
Build a salmon gentrome index from this transcriptome + genome
Run kallisto / kb count on my paired-end reads
Turn my salmon quant.sf files into a gene-level count matrix for DESeq2

→ まず デコイ対応インデックスを構築し、サンプルごとに salmon quant（または kb count）を実行します。次に scripts/build_tx2gene.py と scripts/import_quant.py を実行して tximport の遺伝子マトリックスを作成し、alterlab-pydeseq2 に渡します。

---

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

軽量な定量ツールを使って、FASTQ から転写産物量または遺伝子量を取得するリクエストには、このスキルを使用します。

• 「salmon／kallisto で RNA-seq を定量してください。」
• 「デコイ対応 salmon インデックス（gentrome + decoys.txt）を構築してください。」
• 「--gcBias を使って選択的アラインメントを実行してください。」／「salmon の --sketch モードを使ってください。」
• 「quant.sf ファイルがあります。DESeq2 用の遺伝子レベルのカウントマトリックスを作成してください。」
• 「tx2gene マップを使って tximport／tximeta を設定してください。」
• 「kb-python／kb count でこれらのリードを疑似アラインメントしてください。」

該当しない場合 — 適切な関連スキルに振り分けてください

リクエストの内容…	振り分け先
カウントマトリックスを使った差次的発現解析（DESeq2 の Wald 検定、FDR、ボルケーノプロット）	alterlab-pydeseq2
シングルセル RNA-seq の定量（以前の salmon alevin）	piscem + alevin-fry — references/single_cell_alevin.md を参照。下流解析 → alterlab-scanpy／alterlab-scvi-tools
FASTQ から VCF への生殖細胞系列／体細胞変異のコール	alterlab-nf-core-sarek
16S／ITS アンプリコン（マイクロバイオーム）の FASTQ からフィーチャーテーブルへの変換	alterlab-qiime2-amplicon
空間トランスクリプトミクス（Visium／Xenium）の近傍解析	alterlab-squidpy-spatial
結果のマトリックスを AnnData オブジェクトとして読み込み／操作	alterlab-anndata
BLAST／DIAMOND による配列相同性検索	alterlab-blast
遺伝子／転写産物のID の検索と参照情報の取得（Ensembl／RefSeq）	alterlab-gget
アラインメント済みの BAM の操作、カバレッジ、アラインメントからのリードカウント	alterlab-pysam

このスキルはカウント／量のマトリックスまでを対象とします。DEG の検出、シングルセルのバーコード処理、変異コールのためのゲノムへのアラインメントは行いません。

---

重要な正確性に関する落とし穴 2 つ（定量を始める前にお読みください）

salmon が 2.x に移行した現在（ゼロから書き直された Rust 版。bioconda 2.7.0）、古い RNA-seq 手順でよく見られる間違いは次の 3 つです。いずれも上流の MIGRATION.md で確認されています（references/tool_versions.md を参照）。

1. salmon 2.0 は C++（pufferfish）インデックスを読み込めません。 書き直しに伴いインデックス形式が変更されました。古いインデックスを読み込もうとした場合や、古い salmon で 2.x のインデックスを指定した場合は検出され、明確なエラーで拒否されるため、サイレントに失敗するのではなく、はっきりと失敗します。定量に使うバージョンと同じ salmon バージョンですべてのインデックスを再構築してください。 quant.sf と bootstrap/Gibbs の出力は変更されていないため、tximport、tximeta、fishpond、swish は変更なしで 2.x の出力を引き続き利用できます。

2. --validateMappings は何もしなくなりました。 2.x では選択的アラインメントがデフォルトかつ唯一のアラインメントモードのため、このフラグは解析されて警告がログに出ます。2019 年頃のチュートリアルからコピーしても問題はありませんが、誤解を招くため削除してください。ほかにも、C++ 版のいくつかのフラグはエラーになります。--mimicBT2、--mimicStrictBT2、--minAssignedFrags、--numBiasSamples、--alternativeInitMode です。

3. salmon alevin は削除されました。 シングルセル定量は salmon に含まれなくなりました。salmon alevin … はリダイレクトを表示して終了します。代わりに piscem + alevin-fry エコシステムを使用してください。ユーザーのデータがシングルセル／ドロップレットデータの場合は、上の表に従って案内し、references/single_cell_alevin.md を参照してください。

以前の動作がどうしても必要な場合は、C++ 版の最終リリース（salmon 1.12.0）が上流の cpp ブランチにあり、salmon-cpp として別途パッケージ化されています。

---

パイプライン（デフォルトの経路である salmon）

1. デコイ対応 gentrome インデックスを構築する（参照ごとに 1 回）

デコイ対応インデックスを使うと、トランスクリプトームよりゲノムに適切にアラインするリードを salmon が識別でき、誤った割り当てを減らせます。「gentrome」は、トランスクリプト FASTA とゲノム FASTA を連結したものに相当します。さらに、ゲノム配列名をデコイとして列挙した decoys.txt を用意します。

# 1. decoys.txt = the genome's sequence (chromosome) names, one per line
grep "^>" genome.fa | sed 's/^>//; s/ .*//' > decoys.txt

# 2. gentrome = transcripts FIRST, then genome (order matters)
cat transcripts.fa genome.fa > gentrome.fa

# 3. build the index (rebuild under salmon 2.x — see trap #1)
salmon index \
  -t gentrome.fa \
  -d decoys.txt \
  -i salmon_index \
  -k 31 \
  -p 8

• -k 31 はデフォルトの k-mer です。非常に短いリードの場合に限り、値を小さくしてください。
• ヘルパー scripts/make_decoys.py が decoys.txt と gentrome.fa を作成します。また、ゲノム名がトランスクリプト FASTA に存在しない場合は処理を拒否します（よくある、気づかれにくいミスです）。references/decoy_index.md を参照してください。

2. 各サンプルを定量する

salmon quant \
  -i salmon_index \
  -l A \
  -1 sampleA_R1.fastq.gz -2 sampleA_R2.fastq.gz \
  --gcBias \
  -p 8 \
  -o quants/sampleA

• -l A — ライブラリタイプ（ストランド特異性）を自動検出します。プロトコルが文書化されている場合を除き、salmon に推定させてください。推定されたタイプは lib_format_counts.json で確認します。
• 選択的アラインメントは salmon 2.x のデフォルトです。有効にするフラグはありません（--validateMappings は受け付けられますが、無視されます）。--sketch を指定すると選択的アラインメントを無効にして、より高速なアラインメント不要の疑似アラインメントに切り替わります。定量結果を差次的発現解析に使用する場合はデフォルトを優先し、速度を最優先する大規模スクリーニングでのみ --sketch を使用してください。
• --gcBias — フラグメント単位の GC バイアスを補正します。差次的発現解析で推奨され、低コストで有効にできます。必要に応じて、5'/3' の配列特異的バイアスに対して --seqBias を追加してください。
• --ignoreTxVersion（2.x で追加）— -g/--geneMap と併用すると、tximport の同名オプションと同様に、末尾の .N を無視してトランスクリプト ID を照合します。Ensembl cDNA インデックスを Ensembl GTF と組み合わせる場合に必要です。トランスクリプトが遺伝子マップと一致しない場合、2.x は該当する名前を aux_info/genemap_unmatched_txps.json に書き込み、トランスクリプトごとに警告していた C++ 版 salmon とは異なり、一致しなかった数を添えて 1 回だけ警告します。-g を使った実行後は、このファイルがあるか確認してください。
• シングルエンドリードの場合は、-1/-2 の代わりに -r reads.fastq.gz を渡してください。

各サンプルについて、quants/<sample>/quant.sf（トランスクリプトレベルの推定値）と quants/<sample>/lib_format_counts.json（推定されたライブラリタイプ）が生成されます。フラグの詳細とサンプルごとの QC 確認については、references/salmon_quant.md を参照してください。

3. tximport で遺伝子レベルに集約する

アノテーションから転写産物→遺伝子の対応表（tx2gene）を作成し、サンプルごとの quant.sf ファイルを、pydeseq2 が読み込む遺伝子レベルの行列に要約します。

# tx2gene from a GTF/GFF3 (transcript_id -> gene_id)
uv run python scripts/build_tx2gene.py annotation.gtf --out tx2gene.tsv

# import + summarize to gene level (tximport "lengthScaledTPM" counts)
uv run python scripts/import_quant.py \
  --quants quants \
  --tx2gene tx2gene.tsv \
  --out-counts gene_counts.tsv \
  --out-tpm gene_tpm.tsv

import_quant.py は、alterlab-pydeseq2 が入力として想定する、整数に丸めた遺伝子 × サンプルのカウント行列と、 遺伝子 × サンプルの TPM 行列を生成します。これは転写産物レベルで tximport の makeCountsFromAbundance(..., "lengthScaledTPM") を実装したもので（各転写産物の TPM にサンプルごとの平均実効長を掛けてスケーリングし、各サンプル列を対応するマップ済みリードのライブラリーサイズに戻して再スケーリング）、遺伝子ごとに合計します。つまり、カウントは長さ補正済みでライブラリーサイズにスケーリングされており、NumReads をそのまま合計した値ではありません。標準的な R の方法は、Bioconductor パッケージの tximport／tximeta で countsFromAbundance = "lengthScaledTPM" を指定することです。この Python ヘルパーは同じ計算を再現するため、uv のまま作業できます（異なるのは整数への丸めと、tximeta による来歴情報がない点のみです）。 正確な仕様、tximeta の linkedTxome メタデータオプション、R の方法を選ぶべき状況については、 references/tximport_handoff.md を参照してください。

4. 差次的発現解析に引き継ぐ

gene_counts.tsv とサンプル／条件のシートを alterlab-pydeseq2 に渡します。 このスキルでは DEGs を呼び出しません。これは pydeseq2 の役割です（サイズファクター正規化、 分散、Wald 検定、BH-FDR、ボルケーノプロット／MA プロット）。

---

パイプライン（kallisto、擬似アラインメントの手法）

kallisto（スタンドアロン版 v0.52.0）と kb-python ラッパー（kb）を使うと、 より高速な擬似アラインメントの方法を利用できます。kb-python は kallisto | bustools を実行し、整った形式の出力を書き込みます。

# build a kallisto index from the transcriptome
kallisto index -i kallisto_index.idx transcripts.fa

# quantify a paired-end sample
kallisto quant -i kallisto_index.idx -o quants_kallisto/sampleA \
  sampleA_R1.fastq.gz sampleA_R2.fastq.gz

# OR the kb-python workflow (bulk)
# -f1 is the cDNA FASTA kb WRITES; trailing positionals are genome FASTA THEN GTF
kb ref -i index.idx -g t2g.txt -f1 cdna.fa genome.fa annotation.gtf
kb count -i index.idx -g t2g.txt -x bulk -o quants_kb/sampleA \
  sampleA_R1.fastq.gz sampleA_R2.fastq.gz

• kallisto は abundance.tsv／abundance.h5 を出力します。これらを salmon の quant.sf と同じ方法で tximport に渡してください（type="kallisto"）。
• ロングリード： kb-python では --long フラグ（および k>31 の k-mer）を使って lr-kallisto を利用できます。ONT／PacBio の cDNA に使用してください。詳しくは references/kallisto_kb.md を参照してください。
• kallisto は decoy／gentrome の構築を使用しません。これは salmon 固有のものです。

---

すぐに使える代替手段 — nf-core/rnaseq

トリミング → アラインメント → 定量 → QC を行う、来歴を追跡できるエンドツーエンドのパイプラインとして、 nf-core/rnaseq v3.26.0 はデフォルトで --aligner star_salmon を実行します。STAR がゲノムにマッピングしてトランスクリプトームに投影し、Salmon が定量します。手動でコマンドを実行するより、再現可能な Nextflow パイプラインをユーザーが求めている場合に利用してください。このスキルでは、salmon／kallisto を直接使う方法と tximport への引き継ぎを扱います。詳しくは references/tool_versions.md を参照してください。

---

オフロードに関する注意

インデックス作成とサンプルごとの定量は CPU／IO 負荷が高いものの、完全にオフラインで実行できます。ローカルのワークステーションでは、大きな FASTQ を API セッション経由でストリーミングするよりも、（たとえば夜間に）直接実行するのが適しています。インデックスは一度だけ作成し、ループでサンプルを定量してください。

---

報告前のセルフチェック

• salmon のインデックスを、定量に使用したものと同じ 2.x バイナリで再構築しましたか（落とし穴 #1）？C++／pufferfish インデックスはそのまま拒否されます。
• salmon のインデックスは decoy 対応（gentrome + decoys.txt）になっていますか？ゲノム名が decoys.txt に含まれていることを確認してください。
• -l A で鎖特異性を推定し、推定されたタイプを lib_format_counts.json で妥当性確認しましたか？
• データは実際にシングルセルですか？該当する場合、この方法は使用しないでください。salmon alevin は廃止されています。piscem + alevin-fry に切り替えてください（落とし穴 #2）。
• カウント行列までで処理を止め、ここで DEGs を呼び出すのではなく、差次的発現解析を alterlab-pydeseq2 に引き継ぎましたか？

---

参考資料

• references/tool_versions.md — 固定されたバージョン （salmon 2.x、kallisto v0.52.0、kb-python、nf-core/rnaseq）と、上流のリリースノートに記載された事実（2.0での書き換えとインデックスの互換性破壊、alevinの削除）。
• references/decoy_index.md — デコイ対応gentromeインデックスの構築、注意点、make_decoys.py ヘルパー。
• references/salmon_quant.md — salmon quant のフラグ一覧、ライブラリータイプの推定、サンプルごとのQC。
• references/kallisto_kb.md — kallisto / kb-pythonのワークフロー、--long（lr-kallisto）、出力の処理。
• references/tximport_handoff.md — tximport / tximetaによる集約、tx2gene、countsFromAbundance、pydeseq2への引き継ぎ。
• references/single_cell_alevin.md — salmon alevin が削除された理由と、piscem + alevin-fryによる代替手段。

AlterLab Academic Skillsスイートの一部です。