GPSNorm: Gaussian Process Spatial Normalization for Spatial Transcriptomics

GPSNormの技術要約：空間トランスクリプトミクスにおけるガウス過程空間正規化

問題提起

10x Genomics VisiumやNanoString GeoMxといった空間トランスクリプトミクス（ST）技術は、組織構造を維持しながらゲノムワイドな遺伝子発現測定を可能にする。しかし、STデータは、ライブラリサイズの差異、スライドレベルのバッチ効果、および空間的なアーティファクト（染色性の変化や組織の厚みの違いなど）を含む技術的な変動に大きく左右される。バルクまたはシングルセルRNAシーケンシングとは異なり、STにおけるこれらの技術的要素はしばしば空間的に構造化されており、生物学的に意味のある空間信号と区別することが困難である。

現在の正規化手法は、通常、正規化を逐次的な前処理ステップとして扱う。すなわち、カウントはグローバルまたは空間的に情報を得たスケーリング係数を用いて調整され、その後の解析（例：微分発現解析）は、これらの調整された値を固定された入力として実行される。この二段階のワークフローは、正規化中に導入される不確実性を無視している。技術的成分と生物学的成分が部分的に分離可能な場合、正規化された値を固定値として扱うことは、効果量の推定にバイアスをもたらし、下流の統計解析において反保守的な推論を招く可能性がある。

手法：GPSNormフレームワーク

著者らは、GPSNorm（Gaussian Process Spatial Normalization：ガウス過程空間正規化）を提案する。これは、技術的変動、空間構造、および生物学的信号を統一されたモデル内で同時に推定するベイズ階層モデルである。

1. 階層モデルの仕様:

• データモデル: 遺伝子 $g$ および空間単位 $i$ における観測された転写物カウント $y_{gi}$ は、対数リンク関数を用いた負の二項分布（Negative Binomial distribution）を用いてモデル化される。
• 線形予測子: 対数平均は以下のように分解される：
  • オフセット: 既知のライブラリサイズ／シーケンシング深度 ($L_i$)。
  • 生物学的効果: 遺伝子固有のベースライン ($\alpha_g$) および関心のある共変量に対する回帰係数 ($\beta_g$)（例：疾患状態）。
  • 技術的効果: スライドレベルの効果 ($\gamma_s$) およびスライド固有の空間効果 ($f_{i,s}$)。
• コントロール遺伝子: モデルは、ネガティブコントロール遺伝子（NEG）およびハウスキーピング遺伝子（HK）を利用する。ネガティブコントロールについては、遺伝子固有の生物学的パラメータをゼロに制約することで、生物学的信号の干渉を受けることなく、共有される技術的変動と分散の推定を固定（アンカー）することができる。

2. 空間モデリング:

• 各スライド内の空間依存性は、**固有条件付き自己回帰（ICAR）**事前分布を用いてモデル化される。これは、確率的偏微分方程式（SPDE）による連結を介した、滑らかなガウス過程の疎な離散近似として機能する。
• 空間場は、空間単位の近傍構造から導出されたグラフ・ラプラシアンによって定義されるガウス・マルコフ・ランダム場（GMRF）として表現される。これにより、密なガウス過程共分散行列と比較して計算のスケイラビリティが確保される。

3. 推論:

• 本モデルは潜在ガウスモデル（LGM）のクラスに属する。
• 効率的な近似ベイズ推論のために、**統合ネスト・ラプラス近似（INLA）**が使用される。INLAは、ICAR事前分布の疎な精度行列構造を利用し、ハイパーパラメータに関する決定論的な積分を行うことで、マルコフ連鎖モンテカルロ法（MCMC）の計算コストを回避する。
• 不確実性の伝播: 二段階の手法とは異なり、GPSNormは、微分発現パラメータ ($\beta$) を推定する際に、潜在的な正規化成分 ($u$) の不確実性を統合する。これにより、正規化と生物学的信号が混在する場合でも、より広く、較正された（calibrated）信頼区間が得られる。

4. 正規化出力:
本フレームワークは、事後推定値をカウントスケールに投影する2つの方法を提供する：

• 事後中央値投影（Posterior-Median Projection）: 線形予測子の生物学的成分を用いて、負の二項分布の中央値を算出する。
• パーセンタイル調整カウント（Percentile-Adjusted Counts: PAC）: ランクと分散を保持しつつ、技術的効果を除去するために平均構造を調整しながら、観測されたカウントを変換する。

主な貢献

