心血管疾患に関連する変動パターンと分子プロセスの解明のための教師なしマルチオミクス因子解析の応用

サマリー

複雑なマルチオミクスデータセットの教師なし分析のための柔軟で拡張可能なJupyterラボベースのワークフローを紹介します。これは、さまざまな前処理ステップ、マルチオミクス因子分析モデルの推定、およびいくつかのダウンストリーム分析を組み合わせたものです。

要約

病気のメカニズムは通常複雑で、いくつかの異なる分子プロセスの相互作用によって支配されています。複雑な多次元データセットは、これらのプロセスに関するより多くの洞察を生成するための貴重なリソースですが、このようなデータセットの分析は、たとえば、さまざまな病状、時点、さまざまな解像度でプロセスをキャプチャするオミクスなど、高次元性のために困難な場合があります。

ここでは、急性および慢性冠症候群の免疫応答を捕捉する血液サンプルから生成されたデータセットにマルチオミクス因子分析(MOFA)を適用することにより、このような複雑なマルチオミクスデータセットを教師なしの方法で分析および探索するアプローチを紹介します。このデータセットは、サンプルレベルのサイトカインデータ、血漿プロテオミクスおよび好中球prime-seq、シングルセルRNA-seq(scRNA-seq)データなど、異なる分解能の複数のアッセイで構成されています。患者ごとにいくつかの異なる時点といくつかの患者サブグループを測定することで、さらに複雑さが増します。

分析ワークフローでは、(1)データの前処理と調和、(2)MOFAモデルの推定、(3)ダウンストリーム分析など、いくつかのステップでデータを統合および分析する方法を概説しています。ステップ 1 では、さまざまなデータ型の特徴を処理し、低品質の特徴を除外し、それらを正規化して分布を調和させ、さらに分析する方法について概説します。ステップ2では、MOFAモデルを適用し、すべてのオミクスと特徴にわたってデータセット内の分散の主な原因を調査する方法を示します。ステップ3では、捕捉されたパターンの下流分析のためのいくつかの戦略を示し、それらを疾患の状態とそれらの状態を支配する潜在的な分子プロセスに関連付けます。

全体として、複雑なマルチオミクスデータセットの教師なしデータ探索のワークフローを提示し、他のコンテキストやマルチオミクスデータセット(例示的なユースケースで提示された他のアッセイを含む)にも適用できる、異なる分子特性で構成される変動の主軸の同定を可能にします。

概要

病気のメカニズムは通常複雑で、いくつかの異なる分子プロセスの相互作用によって支配されています。特定の疾患につながる、または疾患の進化を支配する複雑な分子メカニズムを解読することは、疾患の理解と治療のための新たな洞察を明らかにする可能性があるため、医学的関連性が高い課題です。

近年の技術の進歩により、これらのプロセスをより高い分解能(シングルセルレベルなど)で、またさまざまな生物学的層(DNA、mRNA、クロマチンアクセシビリティ、DNAメチル化、プロテオミクスなど)で同時に測定することが可能になりました。これにより、大規模な多次元生物学的データセットの生成が増加し、それらを共同で分析して、基礎となるプロセスについてより多くの洞察を得ることができます。同時に、生物学的に意味のある方法で異なるデータソースを組み合わせて分析することは、依然として困難な課題^です1。

異なる技術的限界、ノイズ、および異なるオミクス間のばらつきの範囲が異なることが、1つの課題を提起します。例えば、シングルセルRNAシーケンシング(scRNA-seq)のデータは非常にまばらで、大きなテクニカルエフェクトやバッチエフェクトの影響を受けやすいです。さらに、フィーチャスペースは多くの場合非常に大きく、数千の測定された遺伝子またはタンパク質にまたがっていますが、サンプルサイズは限られています。これは、いくつかの病状、交絡因子、時点、および解決策を含む可能性のある複雑な設計によってさらに複雑になります。たとえば、提示されたユースケースでは、シングルセルレベルまたはサンプル(バルク)レベルのいずれかで異なるデータ型が使用可能でした。それに加えて、データが不完全である可能性があり、分析されたすべての被験者がすべての測定値を利用できるとは限りません。

これらの課題により、統合解析を行うことでプロセスの全体像が得られるだけでなく、1つのオミクスからの生物学的および技術的なノイズが他のオミクス^3,4によっても補償される可能性があるにもかかわらず、異なるオミクスと含まれる特徴は^{依然として別々にのみ}解析されることがよくあります²。マルチオミクスデータの統合解析を行うためには、ベイズ法、ネットワークベース法^5,6、マルチモーダル深層学習⁷、行列因数分解^8,9による次元削減法など、いくつかの異なる方法が提案されている。後者については、大規模なベンチマーク研究¹⁰の結果により、MOFA⁹(マルチオミクス因子分析)法は、データを臨床アノテーションにリンクする必要がある場合により適したツールの1つであることが示されています。

特に複雑な設定では、教師なし行列因数分解法は、複雑さを軽減し、さまざまなデータソースや特徴から共有シグナルと相補的なシグナルを抽出するための有用なアプローチです。複素空間を低ランクの潜在表現に分解することにより、データ内の分散の主要な原因を迅速に探索し、既知の共変量に関連付けることができます。同じパターンの変動が複数の特徴(遺伝子やタンパク質など)で共有されている場合、ノイズが減少する一方で、これはいくつかの要因に集約される可能性があります。正則化を使用してモデル係数のスパース性を増加させることができるため、このアプローチは、特徴空間が大きく、サンプル数が限られている設定に適しています⁹。

このプロトコルは、MOFAモデルを使用して、複雑なマルチオミクスデータセットを迅速に探索し、このデータセットを特徴付ける主要な変動パターンを抽出する方法を示す柔軟な解析ワークフローを提供します。ワークフローは、主に 3 つのステップで構成されています。最初のステップである「データの前処理と調和」では、さまざまな入力データタイプ(scRNA-seq、プロテオミクス、サイトカイン、臨床データ)に基づくデータ前処理のさまざまな戦略が提示されます。このプロトコルでは、さまざまな入力データセットの特徴を処理し、低品質の特徴を除外し、それらを正規化して分布を調和させる方法について詳しく説明します。また、これらの前処理の決定が下流の結果にどのように影響するかについても説明します。2番目のステップでは、MOFAモデルをデータに適用し、結果として生じる分散分解を使用して、異なるデータセットの統合を評価できます。3番目のステップでは、捕捉した因子を共変量にリンクし、それらの因子を定義する分子プログラムを明らかにする方法を示します。提示されたワークフローにより、冠症候群に罹患している患者のデータセットで臨床共変量に関連するいくつかの潜在因子を抽出し、以前のプロジェクト¹¹から潜在的な基礎となる多細胞免疫プログラムを特定することができました。ここではこのデータセットを使用しますが、このプロトコルは他のオミクスを含む他のコンテキストにも簡単に適用できます。

このデータセットは、安定した慢性冠症候群(CCS)、急性冠症候群(ACS)、および健康な冠動脈症候群(非CCS)の対照群からのサンプルで構成されています(図1)。ACSは、既存のCCSのプラーク破裂によって引き起こされ、心筋への血流の急性の混乱とそれに続く心臓の虚血性損傷につながります。この損傷は、免疫系による炎症反応を引き起こし、その後修復段階を引き起こし、それは急性イベントの数日後まで続きます¹²。ACS患者のこの免疫応答を特徴付けることができるように、血液サンプルは4つの異なる時点で採取されました:急性(TP1);再開通後(14 [± 8] h)(TP2);60 [± 12] h後(TP3);退院前(6.5 [±1.5]日)(TP4)(図1A)。CCSおよび健康な冠状動脈の患者の場合、利用可能なタイムポイントは1つだけでした-(TP0)。すべての患者と時点について、血液サンプルに基づくさまざまなアッセイが測定されました:炎症の臨床マーカー(クレアチンキナーゼ(CK)、CK-MB、トロポニン、C反応性タンパク質(CRP))、末梢血単核細胞(PBMC)のscRNA-seq、サイトカイン分析、血漿プロテオミクス、好中球のprime-seq¹³ データ。

