- 著者: Rahul Satija, Jeffrey A Farrell, David Gennert, Alexander F Schier, Aviv Regev
- Corresponding author: Alexander F Schier (Harvard University); Aviv Regev (Broad Institute of MIT and Harvard / Massachusetts Institute of Technology / Howard Hughes Medical Institute)
- 雑誌: Nature Biotechnology
- 発行年: 2015
- Epub日: 2015-04-13
- Article種別: Original Article
- PMID: 25867923
背景
組織を構成する細胞は、それぞれが固有の空間配置を持って機能している。特に発生過程や腫瘍微小環境などの複雑な生体組織において、細胞の空間的位置情報は細胞運命の決定や細胞間相互作用を理解するための極めて重要な因子である。しかし、ゲノムワイドな遺伝子発現情報を単一細胞レベルで、かつ元の空間解像度を維持したまま測定する技術は当時確立されていなかった。従来のRNA in situ hybridization (ISH) や免疫組織化学染色、蛍光レポーターを用いた手法は、一度に測定できる遺伝子数が少数の候補遺伝子に限定されるため、選択バイアス (selection bias) が生じるという課題があった。一方、近年急速に発展した単一細胞RNAシーケンシング (scRNA-seq) 技術は、ゲノムワイドな遺伝子発現プロファイルを網羅的に測定できるものの、組織を個々の細胞に解離するプロセスにおいて元の空間情報が完全に失われてしまうという致命的な欠点があった。
この課題に対し、いくつかの先行研究が異なるアプローチを試みてきた。Lovatt et al. (2014) は、生体組織内で光活性化を用いて空間的に定義された単一細胞から直接RNAを回収する TIVA (transcriptome in vivo analysis) 法を開発したが、特殊な光活性化プローブと共焦点顕微鏡システムが必要であり、スループットが極めて低いという問題があった。また、Lee et al. (2014) は、固定組織内でマルチプレックスなサブセルラーRNAシーケンシングを行う FISSEQ (fluorescent in situ sequencing) 法を報告したが、独自のローリングサークル増幅プロトコルと専用の画像解析パイプラインが必要であり、標準的なscRNA-seqプロトコルとの互換性がなかった。さらに、Durruthy-Durruthy et al. (2014) は、マウス耳胞の発生解析において主成分分析 (PCA) の第1主成分を空間的な発生軸に近似的に対応させる手法を試みたが、これは明示的な空間参照マップを持たない不完全な近似に過ぎなかった。
このように、既存の実験的手法は独自の複雑なセットアップを必要とする限定的な解決策に留まっており、標準的なscRNA-seqデータと既存のISHデータベースを統合して、細胞の元の空間配置を確率的かつ網羅的に再構成する汎用的な計算枠組みは未確立であった。特に、数十遺伝子の二値化された (binary) 粗いin situパターンという、ノイズが多くバイアスを含み得る参照データから、解離した各単一細胞をどのように高い事後確率 (posterior probability) で空間領域 (bin) に割り当てるかという統計的アプローチは解明されておらず、発生生物学におけるモルフォゲン勾配やシグナル活性化パターンをゲノムワイドに解読するための計算ツールが決定的に不足していた。したがって、標準的なscRNA-seqプロトコルと互換性があり、既存のin situデータベースを最小限の前処理で活用でき、かつマッピングの不確実性を定量化できる確率的フレームワークの開発が強く求められていた。
目的
本研究の目的は、解離した単一細胞のscRNA-seqデータと、既存の少数遺伝子におけるin situハイブリダイゼーションの空間参照マップを統合し、各細胞の元の空間位置を確率的に推定する汎用計算アルゴリズム Seurat を構築することである。さらに、発生運命決定が活発に行われるゼブラフィッシュ (Danio rerio) の50%エピボリー期 (50% epiboly stage) 胚から得た851細胞のscRNA-seqデータを用いて、Seuratの空間再構成精度を実験的・統計的に検証する。具体的には、① 既知の希少細胞集団 (脊索前板前駆細胞や内胚葉前駆細胞など) を正確な空間位置に局在化できるか、② 空間情報が未公開の遺伝子について発現パターンを予測し実験的に検証できるか、③ 未知の希少細胞集団を教師なし解析から発見しその局在を同定できるか、④ 空間情報を考慮した新規マーカー遺伝子の探索 (spatially aware marker discovery) が可能であるか、の4つの機能を実証することを目的とする。
結果
空間マッピングの精度検証: Seuratの空間マッピング精度を検証するため、まず位置が既知である n=28 cells の参照細胞 (reference cells) を用いて、予測位置と実際の回収位置との距離誤差を評価した。その結果、予測された空間重心と実測位置との距離誤差の中央値は 2 bins であり、背腹軸および動物極-植物極軸のそれぞれにおいて誤差中央値は 1 bin と極めて高い精度を示した (Fig 3c-e)。次に、ランダムに解離された n=682 cells のマッピング結果は胚全体に均一に分布したのに対し (Fig 3a)、動物極側を切除して調製したマージン濃縮細胞 n=141 cells では、動物極付近 (bin 6-8) において約 7.0-fold の顕著な細胞減少 (depletion) が観察され、赤道部 (margin) への対称的な局在濃縮が正確に再現された (Fig 3b)。さらに、47個のランドマーク遺伝子をそれぞれ1つずつ参照マップから除外してマッピングを行い、その遺伝子の空間パターンを再構成する leave-one-out 交差検証を実施した。その結果、全ランドマーク遺伝子における空間パターン予測の受信者動作特性 (ROC) 曲線下面積 (AUC) の中央値は 0.96 (median ROC = 0.96) に達し、特に 12/47 遺伝子においては ROC > 0.98 というほぼ完璧な予測精度を達成した (Fig 3f-h)。
空間アーキタイプパターンの抽出と未知遺伝子の予測:
scRNA-seqデータにおいて高い発現変動を示した 290 個の高変動遺伝子について、Seuratによって補完された空間発現プロファイルを用いて階層的クラスタリングを行った。その結果、赤道部から動物極への段階的な勾配や局在を示す 9 つの空間アーキタイプパターン (RM, VM, DEM, DRM, EM, V, DA, VA, A) が同定された (Fig 4a)。このモデルの予測能を検証するため、50%エピボリー期における空間発現パターンがこれまで報告されていなかった14個の遺伝子を選択し、実際に胚を用いた比色in situ hybridizationを実施した。その結果、tbr1b、slc25a33、pkdcca の背側赤道部 (dorsal margin) への限局、arl4ab の背側濃縮、prickle1b および dusp4 の赤道部限局、irx7、ets2、tcf3b の赤道部濃縮、そして nrarpa、id2a、insm1b、prdm12b の腹側への偏りが実験的に確認された (Fig 4b)。特に、irx7 が示す「動物極で低発現、側方赤道部で高発現」という複雑な複合パターンも正確に予測され、検証した14遺伝子のうち12遺伝子で予測と実験結果が完全一致し、残る2遺伝子 (ets2、cpn1) も予測範囲がわずかに広かったものの、発現の局在傾向は完全に一致していた (Fig 4b)。
希少細胞集団の自動同定と局在化:
scRNA-seqデータの主成分分析 (PCA) と k-means クラスタリングを組み合わせた教師なし解析により、赤道部付近に存在する希少な細胞集団を同定した。具体的には、gsc や frzb を高発現する脊索前板 (prechordal plate) 前駆細胞 n=10 cells と、sox32、cxcr4a、gata5 を高発現する内胚葉 (endodermal) 前駆細胞 n=19 cells が明確に分離された (Fig 5a-b)。Seuratは、これら n=10 cells の脊索前板前駆細胞を最も背側の赤道部 (dorsal-most margin) の単一のbinに集中してマッピングし、n=19 cells の内胚葉前駆細胞については赤道部の最下層に散在する「ソルト・アンド・ペッパー (salt-and-pepper)」状の分布パターンとして正確に再構成した (Fig 5c)。さらに、胚全体でわずか ~1/500 の頻度でしか存在しない極めて希少な細胞集団である生殖細胞 (primordial germ cells; PGC) について、ddx4/vasa、nanos3、dnd1 を極めて高く発現する単一の細胞 (n=1 cell) をデータセット内から同定した。Seuratはこの細胞を赤道部の中層 (mid-margin) にマッピングし、既知の発生解剖学的知見と一致する結果を示した。
空間情報を考慮した新規マーカー遺伝子の同定:
従来の空間情報を考慮しない単純な2群比較 (脊索前板前駆細胞 n=10 cells vs その他すべての細胞) では、gsc や nog1 などの既知の特異的マーカーに加えて、赤道部全体で広く発現する osr1 や mixl1 などの非特異的な遺伝子が偽陽性 (false positive) として多数検出される問題が生じた。これに対し、Seuratの空間マッピング情報を活用し、同じ赤道部領域に局在すると予測された細胞群のみを対照群とした比較解析 (spatially aware marker selection) を実施した。その結果、gsc、nog1、klf17、six3b などの既知の脊索前板マーカーのみが特異的に抽出され、さらに新規のマーカー候補として ripply1 および ptf1a が同定された。2色蛍光in situ hybridizationを用いた実験的検証により、ripply1 が gsc 陽性細胞の一部のサブセットにおいて共発現していることが実証され、新規マーカーとしての妥当性が確認された (Fig 5d-e)。
未知の散在性細胞集団 (アポトーシス様細胞) の発見:
全細胞を対象としたPCAにより、第4主成分において特徴的な遺伝子発現を示す n=12 cells からなる新規クラスターが同定された。これらの細胞は、アポトーシス関連遺伝子 (foxo3b、tp53inp1、casp8、ctsh)、細胞ストレスマーカー (isg15、sesn3、mat2al、gadd45aa)、およびシグナル伝達因子 (igf2a、aplnrb) を特異的に強発現していた (Fig 5f-h)。Subramanian et al. ProcNatlAcadSciUSA 2005 の手法に基づく遺伝子セット富裕化解析 (GSEA) を実施した結果、p53シグナル伝達経路の標的遺伝子群に対して極めて有意な濃縮が確認された (FDR < 10⁻⁶)。Seuratは、これらのアポトーシス様細胞が胚全体に確率的に散在し、特に動物極および腹側極にやや偏って分布していると予測した (Fig 5g)。この予測を検証するため、casp8、gadd45aa、igf2a、tp53inp1 のin situ hybridizationを実施した結果、実際に個々の細胞が胚全体にランダムに散在するパターンが確認された (Fig 5i)。さらに、2色蛍光in situ hybridizationにより、同一細胞内において aplnrb と isg15 が共発現していることが実証された (Fig 5j)。この細胞集団は、解析した10個の独立した胚すべてにおいて一貫して検出され、細胞解離時のアーティファクトではなく、初期発生において stochastically に発生する生理的なストレス細胞集団であることが示された。
ランドマーク遺伝子数の最適化: 空間参照マップに必要なランドマーク遺伝子の数を評価するため、入力する遺伝子数を段階的に削減するダウンサンプリング解析を実施した。その結果、Seuratによる空間マッピングの精度は ~30個のランドマーク遺伝子 を境界として安定し、それ以上の遺伝子を追加しても精度の向上は緩やかになることが示された。また、特定の空間パターンを示す遺伝子のみを重複して追加するよりも、同定された9つの空間アーキタイプパターン全体から均等に遺伝子を選択して参照マップを構成する方が、マッピングの安定性と精度において極めて効果的であることが明らかになった。
考察/結論
先行研究との違い: 本研究で開発された Seurat は、解離した単一細胞のゲノムワイドな発現プロファイルと、既存の少数遺伝子のin situ空間情報を確率的に統合する新しい計算フレームワークであり、従来の空間トランスクリプトミクス技術とは本質的に異なる位置付けにある。LovattらのTIVA (transcriptome in vivo analysis) 法やLeeらのFISSEQ (fluorescent in situ sequencing) 法は、組織内での直接的なシーケンシングを試みているものの、特殊な実験装置やプローブを必要とし、検出感度やスループットに大きな制限があった。これらと対照的に、Seuratは標準的なscRNA-seqプロトコルから得られる高感度なデータと、既存の公的データベースに蓄積された簡便な比色in situパターンのみをインプットとして動作するため、実験的な導入障壁が極めて低い。また、Durruthy-Durruthyらが用いたPCAの主成分による空間軸の近似とは異なり、Seuratは明示的な二値化空間参照マップをガウス混合モデルを介して確率的に組み込むため、より解剖学的に正確で信頼性の高い空間再構成が可能である。同時期に報告されたAchimらの研究 (Platynereis 脳におけるFISHデータとの統合) と並び、本研究は「scRNA-seqと空間参照マップの統合」という新しい解析パラダイムを確立した。
新規性: 本研究で初めて、連続的なscRNA-seq発現量と二値化されたin situパターンを統計的に結びつけるガウス混合モデルを導入し、LASSO回帰を用いた遺伝子発現補完 (imputation) によりscRNA-seq特有の技術的ノイズやドロップアウトの影響を劇的に軽減した新規な確率的フレームワークを開発した。さらに、細胞を単一のbinに強制的に割り当てるのではなく、事後確率分布として表現するソフトアサインメント (soft assignment) を採用し、マッピングの不確実性を許容した点、および空間情報を考慮した差分的発現解析 (spatially aware marker selection) を提案し、組織内の局在文脈に依存した新規マーカー同定を可能にした点もこれまで報告されていない。特に、この空間情報を加味した解析アプローチは、従来の単純な細胞群比較における偽陽性問題を解決する極めて強力な手段である。
臨床応用: 本研究はゼブラフィッシュの初期胚を実証モデルとしているが、Seuratのアルゴリズム自体は生物種や組織型を問わない高い汎用性を備えており、bench-to-bedside のトランスレーショナルリサーチ文脈において、本手法は以下のような臨床応用が期待される。(i) がん組織における腫瘍細胞と免疫細胞の空間的相互作用 (腫瘍微小環境) の解明、(ii) 脳腫瘍や転移性病変における解剖学的ニッチ特異的な遺伝子発現プロファイルの同定、(iii) 治療抵抗性を示すがん亜細胞集団 (subclone) の組織内局在パターンの同定。現在、SeuratはscRNA-seq解析における世界的なデファクトスタンダードツールとして広く普及しており、FFPE (ホルマリン固定パラフィン包埋) 組織や臨床生検検体から得られるシングルセルデータと空間情報の統合解析において、疾患病態の理解や個別化医療の推進に大きく貢献している。
残された課題と今後の方向性: 一方で、本手法にはいくつかの制約が存在する。第一に、腫瘍組織のように患者間で空間配置の再現性がない組織や、空間構造が破綻している病変への適用が困難である点。第二に、成体網膜のように、極めて類似した発現プロファイルを持つ多数の細胞型が組織全体に広く均一に分散している場合、十分なマッピング解像度が得られない点。第三に、高精度なマッピングを行うためには、事前に数十個の適切なランドマーク遺伝子を選定し、そのin situパターンをマニュアルでキュレーションする必要がある点である。これらの課題に対し、著者らは今後の検討課題として、(i) scRNA-seqデータから最適なランドマーク遺伝子を反復的に自動探索するアルゴリズムの構築、(ii) クライオセクションを用いた低入力RNA-seq (cryosection RNA-seq) などのバイアスのない空間トランスクリプトミクスデータとの直接的な統合、(iii) 定量的な蛍光in situデータへの対応と混合モデルの多値化拡張、を今後の方向性として提示している。
方法
Seurat アルゴリズムの構築:
Seuratは以下の4つのステップから構成される計算アルゴリズムである。
(i) データ補完 (Imputation): scRNA-seqデータに特有のドロップアウトや技術的ノイズを軽減するため、高変動遺伝子の共発現パターンを利用する。L1正則化を課した LASSO (least absolute shrinkage and selection operator) 回帰モデル (Tibshirani 1996) を用いて、他の高変動遺伝子の発現量から各ランドマーク遺伝子の発現量を予測・補完する。
(ii) 混合モデルによる二値化の対応付け: 補完された連続的な発現値を、空間参照マップの二値情報 (“on” または “off”) に対応付けるため、期待値最大会 (EM) アルゴリズムを用いて二峰性ガウス混合モデル (Gaussian mixture model) をフィッティングする。
(iii) 空間モデルの構築: 各空間binにおけるランドマーク遺伝子の結合確率分布を多変量正規分布 (multivariate normal distribution) としてモデル化する。任意で共分散パラメータを推定する定量的な精緻化ステップ (quantitative refinement step) を導入する。
(iv) 事後確率の計算: ベイズの定理に基づき、各細胞が各空間binに由来する事後確率を算出し、最も確率の高い領域へのマッピングや、複数領域にまたがる確率分布 (soft assignment) を出力する。本手法はRパッケージ Seurat として実装された。
実験プロトコルとデータ取得: 野生型ゼブラフィッシュ胚を50%エピボリー期 (受精後約6時間) まで育成し、1 mg/ml プロナーゼ処理および物理的解離により単一細胞懸濁液を調製した。Random Molecular Tag (RMT, 5塩基のランダム塩基) を含む改良型 SMART (Switching Mechanism at 5’ end of RNA Template) テンプレートスイッチングプロトコルを用いて5’末端scRNA-seqライブラリを調製した。Illuminaシーケンサーを用いて平均 530,000 reads/cell の深度でシーケンシングを行った。1,152個の単一細胞からライブラリを調製し、検出遺伝子数 2,000 個以上の基準を満たした 945 細胞を保持した。さらに、PCAにより包皮層である EVL (enveloping layer) 細胞 (krt18/krt4/cldne陽性) 94 細胞を同定して除外し、最終的に 851 細胞を解析対象とした。
空間参照マップの構築と検証: ZFIN (Zebrafish Model Organism Database) などの公開データベースから、50%エピボリー期における47個のランドマーク遺伝子のカラーin situ hybridization画像を取得した。胚を動物極-植物極軸および背腹軸に沿って128の領域に分割し、左右対称性を考慮して64の空間binに統合した。各binにおける47遺伝子の発現の有無を “on” (1) または “off” (0) としてマニュアルでスコアリングし、空間参照マップを作成した。 実験的検証として、以下の3つの独立したコントロールを用意した。 ① 位置既知の参照細胞 (Reference cells): 胚の特定位置から顕微鏡下で手動で吸引回収した 28 細胞。背腹軸の指標として、sphere stageにおいてdextran-Alexa 488を注入した細胞を移植し、シールド期 (shield stage) における位置を逆算した。 ② マージン濃縮細胞 (Margin-enriched cells): 動物極側をメスで切除し、赤道部 (margin) を濃縮したプロトコルで回収した 141 細胞。 ③ 新規遺伝子のin situ検証: 空間パターンが未解明であった14遺伝子について、予測された発現パターンを検証するため、ジゴキシゲニン (DIG) またはフルオレセインで標識したRNAプローブを用いて比色in situ hybridizationおよび2色蛍光in situ hybridization (FISH) を実施した。 統計的検定および解析の妥当性を検証するため、Rパッケージ内での統計処理として Student t-test や Spearman correlation を用いた相関解析、および多重比較補正として FDR (false discovery rate) 制御を導入した。