1. 結合モデリング: GPSNormは、正規化を単なる前処理ステップではなく、階層ベイズフレームワーク内での結合推定問題として再定義し、正規化の不確実性を下流の推論へと明示的に伝播させる。
2. 空間的技術モデリング: 生物学的な勾配と混同されやすい、滑らかで空間的に構造化された技術的アーティファクトを捉えるために、スライド固有のICAR事前分布を導入した。
3. スケーラビリティ: INLAとGMRFの疎な構造を利用することで、本手法は、密なガウス過程の逆行列計算に伴う立方体の計算量（cubic complexity）を伴わずに、大規模な空間グラフ上での効率的な推論を可能にする。
4. コントロールの統合: ネガティブコントロールおよびハウスキーピング遺伝子を形式的に組み込むことで、技術的変動の推定を固定し、技術的成分と生物学的成分の間の混同を軽減する。

結果

著者らは、経験的なデータ分布に基づいたシミュレーション研究、および3つの実世界のデータセットへの適用を通じて、GPSNormを評価した。

1. シミュレーション研究:

• 技術的回復: GPSNormは、強い生物学的混同（微分発現と空間信号の重複が70%）がある条件下でも、シミュレーションされた空間的技術場を正確に回復した（相関係数 > 0.92）。
• 微分発現（DE）性能: 既存の手法（例：scran, Giotto, SpaNorm, scd）と比較して、GPSNormは以下の特性を示した：
  • 対数フォールドチェンジ推定における最小限のバイアス。
  • 低減された平方根平均二乗誤差（RMSE）。
  • 特に技術的複雑性が高いシナリオにおいて、DE遺伝子のランキング性能（AUC）の向上。
• 不確実性の較正: 真のDE効果がゼロである簡略化されたシミュレーションにおいて、結合推定は、強い空間的混同の下でも公称95%信頼区間カバレッジを維持したが、二段階のアプローチは反保守的なカバレッジ（0%に近い値）を示した。

2. 実データへの適用:

• 10x Visium Human DLPFC: GPSNormは、ライブラリサイズの不均一性によって隠されていた白質領域におけるオリゴデンドロサイトマーカー MOBB の期待される濃縮を正常に回復した。このパターンは、非空間的手法では完全に見落とされていた。また、これらの空間パターンの事後不確実性推定も提供した。
• NanoString GeoMx COVID-19 Lung Damage: 小規模サンプル研究（46 AOI）において、GPSNormは重度の損傷に関連する生物学的に関連のあるケモカイン（CCL18, CXCL10, CCL2）および軽度の損傷に関連する SCGB1A1 を特定した。他の手法は、小規模なコホートにおける正規化の不確実性を伝播できないため、有意な遺伝子を特定できなかった。
• GeoMx Spatial Organ Atlas Kidney: GPSNormは、腎臓の糸球体および近位尿細管コンパートメントの両方について、精査されたマーカー遺伝子の高い優先順位を維持し、解剖学的領域間の期待される転写コントラストを保持する点で他の手法を上回った。

3. 計算のスケーリング:
著者らは、モデル適合の前に関心領域（AOI）をクラスター化することで、高シグナル遺伝子のランキングと効果量の大きさを維持しつつ、実行時間を大幅に短縮（最大5倍の高速化）できることを示した。

意義と主張

本論文は、正規化を決定論的な前処理ステップではなく、潜在的な統計的推論問題として扱うことが、空間トランスクリプトミクス解析の堅牢性と解釈性を大幅に向上させると主張している。

• 堅牢性: 技術的効果と生物学的効果を結合してモデル化することにより、GPSNormは、正規化中に生物学的に意味のある変動を吸収するリスクを軽減し、下流の微分発現における不確実性の過小評価を防ぐ。
• 解釈性: 本フレームワークにより、研究者は空間発現表面および正規化調整の不確実性を定量化することができ、実験的変動のより完全な全体像を提供できる。
• 適用可能性: 本手法は、技術的および生物学的な空間構造が絡み合っている場合や、サンプルサイズが小さい場合（GeoMxのようなターゲット型空間プラットフォームなど）に特に有利である。階層構造により、遺伝子間および空間単位間での情報の共有が可能となる。

著者らは、限界事項として、極めて大規模なデータセットに対する計算コスト（AOIクラスタリングにより緩和可能）や、技術的効果に対する滑らかな空間事前分布の仮定を挙げている。今後の展望として、本フレームワークを空間的に変動する遺伝子の検出と結合させ、シングルセルリファレンスデータを統合することなどを挙げている。