図1:心筋梗塞マルチオミクス入力データセット。 入力データセット: 分析されたデータには、急性冠症候群 (ACS)、慢性冠症候群 (CCS)、および健康な冠動脈 (非 CCS) の患者 (n = 62) の血液サンプルが含まれます。ACS 患者の場合、血液サンプルは 4 つの異なる時点 (TP1-4) に含まれ、CCS 患者と非 CCS 患者の場合は 1 つの時点 (TP0) に含まれていました。各患者と時点の組み合わせは、分析で個別のサンプルとして扱われます。サンプルでは、臨床血液検査(n = 125)、scRNA-seq(n = 121)、血漿プロテオミクス(n = 119)、サイトカインアッセイ(n = 127)、好中球prime-seq(n = 121)など、さまざまなオミクスアッセイがサンプルで測定されました。その後、記載されたプロトコルを適用して、すべてのオミクスにわたるデータを統合し、MOFAモデルを用いて探索し、さらに下流の分析(因子分析、パスウェイエンリッチメント)を行いました。 この図の拡大版を表示するには、ここをクリックしてください。

ここで紹介するワークフローの入力として、例えばscanpy¹⁴ の前処理チュートリアルで概説されているように、cellrangerおよび品質管理(QC)で処理した後のscRNA-seqデータから生のカウントを取ります。細胞タイプのアノテーションには、自動化されたAzimuth¹⁵ パイプラインを使用しました。次に、カウントは、各サンプルとセルタイプのすべてのセルの平均を取ることにより、各セルタイプのサンプルレベルで集計されます(疑似バルク集計)。プラズマプロテオミクスは正規化された強度と中央中心の強度として含まれ、好中球については、prime-seqからumi unique molecular identifier(UMI)エクソンカウントを取得します。サイトカインおよび臨床値については、以前の前処理は適用されていません。(実験的)データ生成に関するさらなる詳細は、対応する原稿¹¹に概説されている。ここで紹介した結果は、scRNA-seqデータ中の細胞型に対して自動Azimuthアノテーションを使用した結果に基づいており、参照された論文で使用されたマーカーベースの戦略と比較した結果であるため、ここで提示された結果は、論文で提示されたものと類似していますが、まったく同じではありません。この原稿では、細胞型アノテーション戦略は分析の主要なパターンと生物学的解釈を変更しないが、モデルから生じる正確な値のわずかな変化は異なる可能性があることを示すことができます。全体として、入力データは、10,000を超える異なる特徴(遺伝子、タンパク質、臨床値)のさまざまな時点と測定レベル(単一細胞とバルク)を含む複雑な多次元データセットでした。MOFA分析に続く厳密な前処理とデータ調和戦略は、データを探索し、関連する免疫プログラムを抽出するための有用で迅速なツールであることが示されています。各時点と患者の組み合わせは、MOFA 分析では独立したサンプルとして扱われます。各データ型とセルの種類は、MOFA 分析では個別のビューと見なされます。

このプロトコルでは、ワークフローの入力データの準備、さまざまなワークフローステップの実行、構成のカスタマイズ、結果の数値の解釈、および解釈に基づく構成の反復調整を行う手順を提供します。プロトコルのさまざまなステップの概要、各ステップで必要な入力データセット、および結果として得られる図とデータセットは、テクニカルワークフローの概要に記載されています(図2)。

図2:テクニカルワークフローの概要。 マルチオミクスデータセットの解析ワークフローの概要。さまざまな要素が、さまざまな色と記号で強調表示されます。データの前処理と調和 (1) のステップに属する Jupyter Notebook は青色で色分けされています。「MOFA Model」(2)ステップに属するJupyter Notebookはオレンジ色で色づいています。「Downstream Analysis」(3)ステップに属するJupyter Notebookは緑色で色付けされています。結果の比較に使用する 1 つの Jupyter Notebook は、黄色で色付けされています。ワークフローの実行パラメータを変更できる設定ファイルは、紫色で強調表示されます。ワークフローの実行に必要な入力データセットは、データセット記号で示され、グレーでハイライト表示されます。ワークフローの実行中に生成されるすべての Figure 出力は、虫眼鏡記号で示されます。ワークフローの実行中に生成されたデータセットは、テーブルとして示されます。一般に、ワークフローは順番に実行されます:(1)データの前処理と調和は、scRNA-seq入力データに基づく疑似バルクテーブルの第一生成(01_Prepare_Pseudobulk)と、その後のこのデータと他のすべてのサンプルレベル(バルク)入力(02_Integrate_and_Normalize_Data)の統合と正規化の2つのステップで構成されます。このステップでは、設定ファイルを使用して、指定された前処理ステップと正規化ステップ(サンプルフィルターなど)をデータセットごとに個別に構成することができます。(2) 「MOFAモデル」:設定ファイル(03_MOFA_configs.csv)で指定された構成を使用して、最初のステップで生成された入力に対してMOFAモデルを実行します (3) 「ダウンストリーム分析」:生成されたMOFA結果に対する洞察を生成し、それらを「Sample Meta Data.csv」ファイルを介して入力として提供されるサンプルメタデータ(共変量)に関連付けるために、互いに独立して実行できる3つの異なるノートブックで構成されています。(4) 「モデル比較」:ステップ2で生成された異なるモデルを比較するために使用できる小さな別個のステップです。 この図の拡大版を表示するには、ここをクリックしてください。

ワークフローは、R と Python で記述された複数の Jupyter Notebook で構成されています (ワークフローの実行に R と Python 言語の知識は必要ありませんが、エラーが発生した場合に役立つ場合があります)。プロトコルのさまざまなステップで、パラメータは設定ファイル(名前に接尾辞「_Configs」を含む「.csv」ファイル)を介して変更されます。プロトコル内では、デフォルト設定から変更する必要があるパラメータのみを概説します。

他のいくつかのパラメータも、たとえば前処理をカスタマイズするために変更できます。これらのパラメータと説明のドキュメントは、ダウンロードしたリポジトリに含まれているファイル 'Documentation_Config_Parameter' に記載されています。

プロトコル

1.準備:技術的なセットアップとインストール

注: このプログラムを実行するには、wget、git、および Apptainer がデバイスにプリインストールされています。さまざまなシステム(Linux、Windows、Mac)にApptainerをインストールするためのガイドは、ここにあります:https://apptainer.org/docs/admin/main/installation.html。git のインストール情報は、https://git-scm.com/book/en/v2/Getting-Started-Installing-Git にあります。さまざまな入力データセットのサイズに応じて、適切なマシン(16 CPU、64 GB メモリ)でワークフローを実行することをお勧めします。提供されたサンプルデータを使用したスモークテストは、ローカルマシンで実行できます。サンプル データでプロトコルを実行した場合の手順と予想される出力は、補足ファイル 1 に記載されています。上記のデータセットで実行されるプロトコルの重要な手順については、 補足ビデオファイル1 を参照してください。

