ここでは、トランスクリプトームデータをmqTransビューに変換し、ダークバイオマーカーの同定を可能にするプロトコルを紹介します。これらのバイオマーカーは、従来のトランスクリプトーム解析では発現差が認められませんが、mqTransでは発現差が見られます。このアプローチは、従来の手法を補完する手法として機能し、これまで見過ごされていたバイオマーカーを明らかにします。
Method Article
ここでは、トランスクリプトームデータをmqTransビューに変換し、ダークバイオマーカーの同定を可能にするプロトコルを紹介します。これらのバイオマーカーは、従来のトランスクリプトーム解析では発現差が認められませんが、mqTransでは発現差が見られます。このアプローチは、従来の手法を補完する手法として機能し、これまで見過ごされていたバイオマーカーを明らかにします。
トランスクリプトームは、サンプル中の多くの遺伝子の発現レベルを表し、生物学研究や臨床現場で広く使用されています。研究者は通常、サンプルの表現型グループとコントロールグループの間で異なる表現を持つトランスクリプトームバイオマーカーに焦点を当てました。本研究では、参照サンプルの複雑な遺伝子間相互作用を学習するためのマルチタスクグラフアテンションネットワーク(GAT)学習フレームワークを提示しました。実証的参照モデルは、健康なサンプル(HealthModel)で事前にトレーニングされており、独立したテストトランスクリプトームのモデルベースの定量的転写調節(mqTrans)ビューを生成するために直接使用できます。生成されたトランスクリプトームの mqTrans ビューは、予測タスクとダークバイオマーカー検出によって実証されました。「ダークバイオマーカー」という造語は、ダークバイオマーカーはmqTransビューでは差異表現を示すが、元の発現レベルでは差異的発現を示さないという定義に由来しています。ダークバイオマーカーは、発現差がないため、従来のバイオマーカー検出研究では常に見過ごされていました。パイプラインHealthModelPipeのソースコードとマニュアルは、http://www.healthinformaticslab.org/supp/resources.php からダウンロードできます。
トランスクリプトームは、サンプル中のすべての遺伝子の発現で構成されており、マイクロアレイやRNA-seq1などのハイスループット技術によってプロファイリングすることができます。データセット内の1つの遺伝子の発現レベルはトランスクリプトーム特徴と呼ばれ、表現型と対照群の間でトランスクリプトームの特徴が異なって表現されることで、この遺伝子はこの表現型のバイオマーカーとして定義されます2,3。トランスクリプトームバイオマーカーは、疾患診断4、生物学的メカニズム5、生存解析6,7などの研究に広く利用されています。
健康な組織の遺伝子活性パターンは、生命に関する重要な情報を運びます8,9。これらのパターンは、非常に貴重な洞察を提供し、良性疾患10,11および致死性疾患12の複雑な発達の軌跡を理解するための理想的な参考資料として機能します。遺伝子は互いに相互作用し、トランスクリプトームは複雑な相互作用の後の最終的な発現レベルを表します。このようなパターンは、転写調節ネットワーク13および代謝ネットワーク14等として定式化される。メッセンジャーRNA(mRNA)の発現は、転写因子(TF)と長鎖遺伝子間ノンコーディングRNA(lincRNA)によって転写制御されている15,16,17。従来の差次的発現解析では、このような複雑な遺伝子相互作用は無視され、特徴間の独立性を仮定していた18,19。
グラフニューラルネットワーク(GNN)の最近の進歩は、がん研究20、例えば共発現モジュールの同定21など、OMICベースのデータから重要な情報を抽出する上で並外れた可能性を示しています。GNNの生来の能力は、遺伝子間の複雑な関係と依存関係をモデル化するのに理想的です22,23。
生物医学研究では、対照群に対する表現型を正確に予測することに重点が置かれることがよくあります。このようなタスクは、一般に二項分類24、25、26として定式化されます。ここで、2 つのクラス ラベルは通常、1 と 0、true と false、または正と負の27 としてエンコードされます。
この研究は、事前学習済みのグラフアテンションネットワーク(GAT)参照モデルに基づいて、トランスクリプトームデータセットの転写制御(mqTrans)ビューを生成するための使いやすいプロトコルを提供することを目的としています。以前に発表された研究26 のマルチタスクGATフレームワークを使用して、トランスクリプトームの特徴をmqTransの特徴に変換しました。カリフォルニア大学サンタクルーズ校(UCSC)のXenaプラットフォーム28 の健康なトランスクリプトームの大規模なデータセットを使用して、制御因子(TFおよびlincRNA)から標的mRNAまでの転写制御を定量的に測定する参照モデル(HealthModel)を事前トレーニングしました。生成されたmqTransビューは、予測モデルの構築やダークバイオマーカーの検出に使用できます。このプロトコルは、例示的な例として、がんゲノムアトラス(TCGA)データベース29 からの結腸腺癌(COAD)患者データセットを利用します。これに関連して、ステージIまたはIIの患者は陰性サンプルに分類され、ステージIIIまたはIVの患者は陽性サンプルと見なされます。また、26種類のTCGAがんにおけるダークバイオマーカーと従来のバイオマーカーの分布も比較しています。
HealthModel パイプラインの説明
このプロトコルで採用されている方法論は、図1に概説されているように、以前に公開されたフレームワーク26に基づいています。まず、ユーザーは入力データセットを準備し、提案されたHealthModelパイプラインにフィードし、mqTrans特徴を取得する必要があります。詳細なデータ準備手順は、プロトコルセクションのセクション2に記載されています。その後、ユーザーは、mqTransの特徴を元のトランスクリプトームの特徴と組み合わせるか、生成されたmqTransの特徴のみを続行するかを選択できます。生成されたデータセットは、特徴選択プロセスにかけられ、ユーザーは分類のためのk分割交差検証でkの優先値を柔軟に選択できます。このプロトコルで使用される主な評価指標は精度です。
HealthModel26は、トランスクリプトームの特徴をTF(転写因子)、lincRNA(長鎖遺伝子間ノンコーディングRNA)、mRNA(メッセンジャーRNA)の3つのグループに分類しています。TFの特徴は、Human ProteinAtlas30,31で利用可能なアノテーションに基づいて定義されています。この研究は、GTExデータセット32からのlincRNAのアノテーションを利用しています。KEGGデータベース33における第3レベルの経路に属する遺伝子は、mRNAの特徴と見なされる。TRRUSTデータベース34に記録されているように、mRNAの特徴が標的遺伝子に対する調節的役割を示す場合、それはTFクラスに再分類されることは注目に値する。
また、このプロトコルでは、規定因子(regulatory_geneIDs.csv)と標的mRNA(target_geneIDs.csv)の遺伝子IDの2つのサンプルファイルを手動で生成します。調節特性(TFおよびlincRNA)間のペアワイズ距離マトリックスは、ピアソン相関係数によって計算され、一般的なツール加重遺伝子共発現ネットワーク解析(WGCNA)36 (adjacent_matrix.csv)によってクラスタリングされます。ユーザーは、HealthModel パイプラインをこれらのサンプル構成ファイルと共に直接利用して、トランスクリプトーム データセットの mqTrans ビューを生成できます。
HealthModelの技術的詳細
HealthModel は、TF と lincRNA の間の複雑な関係をグラフとして表し、入力特徴は V で示される頂点として機能し、頂点間エッジ行列は E として示されます。各サンプルは、VK×1として表されるK個の調節特徴によって特徴付けられます。具体的には、データセットには 425 の TF と 375 の lincRNA が含まれ、サンプルの次元は K = 425 + 375 = 800 でした。エッジ行列 E を確立するために、この作業では一般的なツール WGCNA35 を使用しました。と
で表される
2 つの頂点をつなぐペアワイズ重みは、ピアソン相関係数によって決定されます。遺伝子調節ネットワークは、極めて重要な機能的役割を有するハブ遺伝子の存在を特徴とするスケールフリートポロジー36を呈する。2 つの特徴または頂点 と
の間の相関関係は、
トポロジカル オーバーラップ メジャー (TOM) を使用して次のように計算します。
(1)
(2)
ソフトしきい値βは、WGCNA パッケージの 'pickSoft Threshold' 関数を使用して計算されます。べき乗指数関数 aij が適用され、ここで
、 は i と j を除いた遺伝子を表し、
頂点接続性を表します。WGCNAは、トランスクリプトームの特徴の発現プロファイルを、一般的に採用されている非類似度尺度
を用いて複数のモジュールにクラスタリングする(37.
HealthModelフレームワークは、もともとマルチタスク学習アーキテクチャ26として設計されました。このプロトコルは、トランスクリプトーム mqTrans ビューの構築にモデルの事前トレーニング タスクのみを利用します。ユーザは、追加のタスク特異的トランスクリプトームサンプルを用いて、マルチタスクグラフアテンションネットワークの下で事前訓練されたHealthModelをさらに洗練させることを選択することができる。
機能の選択と分類に関する技術的な詳細
特徴選択プールには、11 個の特徴選択 (FS) アルゴリズムが実装されています。その中には、最大情報量係数(SK_mic)を使用したK個の最良の特徴量の選択、MICのFPRに基づくK個の特徴量の選択(SK_fpr)、MICの偽発見率が最も高いK個の特徴量の選択(SK_fdr)の3つがフィルタベースのFSアルゴリズムです。さらに、3つのツリーベースのFSアルゴリズムは、ジニ指数(DT_gini)、適応型ブースト決定木(AdaBoost)、およびランダムフォレスト(RF_fs)を備えた決定木を使用して個々の特徴を評価します。また、このプールには、線形サポートベクター分類器による再帰的特徴量消去法 (RFE_SVC) とロジスティック回帰分類器 (RFE_LR) による再帰的特徴量消去法の 2 つのラッパー手法も組み込まれています。最後に、最上位の L1 特徴重要度値 (lSVC_L1) を持つ線形 SVC 分類器と、最上位の L1 特徴重要度値 (LR_L1) を持つロジスティック回帰分類器の 2 つの埋め込みアルゴリズムが含まれています。
分類子プールでは、分類モデルを構築するために 7 つの異なる分類子が採用されています。これらの分類器は、線形サポート ベクター マシン (SVC)、ガウス単純ベイズ (GNB)、ロジスティック回帰分類器 (LR)、k 最近傍 (k は既定で 5 に設定されている) (KNN)、XGBoost、ランダム フォレスト (RF)、および決定木 (DT) で構成されます。
データセットをtrain:testサブセットにランダムに分割する方法は、コマンドラインで設定できます。この例では、train: test = 8:2 の比率を使用します。
注: 次のプロトコルでは、主要なモジュールのインフォマティクス解析手順と Python コマンドの詳細について説明します。 図2 は、このプロトコルで利用されるコマンドの例とともに3つの主要なステップを示しており、技術的な詳細については、以前に公開された作品26、38 を参照してください。コンピュータシステムの通常のユーザーアカウントで次のプロトコルを実行し、管理者またはrootアカウントの使用を避けてください。これは計算プロトコルであり、生物医学的危険因子はありません。
1. Python環境の準備
2. 事前学習済みの HealthModel を使用して mqTrans 特徴を生成する
3. mqTrans 機能の選択
トランスクリプトームデータセットのmqTransビューの評価
このテストコードでは、11 個の特徴選択 (FS) アルゴリズムと 7 個の分類器を使用して、トランスクリプトームデータセットの生成された mqTrans ビューが分類タスクにどのように寄与するかを評価します(図 6)。テストデータセットは、The Cancer Genome Atlas(TCGA)データベース29の317の結腸腺癌(COAD)で構成されています。I期またはII期のCOAD患者は陰性サンプルと見なされ、III期またはIV期の患者は陽性サンプルと見なされます。
テストコードには11個のFSアルゴリズムが実装されています。フィルターベースの FS アルゴリズムには、MIC による K 個の最適な特徴量の選択 (SK_mic)、MIC の FPR による K 個の特徴量の選択 (SK_fpr)、MIC の最も高い FDR による K 個の特徴量の選択 (SK_fpr) の 3 つがあります。3つのツリーベースのFSアルゴリズムは、ジニ指数(DT_gini)、適応ブースト決定木(AdaBoost)、ランダムフォレスト(RF_fs)を使用して、それぞれ決定木によって個々の特徴を評価します。テスト コードの FS プールでは、線形サポート ベクター分類器 (SVC)(RFE_SVC) を使用した再帰的特徴消去 (RFE) とロジスティック回帰分類器 (RFE_LR) を使用した RFE、および最上位の L1 特徴重要度値 (lSVC_L1) を持つ線形 SVC 分類器と、最上位の L1 特徴重要度値 (LR_L1) を持つロジスティック回帰分類器の 2 つの埋め込みアルゴリズムも評価されます。
テスト コードでは、線形サポート ベクター マシン (SVC)、ガウス単純ベイズ (GNB)、ロジスティック回帰分類器 (LR)、k 最近傍法、既定で k-5 (KNN)、XGBoost、ランダム フォレスト (RF)、デシジョン ツリー (DT) の 7 つの分類器を使用して分類モデルを構築します。
図 6 は、mqTrans の特徴、元の mRNA の特徴、および各 FS アルゴリズムで推奨されている mRNA と mqTrans の特徴を組み合わせたサブセットの最大テスト精度を示しています。
組み合わされた特徴サブセット(mRNA+mqTrans)は、「SK_fpr」FS法で最高精度0.7656を達成し、個々の特徴タイプmqTrans(0.7188)および元のmRNA(0.7188)よりも優れています。他のFSアルゴリズムでも同様のパターンが見られます。ユーザーは、選択したフィーチャを出力ファイル Output-SelectedFeatures.csvで確認できます。
ダークバイオマーカーの検出
以前の研究では、表現型群と対照群の間で有意に異なるmqTrans値を持つ未分化発現遺伝子の存在が示されました26,38,39。これらの遺伝子は、従来のバイオマーカー検出研究では、その非差別的な発現によって無視されているため、ダークバイオマーカーと呼ばれています。Microsoft Excel の統計分析関数 t.test を使用して、統計的 p 値が 0.05 より小さい場合に差分で表現される特徴量を定義できます。
生成されたmqTrans値を持つ3062の特徴のうち、221のダークバイオマーカーが検出されました(図7)。3位の遺伝子ENSG00000163697(APBB2、Amyloid Beta Precursor Protein Binding Family B Member 2)は、mqTrans値が有意に異なる(mqTrans.P = 2.03 x 10-4)のに対し、元の発現レベルは発現差がない(mRNA.P = 3.80 x 10-1)。キーワードAPBB2は、PubMedデータベース40の27の出版物にヒットしたが、結腸や腸との関連は検出されなかった。
別の遺伝子ENSG00000048052(HDAC9、ヒストン脱アセチル化酵素9)は、表現型群と対照群(mRNA.P = 9.62 x 10-1)の間で実質的に同じ正規分布を維持しながら、差動的に表現されたmqTrans値(mqTrans.P = 6.09 x 10-3)を持っています。キーワードHDAC9は、PubMedデータベースの417の出版物にヒットしました。また、3件の研究では、抄録41,42,43で「結腸」または「腸」というキーワードが言及されている。しかし、結腸がんにおけるHDAC9の役割を調査した研究は行われていませんでした。
このデータは、これらの転写後の活動から、これらのダークバイオマーカー(例えば、翻訳されたタンパク質レベル44,45)をさらに評価する必要性を示唆しています。
代謝関連のダークバイオマーカーと従来のバイオマーカーの汎がん分布
代謝関連の従来のバイオマーカーがスクリーニングされ、TCGAデータセットの26種類のがんについてダークバイオマーカーと比較された38。両方のカテゴリーのバイオマーカーは、早期(ステージIおよびII)および後期(ステージIIIおよびIV)のがん病期にわたる有意レベルを識別するために統計的評価を受けた。この評価では、スチューデントのt検定のp値を使用し、その後、偽発見率(FDR)を使用して多重検定を補正しました。26種類のがんのそれぞれの詳細なデータを 図8に示す。
FDR補正p値が0.05未満の遺伝子は、従来のバイオマーカーに分類された。対照的に、ダークバイオマーカーは、mqTransビューでFDR補正p値が0.05未満であり、同時に発現レベルに統計的に有意な差を示さないバイオマーカーとして定義されました。
図9は、ほとんどのがん種において、従来のバイオマーカーと比較して、ダークバイオマーカーが全般的に不足していることを示しています。注目すべき例外としては、BRCA、MESO、およびTGCTがあり、これらはダークバイオマーカーの有病率が高いことを示しています。転写因子、メチル化パターン、遺伝子変異、環境条件など、さまざまな要因が、これらのダークバイオマーカーの転写異常を調節している可能性があることが明らかになっています。ダークバイオマーカーの発現レベルを混乱させる可能性のあるノンコーディングRNA転写産物が重複しているため、さらに複雑になる可能性があります。いくつかのダークバイオマーカーの転写調節不全は、それらのタンパク質レベルの違いによって裏付けられました44,45。ダークバイオマーカーは、従来の研究では見落とされがちであり、将来のメカニズム研究のための興味深い道筋を示しています。

図1:このプロトコルのHealthModelと機能選択モジュールの概要。 ユーザーが Python プログラミングに精通している場合は、特徴選択プールと分類器プールの特定のアルゴリズムを置き換えます。 この図の拡大版をご覧になるには、ここをクリックしてください。

図2:このプロトコルの完全なコードフロー。 (A) Python環境を準備します。まず、仮想環境を作成し、重要なパッケージをインストールします。包括的な手順については、セクション 1 を参照してください。(b) mqTrans 特徴量を生成します。提供されたコードを段階的に実行して、mqTrans 機能を取得します。詳細な説明はセクション2にあります。(C) mqTrans Featuresを選択します。このセクションでは、mqTrans 機能の評価に焦点を当てます。詳細については、セクション 3 を参照してください。 この図の拡大版をご覧になるには、ここをクリックしてください。

図3:Pythonの環境を整える (A) healthmodel を作成するコマンド。(B) VE作成処理中に y を入力します。(C) VEをアクティブにするための最も一般的なコマンド。(D)トーチ1.13.1をインストールするためのコマンド。(E) トーチジオメトリック パッケージの追加ライブラリをインストールします。(F) トーチジオメトリック パッケージをインストールします。 この図の拡大版をご覧になるには、ここをクリックしてください。

図4:HealthModelを実行してmqTrans機能を取得します。 (A)コードをダウンロードします。(B) データファイルの例。各列には調節因子のすべての値があり、最初の項目は遺伝子IDです。各行には特定のサンプルの値が表示され、最初の項目はサンプル名です。(C) ラベルファイルの例。最初の列にはサンプル名が表示され、各サンプルのクラスラベルは label というタイトルの列に示されます。ラベル列の値 0 は、このサンプルが生きていることを意味し、1 は死んでいることを意味します。(D) mqTrans の出力。 この図の拡大版をご覧になるには、ここをクリックしてください。

図5:mqTrans特徴の特徴選択アルゴリズムを実行します。 特徴選択アルゴリズムの結果がユーザーに表示されます。 この図の拡大版をご覧になるには、ここをクリックしてください。

図6:各特徴選択アルゴリズムの最大テストセット精度。 横軸は特徴選択アルゴリズムをリストし、縦軸は精度の値を示します。ヒストグラムは、mqTrans、mRNA、mRNA+mqTransの3つの設定の実験データを示しています。 この図の拡大版をご覧になるには、ここをクリックしてください。

図7:mqTransビューでp値が最小の上位50のダークバイオマーカー。 「ダークバイオマーカー」の列には、ダークバイオマーカーの名前が表示されます。列「mRNA.P」および「mqTrans.P」は、表現型群と対照群の間の統計的t検定p値です。p 値の背景色は p 値 1.00 (青) と 0.00 (赤) の間で色付けされ、白色は p 値 = 0.05 を表します。 この図の拡大版をご覧になるには、ここをクリックしてください。

図8:The Cancer Genome Atlas(TCGA)に掲載された26のがんの詳細。 列「コホート」と「疾患組織」は、各データセットの患者グループと疾患のある組織を記述します。最後の4つの列は、それぞれ発達段階I、II、III、IVのサンプル数を示しています。 この図の拡大版をご覧になるには、ここをクリックしてください。

図9:26のがんにおけるダークバイオマーカーと従来のバイオマーカーの数。 横軸には、26種類のがんがリストされています。縦軸は、これらのがんタイプのダークバイオマーカーと従来のバイオマーカーの数を示しています。 この図の拡大版をご覧になるには、ここをクリックしてください。
補足コーディングファイル1:HealthModel-mqTrans-v1-00.tar このファイルをダウンロードするには、ここをクリックしてください。
プロトコルのセクション 2 (事前トレーニング済みの HealthModel を使用して mqTrans 特徴を生成する) は、このプロトコル内で最も重要なステップです。セクション 1 で計算作業環境を準備した後、セクション 2 では、事前学習済みの大規模参照モデルに基づいて、トランスクリプトームデータセットの mqTrans ビューを生成します。セクション3は、バイオマーカーの検出と予測タスクのために生成されたmqTrans特徴を選択する実証例です。ユーザーは、独自のツールまたはコードを使用して、このmqTransデータセットに対して他のトランスクリプトーム分析を行うことができます。
元の HealthModel フレームワークでは、26 で説明するように、マルチタスク アーキテクチャを使用して、事前トレーニング済みの HealthModel をさらに改良できます。このプロトコルは、トランスクリプトームデータセットのmqTransビューを生成するための事前トレーニング済み参照モデルの利用に焦点を当てています。
デフォルトの事前トレーニング済み参照モデルは、健康なサンプルで確立されたものであり、原発がんと転移性がんの間の調査など、特定のタスクには適していない可能性があります。また、大規模なトランスクリプトームデータセットでは計算速度も遅くなります。
このプロトコルの重要性は、最も豊富に利用可能なOMICデータタイプ、すなわちトランスクリプトームの補完的なmqTransビューを提供することです。ダークバイオマーカーは、従来のトランスクリプトーム解析では無視されていた未分化発現遺伝子から明らかにすることができます。最近の研究では、合計44の805サンプルからなる3つの独立したコホートに基づいて、転移性結腸がん(mCC)の7つのダークバイオマーカーが検出されました。ダークバイオマーカーは、発現が異なるため、ウェットラボでの調査は限られていました。しかし、検出されたmCCダークバイオマーカーの1つであるYTHDC2は、2を含むタンパク質YTHドメインをコードしており、そのタンパク質レベルは、ヒト胃癌細胞46 および結腸癌47の転移状態と正の相関があることが観察された。ダークバイオマーカーの新たな生物学的知見は、in vitroおよびin vivo技術によって解決されるべきものです。
このプロトコルは、完全にモジュール化するように設計されています。原発がんなど、他の大規模なデータセットで事前にトレーニングされた参照モデルは、腫瘍転移の調査を容易にします。このプロトコルは、植物、菌類、微生物など、他の生命領域への応用についても検討されます。
このプロトコルの計算効率は、並列化とアルゴリズムの最適化によって強化される予定です。
このプロトコルでは、トランスクリプトームデータセットを新しいmqTransビューに変換する手順を説明し、遺伝子の変換されたmqTrans値は、参照サンプルと比較して転写規則の変化を定量的に測定します。既定のモデルは、正常なトランスクリプトームで事前にトレーニングされ、参照 HealthModel としてリリースされました。
生物医学研究者がこのプロトコルを簡単に利用できるように、2つのダウンストリームタスクのソースコードが提供されています。実験データは、変換されたmqTrans特徴が、元の発現レベルのみを使用して予測タスクを改善できることを示しています。また、mqTransビューは、元のトランスクリプトームデータに異なる発現を伴わない、いくつかのダークバイオマーカーの潜在的な表現型接続を明らかにすることもできます。
著者は何も開示していません。
この研究は、シニアおよびジュニア技術革新チーム(20210509055RQ)、貴州省科学技術プロジェクト(ZK2023-297)、貴州省衛生委員会科学技術基金会(gzwkj2023-565)、吉林省教育部科学技術プロジェクト(JJKH20220245KJおよびJJKH20220226SK)、中国国家自然科学基金会(U19A2061)、吉林省ビッグデータインテリジェントコンピューティング重点実験室の支援を受けました(20180622002JC)、およびJLU中央大学基盤研究費。このプロトコルの厳密さと明瞭さを大幅に向上させるのに役立った建設的な批評について、査読編集者と3人の匿名の査読者に心から感謝します。
| Name | Company | Catalog Number | Comments |
|---|---|---|---|
| Anaconda | Anaconda | バージョン 2020.11 | Python プログラミング プラットフォーム |
| コンピューター | N/A | N/A | 任意の汎用コンピューターが要件を満たします |
| GPU カード | N/ | A N/A | CUDA コンピューティング ライブラリを備えた任意の汎用 GPU カード |
| pytorch | バージョン 1.13.1 | ソフトウェア | |
| torch-geometric | Pytorch | バージョン2.2.0 | ソフトウェア |
Request permission to reuse the text or figures of this JoVE article
Request Permission