シングルセル・マルチオミクス(single-cell multiomics/scRNA-seq と scATAC-seq の統合解析)は、RNA 側と ATAC 側をそれぞれ整えてから統合します。この記事では ATAC 側の前処理を、統合(MultiVI / PeakVI)を見据えて進めます。
1. 統合を見据えた ATAC 前処理の全体像
ATAC 側の前処理は、RNA とは指標も特徴も異なりますが、ゴールは同じ ─ 最終的に統合ツール(MultiVI / PeakVI)へ生カウントで渡すところまで持っていきます。ATAC-seq は開いたクロマチンを検出する手法です(Buenrostro et al., Nat. Methods, 2013)。全体の流れを図で示します。
QC(ATAC 特有の指標)→ 特徴の整備(ピーク/タイル)→ ダブレット除去 → 生ピークカウントの退避 → 古典(TF-IDF→LSI)か深層(PeakVI / MultiVI)へ。RNA と同じく生カウントを壊さないのが要点です。実装は muon の ATAC モジュール(Bredikhin et al., 2022)を使いますが、SnapATAC2(Zhang et al., 2024)や R の Signac / ArchR も同じ役割です。ATAC 単独の詳細は scATAC-seq 解析の全体像 にもあります。
2. ATAC 特有の QC ─ TSS 濃縮・ヌクレオソーム
ATAC の QC 指標は RNA と違います。中心は、TSS 濃縮スコア(転写開始点付近にシグナルが集まるか=S/N)、ヌクレオソームシグナル(ヌクレオソーム由来フラグメントの比)、総フラグメント数、そしてピーク内フラグメント率(FRiP)です。
python
import scanpy as sc
import muon as mu
from muon import atac as ac
atac = mdata["atac"]
# ATAC 特有の QC 指標を計算
ac.tl.nucleosome_signal(atac) # ヌクレオソームシグナル
ac.tl.tss_enrichment(mdata, n_tss=2000) # TSS 濃縮 → obs["tss_score"]
sc.pp.calculate_qc_metrics(atac, inplace=True) # 総カウント・特徴数
# フィルタ:TSS 濃縮が高く・ヌクレオソーム比が低い細胞を残す
mu.pp.filter_obs(atac, "tss_score", lambda x: x > 2)
mu.pp.filter_obs(atac, "nucleosome_signal", lambda x: x < 4)
mu.pp.filter_obs(atac, "total_counts", lambda x: (x > 2000) & (x < 40000))
TSS 濃縮:TSS 付近/周辺のフラグメント比。高いほど良い(目安 2 以上)。ヌクレオソームシグナル:モノヌクレオソーム/ヌクレオソームフリーの比。高すぎ(4 超)は品質不良。総フラグメント数:深さ。低すぎ・高すぎ(ダブレット疑い)を除く。RNA の「総カウント・遺伝子数・ミト率」に対応する ATAC 版と考えると整理しやすいです。
3. 特徴の作り方 ─ ピーク vs タイル
ATAC の特徴は、ピーク(アクセス可能領域を可変長で切り出したもの。cellranger や MACS で検出;Zhang et al., Genome Biol., 2008)か、タイル(ゲノムを等間隔で区切った固定幅のビン)で表します。
どちらも行列は cells × features になります。ピークは生物学的に意味のある領域に集中でき、タイルはピークが未定義でも特徴を作れる利点があります。少数の細胞にしか出ない低頻度ピークは、統合の前に落としておくと安定します。
python
# 低頻度ピーク(ごく少数の細胞にしか出ない)を除く
sc.pp.filter_genes(atac, min_cells=int(0.01 * atac.n_obs))
4. ダブレット除去(AMULET)
ATAC のダブレットは、RNA とは別の手法で検出します。代表が AMULET で、1バーコードあたり「2コピーを超えて重なる座位」の多さからダブレットを見つけます(Thibodeau et al., Genome Biol., 2021)。AMULET は別途実行し、得られたラベルで除去します。
python
# ダブレットは ATAC 専用の手法(AMULET)で検出する
# → 得られたラベルで除去(AMULET は別途実行)
atac = atac[~atac.obs["amulet_doublet"]].copy()
RNA と同様、ATAC の QC・ダブレット除去もモダリティ単独で行います。ペアデータでは、このあと両モダリティを通過した細胞にそろえます(合同 QC)。SnapATAC2 や ArchR(Granja et al., 2021)にもダブレット検出があります。
5. 次元削減 ─ TF-IDF → LSI と、深層向けの生カウント保持
ここが RNA の HVG に対応する分岐点です。ATAC の古典的な次元削減は TF-IDF → LSI(TF-IDF 行列の SVD)で、scRNA の PCA に相当します(Cusanovich et al., Science, 2015)。ただし TF-IDF は .X を書き換えるので、深層モデル用に生ピークカウントを先に退避します。
python
from muon import atac as ac
# 生ピークカウントを退避(PeakVI / MultiVI 用)
atac.layers["counts"] = atac.X.copy()
# 古典:TF-IDF → LSI(= TF-IDF + SVD)
ac.pp.tfidf(atac, scale_factor=1e4)
ac.tl.lsi(atac)
# atac.obsm["X_lsi"] に LSI 成分が入る
LSI の第1成分はシーケンス深度と強く相関することが多く、生物学的な違いより技術的な深度を拾いがちです。慣例として第1成分を除いて2成分目以降を使います。深層(PeakVI)ではこの調整は不要で、モデルが深度を扱います。
6. PeakVI / MultiVI に渡す形に整える
深層の道では、生ピークカウントを入力に PeakVI(ATAC 単独の潜在表現。Ashuach et al., Cell Rep. Methods, 2022)を回します。ペアデータの統合では、この ATAC と RNA をまとめて MultiVI に渡します(次のブロック。Ashuach et al., Nat. Methods, 2023)。
python
import scvi
# 深層:生ピークカウントを入力に PeakVI(ATAC 単独の潜在表現)
scvi.model.PEAKVI.setup_anndata(atac, layer="counts", batch_key="batch")
model = scvi.model.PEAKVI(atac)
model.train()
atac.obsm["X_peakvi"] = model.get_latent_representation()
TF-IDF / LSI は可視化・古典解析用です。PeakVI / MultiVI には生ピークカウント(layer=’counts’)を渡します。RNA 側と同じ「可視化は変換後、モデルは生カウント」の使い分けです。
これで RNA・ATAC の両方が統合の準備を終えました。次のブロックからは、いよいよ統合(MultiVI / scGLUE)に入ります。
まとめ
- ATAC 前処理も、ゴールは「生ピークカウントで MultiVI / PeakVI に渡す」。
- QC は ATAC 特有:TSS 濃縮(2以上)・ヌクレオソームシグナル(4未満)・総フラグメント数・FRiP。
- 特徴はピーク(可変長)かタイル(固定幅)。低頻度ピークは落とす。ダブレットは AMULET。
- 古典は TF-IDF → LSI(第1成分は深度と相関 → 除く)。生カウントは先に退避。
- 深層は生カウントで PeakVI、ペア統合は MultiVI。TF-IDF 後の .X は深層に渡さない。
関連記事
- 📖 scRNA-seq 側の前処理 ─ 統合に向けて(前の記事) … RNA 側の前処理
- 📖 scATAC-seq 解析の全体像 … ATAC 単独の詳細
- 📖 scvi-tools の全体像 … PeakVI / MultiVI の土台
- 📖 AnnData のデータ構造 … データ形式の基礎
- 📖 scRNA-seq 前処理の全体像 … RNA 単独の前処理
参考文献
- Buenrostro, J. D., Giresi, P. G., Zaba, L. C., Chang, H. Y., & Greenleaf, W. J. (2013). Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nature Methods, 10(12), 1213–1218. doi:10.1038/nmeth.2688
- Bredikhin, D., Kats, I., & Stegle, O. (2022). MUON: multimodal omics analysis framework. Genome Biology, 23, 42. doi:10.1186/s13059-021-02577-8
- Cusanovich, D. A., Daza, R., Adey, A., et al. (2015). Multiplex single-cell profiling of chromatin accessibility by combinatorial cellular indexing. Science, 348(6237), 910–914. doi:10.1126/science.aab1601
- Zhang, K., Zemke, N. R., Armand, E. J., & Ren, B. (2024). A fast, scalable and versatile tool for analysis of single-cell omics data. Nature Methods, 21(2), 217–227. doi:10.1038/s41592-023-02139-9
- Stuart, T., Srivastava, A., Madad, S., Lareau, C. A., & Satija, R. (2021). Single-cell chromatin state analysis with Signac. Nature Methods, 18(11), 1333–1341. doi:10.1038/s41592-021-01282-5
- Granja, J. M., Corces, M. R., Pierce, S. E., et al. (2021). ArchR is a scalable software package for integrative single-cell chromatin accessibility analysis. Nature Genetics, 53(3), 403–411. doi:10.1038/s41588-021-00790-6
- Zhang, Y., Liu, T., Meyer, C. A., et al. (2008). Model-based Analysis of ChIP-Seq (MACS). Genome Biology, 9(9), R137. doi:10.1186/gb-2008-9-9-r137
- Thibodeau, A., Eroglu, A., McGinnis, C. S., et al. (2021). AMULET: a novel read count-based method for effective multiplet detection from single nucleus ATAC-seq data. Genome Biology, 22, 252. doi:10.1186/s13059-021-02469-x
- Ashuach, T., Reidenbach, D. A., Gayoso, A., & Yosef, N. (2022). PeakVI: A deep generative model for single-cell chromatin accessibility analysis. Cell Reports Methods, 2(3), 100182. doi:10.1016/j.crmeth.2022.100182
- Ashuach, T., Gabitto, M. I., Koodli, R. V., et al. (2023). MultiVI: deep generative model for the integration of multimodal data. Nature Methods, 20(8), 1222–1231. doi:10.1038/s41592-023-01909-9


コメント