1. コンソールを開き、すべての解析コードと出力を保存するフォルダを選択または作成します。ターミナルでコマンドcd path_to_folderを入力して、フォルダに移動します。
2. Github(https://github.com/heiniglab/mofa_workflow)からコードリポジトリをダウンロードまたはクローンするには、ターミナルウィンドウにgit clone https://github.com/heiniglab/mofa_workflow.git と入力します。
3. 必要なすべてのインストールを含むイメージを Zenodo からダウンロードするには、ターミナル ウィンドウに 「wget https://zenodo.org/records/11192947/files/mofa_image.sif」と入力します。
4. ターミナルウィンドウに 「mkdir results 」と入力して、すべての結果データが保存されるフォルダを生成します。
5. 解析で使用するすべての入力データを追加するフォルダを生成するには、ターミナルウィンドウに mkdir input_data と入力します。
6. JupyterLab セッションを開始するコンテナを実行するには、ターミナルで次のコマンドを入力します。 apptainer run mofa_image.sif。コマンドによって返された URL をブラウザーにコピーすると、Jupyter-lab セッションが開きます (Jupyter-lab の詳細については、ソフトウェアのドキュメント¹⁶ を参照してください)。
   注: ワークフローがラップトップ上でローカルに実行される場合は、代わりに コマンド apptainer exec mofa_image.sif jupyter-lab を使用することをお勧めします。これにより、ローカル ホスト アドレスが直接返されます。コンテナがクラスタ化されたコンピューティング環境内で実行される場合、 ssh 経由で実行できるポートフォワーディングの設定が必要になることがあります。

2. 初期化とデータ準備

1. Jupyter-Lab セッションでは、左側のナビゲーション メニューを使用します。input_dataをダブルクリックしてinput_dataフォルダに移動します。
2. 解析の入力として使用するすべてのデータセットを、 ドラッグ&ドロップを使用してinput_dataディレクトリにコピーします。現在配置されているフォルダからファイルをドラッグし、Jupyter-lab セッションの input_data フォルダの下の領域にドロップします。
   注: すべてのデータセットは、.csv または .h5ad (シングルセル データの場合) 形式である必要があります。すべての.csvファイルには、一致する sample_id 列が含まれている必要があります (データセット全体で同一の ID を使用する必要があります)。他のすべての列はフィーチャとして使用されます。h5ad- ファイル内では、セル注釈には、sample_idとcluster_idを指定する 2 つの識別子が含まれている必要があります。これらは、集計と照合に使用されます。他の形式のオミクスデータセットは、使用する前に指定された.csv形式に変換する必要があります(図3)。 .h5seurat 形式で与えられた scRNA-seq データセットは、Jupyter-notebook: 00_Data_Conversion.ipynb を実行して .h5ad に変換できます。
3. フォルダの記号をクリックし、フォルダmofa_workflow、スクリプト、および構成をダブルクリックして、構成フォルダに移動します。フォルダ内で、ファイルData_configs.csvをダブルクリックして開きます。
4. [値] 列に、input_data (data_path) フォルダーと結果 (result_path) フォルダーのフォルダーへのパスを追加します。configuration_nameの値列に、保存されているすべてのファイルにファイル拡張子として追加される名前を追加します (このプロトコルは [Myocardial Infarction version1] MI_v1使用されます) (図 4)。
5. 上部のメニューで[ ファイル]>[CSVファイルを保存 ]をクリックして、変更を保存します。
6. 左側のナビゲーションメニューを使用して、scriptsをクリックしてscriptsフォルダに移動します。00_Configuration_Update.ipynb をダブルクリックして、初期化ノートブックを開きます。スクリプトを実行するには、上部にある[Restart kernel and run all cells]ボタンをクリックし、ポップアップで[Restart]をクリックします(図5)。

図3:データ入力とセットアップ。 ワークフローを実行するには、すべてのデータを指定したinput_dataフォルダーに保存する必要があります。入力データセットごとに、個別のファイルを用意する必要があります。シングルセルデータは、cluster_id上の細胞アノテーション(例えば、以前の細胞タイプのアノテーションステップから得られる)とsample_id列(分析すべき各サンプルを一意に識別する)を含む.h5adとして指定する必要があります。他のすべての入力データセットは、sample_idを指定する1つの列(単一セルデータの対応する列と一致する)と、他のすべての列のMOFA分析で使用される特徴を含む「.csv」形式で指定する必要があります。 この図の拡大版を表示するには、ここをクリックしてください。

図4:Jupyter-lab設定ファイル。 ワークフローの実行中、パラメータの変更(フィルタリングオプションの調整など)は、「.csv」設定ファイルを介して指定されます。クローニングされたリポジトリ内には、各ステップのデフォルトの設定ファイルが含まれています。これらは、スプレッドシートと同様に、jupyter-lab コンソールで直接編集できます。 この図の拡大版を表示するには、ここをクリックしてください。

図5:Jupyter-notebooksスクリプト。 完全なワークフローは、対応する設定ファイルが変更された後に順番に実行される一連の Jupyter Notebook で構成されています。左側のJupyterノートブックをダブルクリックすると、対応するファイルが右側で開きます。ファイルの完全な実行は、上部で強調表示されているボタンから開始できます。 この図の拡大版を表示するには、ここをクリックしてください。

3. データの前処理と調和

1. 前処理 - scデータを擬似バルクに変換します。
   注:この手順は、分析でシングルセルデータを使用する場合にのみ実行する必要があります。
   1. 左側のナビゲーション メニューを使用して、 構成をダブルクリックして構成フォルダーに移動します。ダブルクリックして 01_Preprocessing_SC_Data.csv ファイルを開きます。ファイルに自動的に入力された値を確認し、必要に応じて、[ data_name ] 列の値を、解析に使用される input_data フォルダー内の単一セル データセットのファイル名に対応するように調整します。
      注 : デフォルトでは、入力データ フォルダ内のすべての .h5ad ファイルの名前が初期化スクリプトの設定ファイルに追加されます。一部のデータセットを分析に使用しない場合は、ここで削除できます。
   2. 上部のメニューで[ ファイル]>[CSVファイルを保存 ]をクリックして、変更を保存します。
   3. 左側のナビゲーションメニューを使用して、scriptsをクリックしてscriptsフォルダに移動します。ノートブック 01_Prepare_Pseudobulk.ipynb をダブルクリックして開きます。スクリプトを実行するには、上部にある[カーネルを再起動してすべてのセルを実行]ボタンをクリックし、ポップアップで[再起動]をクリックします。
   4. 左側のナビゲーション メニューを使用して figures フォルダーに移動するには、最初に figures をダブルクリックしてから 01_figuresをダブルクリックします。新しく生成されたプロット FIG01_Amount_of_Cells_overview をダブルクリックして開きます。
      注: ノートブックの実行には数分かかる場合があります。 ノートブックが正常に実行されると、ポップアップが表示され、ファイル FIG01_Amount_of_Cells_Overview ノートブックによって更新されるか、新しく生成されます。 [最終変更日 ] 列は、ファイルが新しいファイルか古いファイルかを評価するために、ファイルが生成された日時を示すことができます。
   5. プロットを調査して、サンプルあたりの細胞数が非常に少ない細胞型クラスターを特定します。これらのcluster_idsの名前を書き留めて、後続の手順で除外します(図6)。
   6. 左側のナビゲーションメニューを使用して、をクリックして構成フォルダーに戻ります...次に、[ 構成]をダブルクリックします。ファイル 02_Preprocessing_Configs_SC.csv をダブルクリックして開きます。
   7. configuration_name列とdata_name列の値を確認し、必要に応じて調整します。
      注 : 初期化スクリプト内では、これらの値には、入力データ フォルダ内のすべての .h5ad ファイルの名前と、以前に Data_Configs.csv ファイル内で設定された configuration_name 値が事前に入力されています。ファイルを分析から除外する必要がある場合、またはファイル名の別の拡張子を使用する必要がある場合は、ここで調整できます。
   8. cell_type_exclusion列の値を調整し、前の手順で除外するように特定されたすべてのcluster_idを「,」で区切って追加します。
   9. 上部のナビゲーションバーで[ ファイル]>[CSVファイルを保存 ]をクリックして、変更を保存します。
2. 前処理 - 他のオミクスデータソースを調和させ、統合します。
   1. ファイル 02_Preprocessing_Configs.csv をダブルクリックして開き、 data_input フォルダーに含まれて保存される各データセットの前処理構成を調整します (データセットごとに 1 行)。
   2. configuration_name列とdata_name列の値を確認し、必要に応じて調整します。
   3. 列内の他のパラメータは、適用する前処理ステップに応じて、それに応じて調整します。
      注: デフォルト値は、input_dataset フォルダ内にあるデータセットごとに追加されますが、データの 1 つのデータ型に固有ではありません。そのため、調整が必要になります。パラメータの詳細なドキュメントは、 Documentation_Config_Parameter.doc ファイルに記載されています。
   4. [ ファイル] > [CSV ファイルの保存] をクリックして、変更を保存します。
   5. 左側のナビゲーションメニューを使用して、scriptsをクリックしてscriptsフォルダに移動します。ノートブック 02_Integrate_and_Normalized_Data_Sources.ipynb をダブルクリックして開きます。スクリプトを実行するには、上部にある[カーネルを再起動してすべてのセルを実行]ボタンをクリックし、ポップアップで[再起動]をクリックします。
   6. 左側のナビゲーション メニューを使用して、生成された 02_results フォルダーに移動するには、 フォルダー 記号をクリックし、 結果 と 02_resultsをダブルクリックします。結合された前処理されたデータ入力ファイルを含むファイル 02_Combined_data_'configuration_name'_Integrated.csv が含まれているかどうかを確認します。

図6:データの前処理と調和 '01_Prepare_Pseudobulk' ステップの出力の 1 つは、プロット 'Fig01_Amount_of_Cells_Overview' です。ここでは、各cluster_id(y軸は前の細胞型アノテーションステップの細胞タイプを示す)について、サンプルあたりの細胞数('sample_id')が与えられます。提示された結果の中で、サンプルあたりの細胞量が少ない細胞タイプは、その後の分析から除外されます(取り消し線で示されています)。 この図の拡大版を表示するには、ここをクリックしてください。

4. MOFAの実行

1. Jupyter-Lab では、左側のナビゲーション メニューを使用して、 フォルダー 記号をクリックしてから mofa_workflow をダブルクリックし、 次にスクリプト と 構成をダブルクリックして、構成フォルダーに移動します。ファイル 03_MOFA_Configs.csv をダブルクリックして開きます。
2. configuration_name列とmofa_result_name列のエントリを確認し、別の名前を使用する必要がある場合はエントリを調整します。
   注: mofa_result_name は、MOFAに基づいて生成されたすべての結果ファイルにファイル拡張子として追加されます。これは、同じ入力データで異なるMOFAセットアップが実行される可能性があるため、 configuration_name 値とは異なる場合があります(このプロトコルは MI_v1_MOFAを使用します)。
3. MOFAモデルで推定する係数の量(amount_of_factors 列)を入力し、ファイル内の値を調整して、重み付けとスケーリングを適用するかどうか(weighting_of_views 列と scale_views 列)を定義します。
4. [ ファイル] > [CSV ファイルの保存] をクリックして、変更を保存します。
5. 左側のナビゲーションメニューを使用して、「scripts」をクリックしてscriptsフォルダに移動します。ノートブック 03_Run_MOFA.ipynb をダブルクリックして開きます。スクリプトを実行するには、上部にある[カーネルを再起動してすべてのセルを実行]ボタンをクリックし、ポップアップで[再起動]をクリックします。
6. 図をダブルクリックしてから03_figuresをクリックして、03_figuresフォルダに移動します。生成されたプロットFIG03_Overview_Variance_Decomposition_を開きmofa_result_nameモデル結果を調査します(図7A)。
7. 左側のナビゲーション メニューを使用して、生成された 03_results フォルダーに移動するには、 フォルダー 記号をクリックし、 results と 03_resultsをダブルクリックします。サンプル係数値ファイル 03_Factor_Data_'mofa_result_name'.csv と特徴因子重みファイル 03_Weight_Data_'mofa_result_name'.csv が含まれているか確認してください。

5. ダウンストリーム解析

1. 因子の解釈。
   1. 左側のナビゲーションメニューを使用して、フォルダ記号をクリックし、次にinput_dataをダブルクリックしてinput_dataフォルダに移動します。
   2. 生成された因子に関連して分析されるサンプルのすべてのメタデータ(共変量)を含む.csvファイル(Prepared_Sample_Meta_Data.csv)を準備します。ドラッグ&ドロップでinput_dataフォルダにファイルをコピーし、input_dataフォルダのフォルダ概要にファイルをドロップします。
      注: 以前に使用したデータと照合するための sample_id 列と、分析する必要がある各機能の追加の列が含まれている必要があります。
   3. Jupyter-Lab では、左側のナビゲーション メニューを使用して、 フォルダー 記号をクリックし、 次に [mofa_workflow] をダブルクリックし、 次に [スクリプト と 構成] をクリックして、構成フォルダーに戻ります。ファイル 04_Factor_Analysis.csv をダブルクリックして開きます。
   4. configuration_nameとmofa_result_nameのエントリに、スクリプトで分析される構成とMOFAの結果の名前が含まれていることを確認し、必要に応じて調整します。
   5. [numeric_covariates] 列に、MOFA 係数に関連して調査される Prepared_Sample_Meta_Data.csv ファイル内のすべての数値列の名前をコンマで区切って追加します (このプロトコルでは CRP,CK を使用します)。
   6. categorical_covariates' 列に、MOFA 要因に関連して調査される Prepared_Sample_Meta_Data.csv ファイル内のすべてのカテゴリ列の名前をコンマで区切って追加します (このプロトコルでは測定を使用します)。
   7. [ ファイル]>[CSVファイルを保存]をクリックして変更を保存します。
   8. 左側のナビゲーション メニューを使用して、 スクリプトをクリックして 'scripts' フォルダーに移動します。ノートブック 04_Downstream_Factor_Analysis.ipynb をダブルクリックして開きます。上部にある [カーネルを再起動してすべてのセルを実行する ]ボタンをクリックし、ポップアップで[再起動]をクリックして、スクリプトを実行します。
   9. 左側のナビゲーションメニューを使用して、図をダブルクリックしてから04_figuresし、04_figuresフォルダに移動します。生成されたプロットをダブルクリックして開き、興味深いパターンと関連性の要因を調査します。FIG04_Factor_Association_with_numeric_features_
      'mofa_result_name.pdf (図 7B)。FIG04_Factor_Association_
      with_categorical_features_mofa_result_name.pdf (図 7C)。FIG04_Top_Feature_Overview_per_Factor _'mofa_result_name.pdf (図 8A)。
2. 特徴解析
   1. 左側のナビゲーションメニューを使用して、... 次に、[ 構成]をダブルクリックします。ファイル 05_Feature_Analysis_Configs.csv をダブルクリックして開きます。
   2. configuration_name列とmofa_result_name列のエントリが、ダウンストリーム解析に使用されるコンフィギュレーションと生成されたMOFA結果の名前に対応していることを確認し、必要に応じて調整します。
   3. [係数] 列に、次のスクリプト内で上位の特徴がプロットされる係数を追加します。
   4. 列 faceting_variableで、プロット内のサンプルをグループ化するために使用される Prepared_Sample_Meta_Data.csvカテゴリ列の列名を追加します(このプロトコルは 測定を使用します)
   5. [ ファイル]>[CSVファイルを保存]をクリックして変更を保存します。
   6. 左側のナビゲーション メニューを使用して、スクリプト フォルダーに移動し 、 スクリプトをクリックします。ノートブック 05_Downstream_Investigate_Features_Heatmap.ipynb をダブルクリックして開きます。スクリプトを実行するには、上部にある [カーネルを再起動してすべてのセルを実行 ]ボタンをクリックし、ポップアップで [再起動 ]をクリックします。
   7. 左側のナビゲーション メニューを使用して 05_figures フォルダーに移動するには、最初に figures をダブルクリックしてから 05_figures をダブルクリックします。'mofa_result_name' FIG05_Heatmap_Feature_Overview__生成されたプロットを開いて調査します。pdf ファイルをダブルクリックして貼り付けます (図 8B)。
      注: プロットに表示されるフィーチャの量によっては、05_Feature_Analysis_Configs.csv内のパラメータ plot_width と plot_height を調整し、スクリプトを再実行して、すべてがプロットに収まるようにする必要がある場合があります。
3. パスウェイ解析
   1. 左側のナビゲーションメニューを使用して、フォルダ記号をクリックし、次にinput_dataをダブルクリックしてinput_dataフォルダに移動します。
   2. エンリッチメントについてテストするパスウェイのリストを含む .csv ファイル (Prepared_Pathway_Data.csv) を準備します。ドラッグ&ドロップでinput_dataフォルダにファイルをコピーし、input_dataフォルダのフォルダ概要にファイルをドロップします。
      注: ID (パスウェイの一意の識別子)、 gene (パスウェイに属する遺伝子名(SYMBOL)によって与えられる遺伝子、遺伝子ごとに1行)、 pathway_name (パスウェイの名前/テキストによる説明)の3つの列を含める必要があります。
   3. Jupyter-Lab セッションでは、左側のナビゲーション メニューを使用して、...次に、[ 構成]をダブルクリックします。ファイル 06_Pathway_Configs.csv をダブルクリックして開きます。
   4. mofa_result_name列のエントリを確認し、パスウェイエンリッチメントの計算に使用される生成されたMOFA結果の名前に対応していることを確認します。
   5. types 列のエントリを確認し、Prepared_Pathway_Data.csv ファイルの gene 列に一致する特徴を含まない types 列内のエントリを削除します。
      注: デフォルトでは、MOFA モデル内で使用されているすべての異なるビューは、ワークフローの実行中にこのファイルに追加されます。少なくとも 1 つの経路に一致するフィーチャが含まれていないビューがある場合は、それらを削除する必要があります。そうしないと、実行は失敗します。たとえば、パスウェイファイルには遺伝子のパスウェイアノテーションのみが含まれていますが、タンパク質名を含むビューがあります。
   6. [ ファイル]>[CSVファイルを保存]をクリックして変更を保存します。
   7. ナビゲーションメニュー を使用して、 scriptsフォルダに移動するには、 scriptsをクリックします。ノートブック 06_Downstream_Pathways.ipynb をダブルクリックして開きます。上部にある [カーネルを再起動してすべてのセルを実行する ]ボタンをクリックし、ポップアップで [再起動 ]をクリックして、スクリプトを実行します。
   8. 左側のナビゲーションメニューを使用して 06_figures フォルダに移動するには、最初に 図 をダブルクリックしてから 06_figuresをダブルクリックします。生成されたプロット FIG06_Pathways_and_Genes_ 'mofa_result_name をダブルクリックして開き、視覚化された経路を調査します (図 8C)。
      注:視覚化された経路の選択方法は、設定ファイルを使用して構成できます。詳細については、パラメーターのドキュメントを参照してください。
   9. 左側のナビゲーション メニューを使用して、生成された 06_results フォルダーに移動するには、 フォルダー 記号をクリックし、 結果 と 06_resultsをダブルクリックします。エンリッチされたパス ウェイ 06_Pathway_enrichment__'mofa_result_name' を含むファイルが含まれているかどうかを確認します。

6. 異なる構成とバージョンの比較(補足図1、補足図2、補足図3、補足図4)

1. ワークフロー全体で異なるパラメーター/構成を使用した場合の影響を比較するには、セクション 3 から 5 を再実行し、構成ファイル内のパラメーターを変更し、異なる configuration_name 識別子と mofa_result_name 識別子を使用します。
   注:新しい結果は、異なる実行の比較に使用されるために、これらの名前で保存されます。
2. Jupyter-Lab で、左側のナビゲーション メニューを使用して configurations フォルダーに移動します。ファイル 07_Comparison_Configs.csv をダブルクリックして開きます。
3. [mofa_result_name] 列に、比較する以前のすべての MOFA 実行の名前を追加します (名前/構成ごとに 1 行 (例: MI_v1_MOFA、MI_v2_MOFA)。
4. 「compare_factors」列に、モデル間で比較する因子を追加します。デフォルトでは、Factor1,Factor2,Factor3 です。(補足図2A)。
   注:このスクリプトでは、さまざまなモデルの特徴値と因子値を関連付けることによって比較されます。これは、同じサンプル( sample_idで示される)と同じフィーチャセットに基づくモデルでのみ機能します。サンプルまたは特徴が比較バージョン間で一致しない場合、それらは比較から除外されます。
5. [ ファイル]>[CSVファイルを保存]をクリックして変更を保存します。
6. 左側のメニューを使用して、スクリプトをクリックしてスクリプトフォルダに移動します。ノートブック 07_Compare_Models.ipynb をダブルクリックして開きます。スクリプトを実行するには、[カーネルを再起動してすべてのセルを実行]ボタンをクリックし、ポップアップで[再起動]をクリックします。
7. 左側のメニューを使用して、最初に図をダブルクリックしてから06_figuresをクリックして06_figuresフォルダーに移動します。生成されたプロットを開くには、ファイルをダブルクリックして、異なるバージョンの類似性を分析します。
   FIG07_Variance_Model_Comparison.pdf(補足図2B)
   FIG07_Factor_Correlations.pdf
   (補足図2C)
   FIG07_Feature_Correlations.pdf
   (補足図3C)

7. ワークフローの拡張:他のパラメーターと設定の追加

注: 設定ファイルで現在設定可能なパラメータの他に、コードや他のパラメータに他の調整が含まれている場合があります。例えば、MOFAモデル自体は、コード内で直接変更することも、設定ファイルを通じて調整することもできる他のいくつかの訓練パラメータ¹⁷ を提供する。プロトコルの次のセクションでは、追加のMOFAモデルのトレーニングパラメータに対してこれを行う方法の例を概説します。この部分では、Rプログラミングの知識が必要です。

1. Jupyter-Lab で、左側のナビゲーション メニューを使用して scripts フォルダーに移動します。ノートブック 03_Run_MOFA.ipynb をダブルクリックして開きます。
2. 左側の [目次 ]タブをクリックし、[ 4.3 MOFAトレーニングオプションの設定 ]サブセクションに移動し、クリックしてモデルトレーニングを実行します。下にスクロールして、ノートブック内の構成可能なパラメーターのMOFAモデルの印刷出力を確認します。
3. 見出しの下のコードの R for ループ内では、MOFA データ、モデル、およびトレーニング オプションがすべて設定されます。 model_opts$num_factors = mofa_configs$amount_of_factors[i] の行の下に、次のコードを含む別の行を追加します
   model_opts$likelihoods['data_type'] = 'ポアソン'.
   注: これにより、すべての MOFA 実行の data_type という名前で指定されたビューの入力としてモデルが受け取る分布が変更されます。データ型にポアソンを指定すると、このデータ型の特徴が整数 (RNA-seq からの読み取りカウントなど) の場合にのみモデルが実行されます。MOFAデータ、トレーニング、およびモデルオプションの詳細については、MOFAのチュートリアルとドキュメント¹⁷も参照してください。
4. ノートブックの変更を保存するには、上部にある [保存 ] ボタンをクリックします。
5. .csv設定ファイルを介して新しいパラメータを引き渡すには、左側のナビゲーションを使用して、configurationsをダブルクリックしてconfigurationsフォルダに移動し、ダブルクリックしてファイル03_MOFA_Configs.csvを開きます。
   1. パラメータ名を指定する新しい列( 例:number_iterations )を追加し、 値(例:1000)を入力します。[ ファイル]>[CSVファイルの保存]をクリックして変更を保存します。
   2. ナビゲーションメニューを使用して、 scriptsをクリックしてscriptsフォルダをナビゲートします。ノートブック 03_Run_Mofa.ipynb をダブルクリックして開きます。左側の [目次 ]タブをクリックし、[ 4.3 MOFAトレーニングオプションの設定 ]サブセクションに移動し、クリックしてモデルトレーニングを実行します。
   3. 行 train_opts$maxiter = 50000 を train_opts$maxiter = mofa_configs$column_name[i] に置き換えます (追加された列の名前が number_of_iterations の場合は train_opts$maxiter = mofa_configs$number_of_iterations[i]]) です。
      注: このノートブックで03_MOFA_Configs.csv設定ファイルは、ノートブックの先頭 (サブセクション: 前提条件 設定とパラメータ) でセッションの data.frame として読み取られmofa_configしたがって、このコード行では、このオブジェクトと対応する新しく生成された列が参照されます。複数の構成を同時に実行できるため、i は data.frame の行を識別し、モデル推定は .csv ファイル内のすべての異なる行に対して for ループで実行されます。ノートブックの冒頭にある「前提条件、構成とパラメータ」のセクションにある設定ファイルを読み込む原則は、すべてのノートブックで同じであり、このようにさらに変更を加えることができます。
   4. [ 保存 ] ボタンをクリックして、ノートブックの変更を保存します。

結果

ワークフローが正常に実行されると、 図 2 に示すように、いくつかのテーブルと図が生成されます。図は /figures フォルダ (図 6、 図 7、 図 8、 補足図 1、 補足図 2、 補足図 3、 補足図 4) に配置され、テーブルは指定された /results フォルダに配置されます。

ワークフローの実行が成功しない場合、これは主に、メモリ不足(特に大規模な単一セルデータセットがロードされる最初のステップ)によって引き起こされる技術的なエラー、データの形式が正しくない(データセット間のsample_id列が一致しないなど)、または設定ファイルの仕様が正しくない(多くの機能に除外されるなど)が原因である可能性があります。この場合、通常、実行中にJupyter-notebookスクリプト内でエラーメッセージが表示され、プロットとデータは生成されません。スクリプトの実行中に生成されたデフォルトの設定ファイルを使用し、プロトコルで説明されている特定のパラメータのみを変更することをお勧めします。

成功した実行は、結果のプロットとテーブルの生成によって示され、各ステップでは、データとそれに固有の分散の主なパターンに関する追加情報が明らかになります。それでも、必ずしも各実行が生物学的に有用で解釈可能な結果をもたらすわけではありません。多くの場合、データは大きな技術的影響と異なる分布によって特徴付けられますが、これらは「データの前処理と調和」ステップまたは「MOFA⁹ モデル」(入力データタイプに異なる分布を指定することもできます)で考慮する必要があります。これにより、基礎となる生物学的プロセスを反映したデータの変動を抽出できます。

提示されたワークフロー内では、異なるマルチオミクスデータセットを入力として使用できます。現在、ワークフローでは、シングルセルデータ用の一般的な .h5ad ファイル形式と、他のすべてのデータセット用の非常に一般的な .csv ファイル形式を入力として受け入れています(図 3)。オミクスデータセットが異なれば、ファイル形式も大きく異なるのが一般的です。ワークフローの実行を特定のファイル形式に制限しないように、.csv は非常に一般的な形式として使用されます。したがって、あらゆる種類の異なるオミクスデータセットをワークフローの入力として使用できますが、このワークフローで使用する前に、まず図3に示すように対応する.csv形式に変換する必要があります。これは、スプレッドシートまたはオミクス専用のソフトウェアを使用して準備できます。異なるオミクスデータセットを前処理するために、ワークフロー内でいくつかのオプションが利用可能で、02_Pre_Processing_Configs.csvファイルと02_Pre_Processing_Configs_SC.csvファイルを設定することで、異なる入力データセットに異なる前処理および正規化ステップ(ライブラリサイズ調整、対数変換、サンプル分位点の正規化など)を適用することができます(図2).それにもかかわらず、ここで利用可能なオプションは、主にここに示すデータセットで利用可能な特定の入力データ(scRNA-seq、サイトカインアッセイ、プロテオミクス、prime-seq)に基づいています。他のオミクス/データタイプを使用する場合は、既存のベストプラクティスに従って、追加のオミクス固有の正規化手順を適用する必要がある場合があります。この場合、データはすでに前処理された形式でワークフローに引き渡すことができ、追加の前処理手順を適用せずに他のデータセットと一緒に統合されます。多くの場合、[Feature Wise Quantile Normalization] ステップを適用すると、すべてのデータ タイプの分布を正規分布に揃え、異なる入力フィーチャ間のダウンストリーム解析をガウス ノイズのモデル仕様との比較と互換性を高めるのに役立ちます。

ワークフローの実行中に、データ統合とその後の生物学的下流の解釈のプロセスをサポートするいくつかのプロットと出力が生成されます。scRNA-seqデータの場合、 FIG01_Amount_of_Cells_Overview のプロット(図6)は、サンプルあたりの細胞数が少なすぎる可能性のある細胞タイプと、遺伝子発現シグナルを確実に測定するための細胞タイプを示しています。その後の解析では、サンプルあたりの細胞タイプの全細胞の平均値を発現推定値として使用します(psedobulk-アプローチ)。このユースケースでは、大部分のサンプルで細胞が3つ未満の細胞タイプを除外します。

分散分解プロットFIG03_Overview_Variance_Decomposition(図7、補足図1)は、さまざまなデータソースがどの程度うまく統合されているか、およびさまざまなデータソースの分散が各データソースに共有され、一意であることを示すことができます。たとえば、ここで使用したデータセットでさまざまな前処理戦略をテストすると、たとえば、前処理から特徴量分位数の正規化ステップを削除すると、特定のデータ ビューに重点を置く潜在的な要因が発生し、プロテオミクス データと他のデータ ソースとの統合が減少することがわかります。これは、説明された分散の量が減少していることに見ることができます(補足図1B)。特徴のフィルタリングや正規化を行わずにMOFAモデルを実行すると、潜在因子によってキャプチャされた異なるビュー間の共有分散が少なくなります(補足図1C)。これは、潜在的な要因が主にデータ型固有の技術的影響を反映していることを示しています。それに加えて、MOFA⁹モデル自体も、データの前処理が不十分な場合に警告を返す可能性があります。このような警告の例を、代替の前処理構成 (MI_v2 と MI_v3) の補足図 1 に示します (特定の構成ファイルの例は、クローンされた GitHub リポジトリの config_examples フォルダーに格納されています)。

さらに、MOFAモデルを実行した後、因子をサンプルに関する既知の生物学的メタ情報、および技術的およびその他の交絡共変量(04_Downstream_Factor_Analysis)と関連付けて、因子によって捕捉された変動の考えられる原因を特定することにより、いくつかのダウンストリーム分析で結果を評価できます。たとえば、MOFA モデルの因子の 1 つが技術共変量の 1 つ (バッチ情報など) と強く関連している場合、この因子は生物学的な変動ではなく、データ内の技術的な変動を捉えていることを示している可能性があります。

ダウンストリーム解析部分の生物学的解釈を絞り込むために、入力データセットに基づくいくつかの知見(より洗練された解釈は元の出版物¹¹に記載されています)をここで概説します。最初のステップでは、適用された前処理戦略により、複数の細胞タイプだけでなく、他のオミクスデータタイプ間での分散を捉えるいくつかの因子が見つかることを観察できます(図7A)。例えば、Factor 2は、scRNA-seqデータセットの臨床インプット特徴といくつかの細胞タイプの分散を捉えます。最初の3つの因子を「CRP」や「CK」などの関連する臨床共変量(図7B)と関連付け、異なる患者サブグループ(「コントロール(CCSおよび非CCSを含む)」と異なる時点(TP1-TP4)で測定された「ACS」(図7C)の因子値の違いを調査すると、Factor2は「CK」値と有意に関連し、Factor3は「CRP」値と有意に関連していることもわかりました。同時に、TP1 および TP2 の「ACS」サンプル (心筋梗塞 (MI) に対する免疫応答の急性期を反映) は、「コントロール」およびそれ以降の時点のサンプル (TP3/TP4) と比較して因子値の増加を示しています。CKは心筋損傷の既知のマーカーであり、通常、Factor2によって捕捉されるパターンと同様に、TP1/TP2での値の増加によって特徴付けられます。

Factor2を形作る生物学的プロセスに関する洞察を得るために、モデルによって生成された特徴重みテーブル(03_Weight_Data.csv)を見て、因子の上位の特徴を評価します。因子の絶対的な重みが最も高い特徴の上位1%を分析すると、主にCD4が見つかります。TCMおよびCD14。モノ由来の特徴は、入力特徴の総数と比較して過大評価されており(図8A)、これらの細胞タイプがMI後の炎症過程に非常に関連性があることを示しています(注:前処理で特徴ごとの分位点の正規化が適用されなかった場合、特徴の異なる分布もこの結果に影響を与える可能性があり、評価はデータタイプごとに個別に行う必要があります)。CD4の上位機能を分析します。TCM細胞タイプを因子として、堅牢なT細胞の活性化に必要なEIF3E¹⁸ や、T細胞の増殖と活性化を促進するHMGB1¹⁹など、いくつかの興味深い遺伝子が見つかります(図8B)。次に、REACTOME²⁰ データベースの免疫経路をパスウェイセット(Prepared_Pathway_Data.csv)として、パスウェイエンリッチメント解析を行います。「インターロイキン-6」シグナル伝達を含むいくつかの「インターロイキン」経路の濃縮が見出されます。scRNA-seqデータの異なる細胞タイプにおけるいくつかの遺伝子の発現レベルと、サイトカインアッセイによって測定された「IL6」サイトカイン値が、この結果に寄与しました(図8C)。これらのデータタイプ間で共有されるパターンを特定することで、統合分析の付加価値が強調されます。全体として、このアプローチは、対応する出版物¹¹でより詳細に説明されているように、疾患状態または関連する治療結果と基礎となる多細胞免疫プログラムを反映する他のいくつかの要因も特定することができる。

複数のオミクスにわたる統合解析の利点をさらに強調するために、プロテオミクスの入力データのみを含めて同じワークフローも実行しました(補足図4)。結果として得られる因子を分析すると、統合分析と同様に、「CRP」値と強く相関する因子(Factor1)が見つかります。このパターンは、プロテオミクスデータ内の変動の主な原因を説明し、統合解析で「Factor3」によってキャプチャされた他のデータセットの変動の一部とも一致しています(図7C)。しかし、Factor2が示したような、統合解析で炎症の時間経過を捉えるパターンは、プロテオミクスのデータだけでは特定できません。

導入されたワークフローとMOFA⁹ モデル自体は、多くの調整可能なパラメーターで高度にカスタマイズ可能です。したがって、さまざまな構成によって生成された結果を視覚化し、体系的に比較することが重要です。このタスクを容易にするために、ワークフローで生成できる最終的な出力は、前処理とモデル推定の異なるパラメーターを持つパイプラインのさまざまな名前付き実行の比較です。例えば、外務省モデルは、異なる数の潜在因子で推定したり(補足図2A)、特徴の数が少ないビューに重み付けをしたり(補足図3A)することができます。ワークフロー '07_Compare_Models' の最後のスクリプトを構成して実行すると、異なるパイプライン実行間の類似性を評価するためのいくつかのプロットが生成されます。FIG07_Variance_Model_Comparison(補足図2B、補足図3B)は、異なる実行について各ビューについて説明された合計分散の比較を示しています。異なる実行間の因子値と特徴因子の重みの相関関係は、特定のパラメータを変更したときに結果がどの程度変化するかを示すことができます(補足図2C、補足図3C)。ここで、因子の数を変更しても、推定因子値と特徴量の重みがわずかに変更されるだけです(補足図2C)。データビューの重み付けを変更すると、特徴の数が少ないビュー(例:「臨床」ビュー(補足図3B))で説明分散がはるかに高くなります。それにもかかわらず、最初の3つの因子内の関連する特徴は、重み付けされていないバージョンで推論された特徴と依然として高い相関があります(補足図3C)。

生成されたモデル出力.csvファイル(推定係数や特徴重量など)を使用して、さらに個々の下流分析を行うことができます。すべてのコードと必要な設定ファイル (ドキュメントを含む) は、GitHub の https://github.com/heiniglab/mofa_workflow で入手できます。分析に必要なcondaパッケージを簡単にインストールできるようにするために作成されたシンギュラリティ・イメージは、https://doi.org/10.5281/zenodo.10815146 からダウンロードできます。パイプラインの初期テストを実行するために使用できる小さなサンプル データセットは、同じ zenodo レコードからもダウンロードできます。

図7:MOFAの出力分析。 MOFAモデル(03_Run_MOFA.ipynb)の実行と因子値の下流分析(04_Downstream_Factor_Analysis.ipynb)の後、いくつかのプロットが生成されます:(A)FIG03_Overview_Variance_Decomposition:異なるビュー内で推定されたMOFA因子の説明された分散の視覚化を返します。ヒートマップ (左): 各ビューの係数によってキャプチャされたビューの合計分散の割合を示します。棒グラフ(右):各ビューのすべての因子によって捕捉された分散の合計パーセンテージを示します。(B)FIG04_Factor_Association_Numerical_Features:因子値と選択された数値サンプル共変量(ここでは臨床変数(CRP、CK)とのピアソン相関を示します。(C) FIG04_Factor_Association_Categorical_Features:カテゴリ別サンプル共変量の因子値の差を箱ひげ図として表示します。ここでは、ACS患者と対照患者の各時点の因子1-3の因子値を比較します。 この図の拡大版を表示するには、ここをクリックしてください。

図8:外務省の特徴分析。 ダウンストリーム解析(04_Downstream_Factor_Analysis.ipynb、05_Downstream_Investigate_Features.ipynb)の実行後、いくつかのプロットが生成されます。ここにあるすべてのプロットは、MOFAファクター2を視覚化しています:(A)FIG04_Top_Feature_Overview_per_Factor:ヒートマップ(左)は、選択したファクターによってキャプチャされた分散の割合を各ビューに示します。棒グラフ(右)は、因子のさまざまなビューの特徴の関連性を示します。左側には、特定のビューの特徴量の合計が、その要素に関するビュー全体で上位 1% の最高ランクの特徴量に表示されます。右側には、上位 1% の合計数をそのビューの機能の合計数で割ったパーセンテージが表示されます。(B)FIG05_Heatmap_Feature_Overview:ヒートマップ(左)は、CD4の最高ランクの1%の特徴を示しています。TCM細胞型:各サンプルの正規化された発現値で、「対照」グループの患者(CCSおよび非CCS)を「ACS」患者の異なる時点と比較します。棒グラフ (右) は、特徴の重みを示しています。重みの符号の方向は、セルタイプ名の前に左側に示されています:「+」正の因子の重み。'-' 負の係数の重み。(C)FIG06_Pathway_and_Genes:濃縮インターロイキン経路に属する因子の上位25%ランクの遺伝子の重みを示します。上部のヒートマップではビュー間で平均化され、下部のヒートマップではビューごとに平均化されます。 この図の拡大版を表示するには、ここをクリックしてください。

補足図1:データ調和効果。 この図は、いくつかの異なるデータ前処理構成のFIG03_Overview_Variance_Decompositionを示しています。つまり、さまざまなビュー内で推定されたMOFA因子の説明された分散の視覚化です。ヒートマップ (左): 各ビューについて、ファクターによってキャプチャされたビューの合計分散の割合を示します。棒グラフ(右):各ビューについて、すべての因子によって捕捉された分散の合計パーセンテージを示します。(A)以前の図で生物学的ダウンストリームの結果が分析された基礎となる構成(「MI_v1」)(クローンリポジトリのデフォルトの構成ファイルで設定されたパラメータ)。(B) 「MI_v1」と同じ前処理構成に、 機能ごとの分位数正規化を適用しない という変更を加えたもの(リポジトリの「config_examples」フォルダ内の例示的な設定ファイルのように設定されたパラメータ)。この構成のMOFAモデル出力警告のスクリーンショットを以下のプロットに追加しています。(C) 前処理ステップが適用されず、 すべてのデータが前処理または特徴のフィルタリングなしで入力として使用された場合に生じる分散分解(リポジトリの「config_examples」フォルダ内の例示的な構成ファイルで設定されたパラメータ)。この構成のMOFAモデル出力警告のスクリーンショットを以下のプロットに追加しています。 このファイルをダウンロードするには、ここをクリックしてください。

補足図2:外務省の構成-ファクター量の影響。 MOFA モデルを実行するためにいくつかの異なる構成を使用して '07_Compare_Models.ipynb' スクリプトによって生成された結果の数値。(A) '03_MOFA_configs.csv': いくつかの異なる量の要素 (10,15,20,25) を指定する '03_Run_MOFA.ipynb' スクリプトを実行するために使用されるさまざまな構成の例。'07_Comparison_configs.csv': スクリプト '07_Compare_Models.ipynb' を実行するための設定入力ファイルを指定する方法の例。(B)「FIG07_Variance_Model_Comparison」は、モデルで指定されたすべての要因にわたるさまざまなモデルの各ビュー(y軸)について説明された分散の合計を示します。(C)「FIG07_Factor_Correlations」は、異なる構成間の因子サンプル値の相関を示しています。 このファイルをダウンロードするには、ここをクリックしてください。

補足図3:MOFAの構成 - 重み付けビューの影響。 MOFA モデルを実行するためにいくつかの異なる構成を使用して '07_Compare_Models.ipynb' スクリプトによって生成された結果の数値。(A) '03_MOFA_configs.csv': '03_Run_MOFA.ipynb' スクリプトの実行に使用されるさまざまな設定の例で、'weighting_of_views' パラメータを 'TRUE' (MI_v1_MOFA_weighted) または 'FALSE' (MI_v1_MOFA) に指定します。'07_Comparison_configs.csv': スクリプト '07_Compare_Models.ipynb' を実行するための設定入力ファイルを指定する方法の例。(B)「FIG07_Variance_Model_Comparison」は、モデルで指定されたすべての要因にわたるさまざまなモデルの各ビュー(y軸)について説明された分散の合計を示します。(C)「FIG07_Feature_Correlations」は、異なる構成間の特徴因子の重みの相関を示す。 このファイルをダウンロードするには、ここをクリックしてください。

補足図4:マルチオミクス統合効果 - プロテオミクスデータのみを使用。 結果として得られるパターンは、プロテオミクスデータを入力としてのみ使用した場合に潜在因子によって捕捉されます。(A)FIG04_Factor_Association_Numerical_Features:因子値と臨床変数(CRP、CK)のピアソン相関。(B) FIG04_Factor_Association_Categorical_Features: ACS 患者と対照患者の各時点の因子値の箱ひげ図比較。 このファイルをダウンロードするには、ここをクリックしてください。

補足ファイル1:Supplementary_File_
Running_Pipeline_with_Exemplary_Data。 サンプル データでパイプラインを実行する方法と予想される出力については、追加で提供される補足ファイルに記載されています。 このファイルをダウンロードするには、ここをクリックしてください。

補足ビデオファイル1:プロトコルのスクリーンキャプチャビデオ。このファイルをダウンロードするには、ここをクリックしてください。

ディスカッション

概説されたプロトコルでは、複雑なマルチオミクスデータセットを迅速に探索するために使用できる、モジュール式で拡張可能なJupyter-notebookベースのワークフローが提示されます。ワークフローの主要な部分は、前処理とデータ調和の部分(データのフィルタリングと正規化のためのさまざまな標準ステップを提供)、MOFA⁹ モデルの推定、およびいくつかの模範的なダウンストリーム分析で構成されています。主な重要なステップの1つは、さまざまなオミクスデータセットを前処理して統合し、調和させることです。ここでは、scRNA-seqデータ、prime-seqバルクRNA、サイトカインアッセイ、血漿プロテオミクス、および臨床値を含むデータセットの戦略を提示し、その結果、MI¹¹の関連する生物学的プロセスを特定するために使用できる統合的で調和のとれたデータセットが得られました。他のオミクスデータセットで追加のデータ前処理戦略が必要な場合は、外務省の分析前に実行する必要があります。その後、前処理されたデータを現在のワークフローへの入力として使用できます。モデルの出力を使用して、さまざまなデータセットの統合の品質や前処理の影響を評価し、潜在的な技術的影響を特定できます。たとえば、最小限の前処理と正規化の手順のみを適用することから始めて、その後、さらに正規化の効果を評価できます。このアプリケーションでは、「機能ごとの分位点の正規化」を追加することで、プロテオミクスと残りのアッセイとのより良い統合につながることが観察されました。

ワークフローのもう1つの大きな部分は、サンプルレベルで集計されたデータに対するMOFA⁹ モデルの推定です。MOFAモデルの他の拡張は、次のようなものがあります
MOFA+²¹(特にシングルセルデータ用)、MEFISTO²² (特に時間成分を含むデータ用)、およびMuVI²³ (因子分析アプローチにおけるドメイン知識の統合)。これらは、特定のタイプのデータセットまたは設定に対する方法の非常に便利な拡張ですが、それらを適用するには特定の要件が伴います。例えば、今回のケースでMOFA+²¹ をシングルセルデータに限定すると、他のオミクスデータのようなサンプルレベルのデータを簡単に利用できなくなることになります。MEFISTO²² 法を使用すると、時間経過を具体的にモデル化することができますが、複数の時間点が利用できない設定や、この場合、1つの時点での測定値のみを持つ「制御」サンプルと「ACS」サンプルの時間分解データを組み合わせる場合には適用されません。したがって、この非常に一般的なワークフローでは、MOFA⁹の方法を使用することを選択します。これは、多くの要件をあまり必要とせずにあらゆる種類のデータセットを迅速に探索するための最も柔軟な方法であるためです。注目すべきは、ワークフローで生成された前処理されたデータをこれらの方法で分析することも、データセットがメソッドの要件と互換性がある場合は、現在のスクリプトを拡張してこれらのメソッドを統合することもできます。その他のユースケース、チュートリアル、およびドキュメントについては、MOFA developers²⁴ の GitHub リポジトリを参照してください。

主成分分析(PCA)のようなさらに一般的な次元削減方法と比較して、MOFA⁹モデルは、特にマルチオミクス設定において、いくつかの利点を提供します。たとえば、分散分解は入力ビューごとに簡単に分析でき、重みをさまざまなビューに割り当てて、さまざまな特徴数を考慮に入れることができます。このモデルではスパース性が促進され、結果の解釈可能性が向上します。さらに、オミクスの1つに欠損データがあるサンプルを除外する必要はなく、モデルはスパース特徴因子の重みの学習に焦点を当てるためにいくつかの特徴を実装します。また、MOFAモデルには、異なる分布によって特徴付けられたデータを統合するためのいくつかの設定が用意されています(「ガウス」、「ベルヌーイ」、または「ポアソン」の尤度のいずれかをモデル化)。このワークフローでは、連続データのみを統合し、それを「ガウス」分布に従うように正規化しますが、必要に応じて既存のコードを拡張して他のデータ型を統合することもできます。scITD²⁵など、特にscRNA-seqデータには、他の分解ベースの方法も存在しますが、他のオミクスと連携するように設計されていないという制限があります。

提示されたワークフローは、大規模で複雑なマルチオミクスデータセットの教師なし探索のプロトコルを示しており、データ内の変動を引き起こす潜在的な生物学的プロセスやその他の特徴を特定します。これは、病気やその他の生物学的または技術的な摂動によって引き起こされるデータの変動を調査する必要があるほぼすべての設定に適用できます。結果として得られるダウンストリーム解析と、異なるオミクス間の分散を駆動する特徴セットにより、特定の状況でさらに調査するための関連する生物学的プロセスを明らかにすることができます。例えば、疾患の文脈では、新しい診断マーカーや治療標的の同定が可能になるかもしれません。

材料

この記事で使用された材料の一覧

名前	会社	カタログ番号	コメント
Apptainer	NA		NA https://apptainer.org/docs/admin/main/installation.html
Compute Server or Workstation or cloud (Linux、Mac、またはWindows環境). さまざまな入力データセットのサイズに応じて、適切なマシンでワークフローを実行することをお勧めします（私たちの設定では、16 CPU、64GBメモリを使用します。	どのメーカーの	16 CPU、64GBメモリ	大容量メモリは、生のシングルセルデータの処理にのみ必要です。前処理後、後の解析手順は、通常のデスクトップまたはラップトップ コンピューター
git	NA		https://git-scm.com/book/en/v2/Getting-Started-Installing-Git
GitHub	GitHub	NA	でも実行できます https://github.com/heiniglab/mofa_workflow