JP6883584B2

JP6883584B2 - マルチオミック癌プロファイルを用いて機能的な患者固有の体細胞異常を識別するための統合された方法及びシステム

Info

Publication number: JP6883584B2
Application number: JP2018530190A
Authority: JP
Inventors: ラジ，アボルファズル; ヴァラダン，ヴィナイ; ディミトロヴァ，ネヴェンカ; バネルジェー，ニランジャナ
Original assignee: Koninklijke Philips NV
Current assignee: Koninklijke Philips NV
Priority date: 2015-08-27
Filing date: 2016-08-26
Publication date: 2021-06-09
Anticipated expiration: 2036-08-26
Also published as: WO2017033154A1; CN108292326B; EP3341875A1; US20180247010A1; CN108292326A; JP2018532214A

Description

関連出願
本出願は、2015年8月27日に出願された米国仮出願第62/210,502号に対する優先権を主張し、その全体が参照により本明細書に具体的に組み込まれる。

技術分野
本発明は、RNAseqベースの発現データ、コピー数多型（CNV：copy number variation）データ及びDNAメチル化データのように共同キュレートされた生物学的経路ネットワーク情報及びオミックスデータを組み込む、遺伝子間調節影響ネットワーク（gene-gene regulatory influence network）を構築し、RNAseqベースの遺伝子発現、アレイベースのDNAメチル化（エピジェネティック）及びSNPアレイベースの体細胞コピー数変異（sCNA：somatic copy-number alterations）を含む、マルチオミック（multi-omic）患者固有測定値と比較することにより、患者固有の遺伝子発現予測を提供するためのデータ駆動型の統合システム及び方法に関する。より具体的には、患者固有の遺伝子発現予測は、癌及び癌治療に関連する予測情報を提供するために、個々の患者サンプルにおいて期待されるレベルから、遺伝子発現レベルの有意な偏差及び不一致を識別するために使用される。

癌の病理生物学は、正常細胞の増殖及び分化を支配する自然の複雑な生物学的プロセス内の有意な異常に関連付けられる。しかしながら、同じ組織タイプ内で発症する癌の中でも、正常なシグナリングネットワークが病理的に変更される可能性がある複数の方法を反映して、有意な異種性が存在する。この異種性は、診断及びセラノスティックのバイオマーカーの開発並びに腫瘍学における潜在的な治療介入において生じる重要な課題の基礎となり、癌の原因及び進行のシステムレベルの理解の必要性を指している。

例えば受容体チロシンキナーゼの上皮成長因子（EGF：epidermal growth factor）受容体ファミリーのメンバーをエンコードし、細胞増殖において重要な役割を果たすERBB2遺伝子は、複数の癌、特に乳癌、胃腸癌及び卵巣癌において高度に過剰発現される。この遺伝子は、乳癌の約20％で調節解除（deregulate）されており、ほとんどの場合、その過剰発現はコピー数重複（copy number amplifications）と関連付けられ、この遺伝子の後に命名される乳癌の特定のサブタイプ、すなわち、HER2陽性乳癌の定義をもたらすことになった。乳癌のこの特定のサブタイプを対象とする治療介入、すなわち、ハーセプチンが利用可能であるにもかかわらず、この治療に対する乳癌患者の反応率は50-55%の範囲にとどまっている。反応におけるこの異種性は、腫瘍進行の他の遺伝子変調因子の存在を指している。実際、腫瘍サプレッサー遺伝子PTENの欠失（deletions）及びPIK3CA遺伝子の突然変異といったAKT/PI3K経路における異常は、ハーセプチンに対する耐性を生じることが分かっている。しかしながら、これらの因子のすべてを治療耐性についての単一の統合的バイオマーカーに統合することができるシステムレベルの経路モデルは現在のところ存在しない。

既知の癌ドライバー遺伝子における固有の反復突然変異の腫瘍化効果は十分に特徴付けられているが、癌全体にわたって観察される大多数の反復突然変異の機能的関連性についてはあまり知られていない。突然変異の機能的関連性を評価するための計算方法は、タンパク質構造に対するそれらの影響を推定することに大きく依存するか、あるいはバックグラウンド突然変異プロセスと比較して、それらの発生の相対的頻度に基づく。下流の細胞プロセスに対する突然変異の潜在的な影響を明らかにするために、近年のアプローチでは、癌サンプル内のマルチオミックス測定値を、共同キュレートされた生物学的経路ネットワークと統合することによって、ゲノム異常の機能的影響を識別することを試みている。しかしながら、大部分のこれらのアプローチでは、下流の遺伝子転写及び経路相互作用の組織特異性に対する複数の調節因子の不均一の、場合によっては非線形の影響を含む、重要な生物学的考察を見落とす傾向がある。

癌サンプルにおける突然変異又はゲノム異常の機能的有意性を評価するために、いくつかの計算フレームワークが開発されている。タンパク質構造に対する突然変異の効果についての推論に基づく方法が社会で広く用いられているが、最近の研究では、バックグラウンドの突然変異プロセスと比べて、突然変異している遺伝子の相対的頻度を評価することによって、遺伝子のドライバー突然変異を決定することに焦点を当てている。サイレント突然変異は典型的にいずれかの候補遺伝子については珍しく、その結果、不正確なバックグラウンドの突然変異率推定値につながる可能性があることを認識し、MutSigCVは、バックグラウンド突然変異率推定値を改善するために、候補遺伝子に対して類似するゲノム特性を有する遺伝子を利用しようと試みた。他の方法は、所与の癌サブタイプ内の体細胞突然変異によって頻繁に当てられるサブネットワークを識別することを目的とする。しかしながら、これらのアプローチは、体細胞異常の下流の調節解除又はシグナリング効果への機構的洞察を提供しない。これらの欠点は、細胞エンティティ（例えば遺伝子、RNA、タンパク質、タンパク質複合体及びmiRNA）間の十分にキュレートされた生物学的相互作用が、経路ネットワークに関してモデルに組み込まれるネットワークベースの方法につながっている。他の研究は、癌の臨床転帰と、遺伝子やタンパク質発現レベルといった分子エンティティの活性化レベルとの間の関連にのみ焦点を当てているが、癌生物学における突然変異の機能的効果を明示的にモデル化していない。最近、経路調節ネットワークをマルチオミックスデータと統合して、経路内の個々のノードの活性に対する体細胞突然変異の機能的影響をモデル化するPARADIGM-SHIFTが提案された。任意の所与のタンパク質における体細胞異常の機能的効果は、その上流調節ネットワークから一旦取得され、かつその下流の標的ノードから再び取得される、対応するノードの活性の差に基づいて推論される。

発展は異なるが、これらの方法には、生物経路ネットワークへの絶対的な依拠という、共有される共通の欠点が存在し、したがって、これらの方法の利用は、十分にキュレートされた経路ネットワークに限定されるべきであり、異なる組織状況から導出された、部分的に検証されたネットワーク又は分子ネットワークについては推奨されない。より重要なことに、これらの技術は、典型的に、すべての親遺伝子が対応する相互作用に等しく寄与していると推定し、ネットワークノード間の相互作用の間の影響の強さの変動の可能性を見落すことになる。例えば複数の遺伝子が特定の標的遺伝子の転写調節因子として現れる場合、これらの遺伝子は生物学的に疑わしい標的遺伝子の発現レベルに等しく寄与すると考えられる。実際には、隣接するノード間のペアワイズの影響は非常に異なる可能性がある。HotNetアルゴリズムではリンク間の異種性が考慮されており、HotNetアルゴリズムは、ネットワークトポロジに基づく遺伝子ペア間のペアワイズ影響の尺度を定義することを通して、この異種性を発見することを意図する。しかしながら、基礎となる複雑な調節相互作用から生じる実際のペアワイズの影響異種性は、推定経路ネットワークトポロジから完全には抽出することができない。

経路レベルの異常は、体細胞突然変異、コピー数変異、エピジェネティック変動（epigenetic variation）及び調節遺伝子発現の変化といった複数の源から生じる可能性があるので、これらの変動源を一緒にモデル化することが、腫瘍学での用途について包括的な経路ベースの予測モデルを開発するために不可欠である。さらに、分子生物学における低コストゲノムワイドデータ取得技術の最近の進歩により、様々な変動源の測定がますます利用可能になってきている。しかしながら、これらのマルチオミックスプロファイルに存在する情報を十分に活用することができるモデリングフレームワークは、研究コミュニティと診断コミュニティの双方で欠けている。したがって、臨床的に有用なバイオマーカーを見つける目的では、RNA発現レベル、コピー数多型、DNAメチル化パターン及び体細胞突然変異を含む様々なデータ源を統合するための計算フレームワークの開発は、腫瘍学のコミュニティにおいて本質的に必要である。

最近、早期癌診断、臨床結果予測及びより関連性の高い治療介入を容易にするために、様々な情報源を統一されたフレームワークに組み込むためのいくつかの統合アプローチが提案された。これらのアプローチの大部分は、いずれかの２つの極端な視点：すなわち、ｉ）概念的な生物学的情報を完全に無視し、データ駆動技術のみに依存すること、又はｉｉ）相互作用する分子エンティティのネットワークを組み込むことにより概念的な生物学的情報を完全に信頼すること、という２つの視点のうちの１つをとる。第１のアプローチでは、細胞分子エンティティ（例えば遺伝子及びタンパク質）の間の生物学的相互作用を無視することは、データの過剰適合の可能性に起因して、有意な集合的予測力を有するエンティティの生物学的に関連するサブセットを見つける際に非常に非効率的である。実際、任意の所与の研究では癌サンプルの数が、測定された分子特徴の数よりも数桁少ない傾向があるので、この問題は、癌研究において特に顕著である。一方、記述的な生物ネットワークへの完全な依存は、それらのネットワークの限界を無視する：経路ネットワークは、典型的に、必ずしも他の組織及び病理学的状況に翻訳可能ではない特有の細胞状況における経験的証拠に基づいて構築される。

本発明内では、ハイブリッドアプローチをとり、測定ベースのオミックスデータと部分的に信頼できる経路情報の双方を統一フレームワークに組み込んで、遺伝子間の影響ネットワークを構築する。このような遺伝子間影響ネットワークは、調節ネットワーク状態を所与として特定の遺伝子発現レベルを予測することができる。このフレームワークは、組織固有のタンパク質間の相互作用についての我々の知識を洗練及び拡大するだけでなく、ネットワークエンティティ（例えば遺伝子）の患者固有の予測及び条件付き分布も提供する。これらの患者固有の遺伝子発現予測を利用して、個々の患者サンプル内において、期待されるレベルからの遺伝子発現レベルの有意な偏差及び不一致を見つけ、したがって、治療反応及び予後のような表現型との潜在的関連の発見を可能にすることができる。

本発明は、生物学的情報及び様々な分子測定データ源を統一されたネットワークベースの計算フレームワークに統合する際のいくつかの重要な限界を克服する。これは、より関連性の高い患者固有の機能不全遺伝子及び攪乱された生物学的プロセス（perturbed biological processes）を明らかにすることにつながる。

例えば本発明の方法は、生物学的情報を組み込み、基礎となるネットワークベースの予測及び患者固有の測定と有意な不一致を示す遺伝子のみを報告する。したがって、このアプローチは、考慮すべき表現型に関連付けられる最も機能的に関連する遺伝子を識別する際のより高い特定及び感度をもたらす。

また、現在のセットベースの方法は、事前の生物学的知識に基づいて特定の表現型又は細胞／生物学的プロセスに一緒に関連付けられる遺伝子のセットに最初に注釈を付けることによって、生物学的情報を考慮に入れる。しかしながら、セットベースの方法は適応的な統合を行うことができず、ユーザは、潜在的により関連性のある遺伝子セットを形成することにより手動で生物学的情報を含めることを必要とされる。対照的に、本発明では、癌生物学に関する事前情報をなんら必要としない。この方法は、経路ネットワークの注釈から各遺伝子の遺伝子調節ネットワークを開発する。結果として得られる表現型に関連付けられる経路サブネットワークは、ロバストなバイオマーカーとともに機能的な洞察を提供し、したがって、癌全体にわたって広く適用可能である。

Paradigm、Pathologist、SPIAといった現在利用可能なネットワークベースの方法は、攪乱された経路やネットワークから得られた予測からの有意な偏差を示す遺伝子を識別するために、経路情報を測定データと統合することを目的とする。これらのアプローチには２つの重要な欠点がある。第１に、これらのアプローチは、経路ネットワーク接続性における組織固有の変異の可能性を許容することなく、生物学的経路ネットワーク関係を完全に信頼する。第２の更に重要な問題は、これらの技術がネットワーク内の相互作用リンクの間の機能的な異種性の可能性を見落とすことである。実際には、いくつかの調節親遺伝子の影響は、他の親遺伝子よりも有意に高いことがあるが、これらは直接の親ノードすべてについて等しい影響を前提とする。

方法及びシステム内では、経路ネットワークに完全に依拠するのではなく、マルチオミックスデータから学習されるネットワークエッジに異なる係数を割り当てることによって影響ネットワークを精緻化する。例えばテーブル２及びテーブル３を参照されたい。上流調節因子を表すネットワークエッジは、祖先の係数を使用して捕捉され、シス調節の影響（cis-regulatory influences）はCNV及びメチル化係数として捕捉される。さらに、緩く接続されたリンクは除去される。したがって、我々の方法は、ネットワークノード（例えば遺伝子、RNA、タンパク質）間の異種の関係を強調し、発見する。

さらに対照的に、我々の方法は、生物学的経路とマルチオミックス測定データの双方を使用して、トポロジだけでなく、上述のようなネットワーク内のノード間の影響の強さも捕捉する。したがって、ネットワークノード間でより正確かつ現実的な影響を提供する。第２に、方法は、体細胞突然変異によって頻繁に影響される経路を見つけることに限定されず、機能不全のノードも見つける。

これらの問題に対処するために、突然変異によって影響を受ける情報フロー（「InFlo-Mut」：Information Flow impacted by Mutations）と呼ぶ本発明のプロセスは、RNAseqベースの遺伝子発現、アレイベースのDNAメチル化（エピジェネティック）及びSNPアレイベースの体細胞コピー数変異（sCNA）を含むマルチオミックス測定値及び生物学的経路ネットワーク情報を組み込み、遺伝子間調節影響ネットワークを構築する。InFlo-Mutは、正常及び癌組織の分子プロファイルからの標的遺伝子に対する調節ノードのペアワイズ影響を学習する。新たなサンプルのノードの活性を推論するために、InFlo-Mutは、トレーニングデータセットから既に学習されたネットワーク係数を使用する。これは、非線形ベイジアンモデルを学習することを通して実現され、生物学的経路ネットワークから推論される上流調節の影響とともに、それ自体のsCNA及びメチル化プロファイルを使用して任意の所与の遺伝子の発現レベルを予測する。このアプローチは、異種ペアワイズ影響係数を捕捉することにより不均一な親ノードの寄与の問題を解決するだけでなく、ノード間の非線形関係を学習することも可能である。InFlo-Mutは、突然変異が標的遺伝子の調節不全に与える影響がより高い遺伝子のサブセットを明らかにするために、下流の標的遺伝子との体細胞突然変異の間の関連の評価も可能にする。InFlo-Mutのロバスト性と生物学的妥当性を乳癌と結腸癌の２つの大きなマルチオミックスデータセットに適用することにより、そのロバスト性と生物学的妥当性を実証し、これらの疾患の主要な発癌経路における遺伝子間の突然変異の潜在的な調整効果を発見する。

特に、本発明の目的は、マルチオミック生物学的情報及び様々な分子測定データ源を伴うキュレートされた経路ネットワークを、統一されたネットワークベースの計算フレームワークに統合して体細胞突然変異の影響を識別することにより、従来技術の上記の問題を解決するシステム及び方法を提供することである。本発明の目的はまた、患者固有の遺伝子発現予測を提供し、予測されたレベルからの患者遺伝子発現レベルの有意な偏差及び不一致を識別し、より関連する機能不全遺伝子及び攪乱された生物学的プロセスを識別するためのシステム及び方法を提供することである。本発明の更なる目的は、治療応答及び予後等の表現型との潜在的関連を識別することである。本発明の目的はまた、先行技術の代替物を提供することである。

したがって、上記の目的及び他のいくつかの目的は、調節不全遺伝子を引き起こす潜在的な体細胞異常を識別及び報告するためのシステム及び方法を提供することにより、本発明の第１の態様において得られるよう意図されている。そのような方法は：
十分にキュレートされた公的に利用可能な経路ネットワークから生物学的ネットワーク経路情報を取得し、該経路情報を受け取るように構成されたプロセッサ上へその経路情報を入力することにより、関心のある各特定の標的遺伝子についての上流調節親遺伝子情報（upstream regulatory parent gene information）の第１データセットを決定するステップと；
適用することにより、各々の特定の標的遺伝子について、該遺伝子の発現レベルとそれ自身のゲノム及びエピジェネティック状態（epigenetic status）、並びにその上流転写調節因子（upstream transcriptional regulator）の間の関係を捕捉する調節ツリー（regulatory tree）を決定するステップであって、関心のある遺伝子は、ルートノードに存在し、ツリーの葉は、直接又は中間のシグナリングパートナーを介して間接的にその転写を潜在的に調節する遺伝子のすべてを表す、ステップと；
RNAseq発現データ、コピー数多型データ及びDNAメチル化データのような測定ベースのオミックスデータ（measurement-based omics data）の第２データセットを決定し、測定ベースのオミックスデータを、当該データを受け取るように構成されたプロセッサ上へ入力するステップと；
コンピュータによって演算技術を適用し、特定の遺伝子発現レベルを調節ツリーの葉に関連する測定値に関連付けるために、遺伝子のエピジェネティック情報及び調節ネットワーク状態に基づいて、関心のある各遺伝子について非線形関数を学習するステップであって、非線形関数のパラメータは、ツリー内のルートノードにより近いノードの潜在的により強い調節影響を捕捉するために新規な深さペナリゼーション機構（novel depth penalization mechanism）を組み込んだベイジアン推論法を使用して推定される、ステップと；
コンピュータによって分析技術を適用して、関心のある各遺伝子の発現レベルを予測するステップと；
所望の標的遺伝子の観察された発現レベルに関連する患者固有の情報を決定し、該患者固有の情報を第３データセットに入力するステップであって、患者固有の情報は、RNA発現データ、CNVデータ、メチル化データ及び体細胞突然異常データ等の新たな癌サンプルデータを含む、ステップと；
患者固有の情報及び予測発現レベル情報を使用して、所与のサンプル内の所望の標的遺伝子について予測された発現レベルと観察された発現レベルとの間の相対的な患者固有の不一致スコアを計算するステップと；
すべてのテストサンプルについて取得された活性化及び不一致スコアを評価して、標的遺伝子の発現レベルの不一致とその特定の遺伝子の上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連を発見するステップと；
を含む。

本発明の第２の態様によれば、個々の患者サンプルにおける標的遺伝子発現レベルにおける不一致と上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連を利用して、患者固有のバイオマーカーを識別するためのシステムが提供され、当該システムは、遺伝子発現レベルの有意な偏差及び不一致を識別するために統合された統一ネットワークを備え；
十分にキュレートされた生物学的ネットワーク経路情報から得られる関心のある各特定の標的遺伝子についての上流調節親遺伝子情報の第１データセットであって、そのような経路情報を受け取るように構成されたプロセッサに含まれる第１データセットと；
各々の特定の標的遺伝子について、該標的遺伝子の発現レベルとその標的遺伝子自体のゲノム及びエピジェネティック状態、並びにその上流転写調節因子の間の関係を捕捉する調節ツリーであって、関心のある遺伝子が、ルートノードに存在し、ツリーの葉は、直接又は中間のシグナリングパートナーを介して間接的にその転写を潜在的に調節する遺伝子のすべてを表す、調節ツリーと；
RNAseq発現データ、コピー数多型データ及びDNAメチル化データのような測定ベースのオミックスデータの第２データセットであって、そのようなデータを受け取るように構成されたプロセッサ上にも配置される第２データセットと；
標的遺伝子のエピジェネティック情報及び調節ネットワーク状態から決定された各標的遺伝子について学習される非線形関数であって、該非線形関数は、その特定の標的遺伝子の発現レベルを調節ツリーに関連する測定値に関連付け、非線形関数のパラメータは、ツリーのルートノードにより近いノードの潜在的により強い調節影響を捕捉するための新規な深さペナリゼーション機構を組み込んだベイジアン推論法を使用して推定される、非線形関数と；
標的遺伝子の観察された発現レベルに関連する患者固有の情報の第３データセットであって、患者固有の情報は、RNA発現データ、CNVデータ、メチル化データ及び体細胞突然変異データ等の新しい癌サンプルデータを含む、第３データセットと；
を備え、
標的遺伝子の発現レベルは非線形関数を利用して決定され、所与のサンプル内の標的遺伝子について予測された発現レベルと観察された発現レベルとの間の相対的な患者固有の不一致スコアを決定し；
標的遺伝子の観察された発現レベルに関連する患者固有の情報の第３データセットについて、活性化及び不一致スコアを決定し、患者固有の情報は、RNA発現データ、CNVデータ、メチル化データ及び体細胞突然変異データ等の新しい癌サンプルデータを含み、
標的遺伝子の発現レベルは非線形関数を利用して決定され、所与のサンプル内の標的遺伝子について予測された発現レベルと観察された発現レベルとの間の相対的な患者固有の不一致スコアを決定し；
すべてのテストサンプルについて取得された活性化及び不一致スコアを評価して、これにより、標的遺伝子の発現レベルの不一致とその特定の遺伝子の上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連を識別する。

本発明による方法を、添付の図面に関してより詳細に説明する。図面は、本発明を実施する方法を示しており、添付の特許請求の範囲内に入る他の可能な実施形態に限定するものとして解釈されるべきではない。

患者固有の遺伝子発現予測を提供するために、遺伝子調節及び／又はシグナリング経路ネットワークを、測定ベースのオミックスデータと統合するステップの経路を説明する方法の概要を示す図である。本発明のこの態様のステップは、ｉ）単離されていない標的遺伝子についての調節ツリーを抽出すること、ｉｉ）トレーニングデータセットを使用して各標的遺伝子について非線形関数を学習すること、ｉｉｉ）標的遺伝子の遺伝子発現値を予測し、活性化及び一貫性スコアを計算すること、そしてｉｖ）機能的突然変異の影響分析である。サンプル遺伝子PPP3CAについて、経路データベースから導出される調節相互作用を使用して生成された調節ツリーを示す図である。経路ネットワーク内のすべての遺伝子についてレベル２までの祖先の数の分布を示しており、ほとんどの遺伝子が１０個〜５０個の間のどこかの調節因子を有することを例示する、遺伝子の祖先カウントのヒストグラムを示す図である。２つの潜在的な非線形効果を捕捉するために、中央シグモイド（centered sigmoid）及びソフト閾値処理（soft thresholding）を含む非線形関数のグラフを示す図であり、２つの潜在的な非線形効果とは、ｉ）近平均感度（near mean-sensitivity）及びｉｉ）近平均無視（near-mean ignorance）であり、ｘ軸は測定されたコピー数又はDNAメチル化レベルを示し、ｙ軸は遺伝子発現に対する測定の影響の程度を示す。近平均感度の場合、平均値近くの測定されたDNAメチル化の小さな変化は、遺伝子発現の大きな偏差をもたらす。しかし、近平均無視では、平均近くのコピー数の小さな変化は遺伝子発現における大きな変化をもたらさない。 CRC正常及び腫瘍サンプルの観察についてのJUN遺伝子発現レベルの予測を示す図である。癌サンプル（*赤）は、正常サンプル（*青）と比べて広範な不一致を示す。方法予測は、誤差のバー┬によって示される最大３つの標準偏差までの事後平均（posterior mean）（o）及び信頼区間に関して提供される。 BRC及びCRC腫瘍サンプルについてのすべての遺伝子の不一致スコアを示す図である。 BRC及びCRC腫瘍サンプルについてのすべての遺伝子の不一致スコアを示す図である。 BRC及びCRC腫瘍サンプルについてのすべての遺伝子の不一致スコアを示す図である。 BRC及びCRC腫瘍サンプルについてのすべての遺伝子の不一致スコアを示す図である。ネットワークベースの予測と患者固有の測定との間の有意な不一致に基づいて、患者固有の機能不全遺伝子を識別するための本発明の方法を要約したフローチャートである。結腸癌サンプルにおける標的遺伝子発現に対する体細胞突然変異の影響を示す本発明の方法の結果のグラフ表示である。遺伝子PTENについてのRNA発現のヒストグラムである。試料サンプルMYB、GATA3、PTEN及びERBB2についての予測対観察を示す図である。遺伝子ERBB2についてのRNA発現レベル対コピー数多型CNVを示す図である。遺伝子発現の不一致に対する、PTENの上流調節サブネットワークにおける体細胞突然変異の影響を示す図である。

本発明は、患者固有の遺伝子発現予測を提供し、期待されるレベルからの遺伝子発現レベルの有意な偏差及び不一致を識別するために、マルチオミック生物学的情報及び様々な分子測定データソースを統合ネットワークベースの計算方法に統合するためのシステム及び方法を提供する。本発明は、図１〜図１２を参照して以下に更に詳細に説明される。

本発明の一実施形態によると、患者固有の遺伝子発現予測を提供し、期待されるレベルからの遺伝子発現レベルの有意な偏差及び不一致を識別し、患者固有のバイオマーカーを報告するための方法の全体的なブロック図を提示するフローチャートが、図１に概説されるステップ又はモジュールによって説明される。図１に図示されるように、この方法は、調節不全遺伝子（dysregulated genes）を引き起こす潜在的な体細胞異常を識別して報告するための４つの主要な逐次的ステップ又はモジュールからなる。第１ステップのモジュール１では、調節ツリーが、関心のある各遺伝子について、それ自体のゲノム及びエピジェネティック状態を有する遺伝子の発現レベル並びにその上流転写調節因子の間の関係を捕捉する経路ネットワークから抽出される。関心のある遺伝子はツリーのルートノードに存在し、ツリーは遺伝子の転写の上流調節因子のネットワークを表す。ツリーの葉は、直接的又は間接的に中間シグナリングパートナーを通じて遺伝子の転写を潜在的に調節する遺伝子のすべてを表す。我々は、「祖先遺伝子」又は単に「祖先」という用語を使用してこれらの遺伝子を参照する。

モジュール２である第２ステップでは、特定の遺伝子発現レベルを調節ツリーの葉に関連する測定値に関連付けるために、各遺伝子について非線形関数を決定する。したがって、各ツリーサブネットワークは、非線形関数を学習して、対応する遺伝子発現レベルを、その遺伝子自体のエピジェネティック情報（例えばDNAメチル化及びコピー数）及びその調節祖先遺伝子発現レベルから予測するために使用される。非線形関数のパラメータは、ツリーのルートノードにより近いノードの潜在的により強い調節影響を捕捉するための新規な深さペナリゼーション機構を組み込んだベイジアン推論法を使用して推定される。これは、特定の組織タイプの状況において特定の遺伝子に各々対応する機能のバンク（a bank of functions）を提供する。この機能データベースを一度学習し、モジュール３及び４によって実行される２つの後続のステップにおいて患者固有の分析のために使用することができる。

ステップ３において、モジュール３は、所与のサンプル内の所望の標的遺伝子についての予測発現レベルと観察発現レベルとの間の相対的な患者固有の不一致スコアを計算する。すなわち、モジュール３は、所与の患者についての情報を受け取り、機能バンクを使用して調節ネットワーク内のすべての遺伝子について遺伝子発現レベルの予測を実行する。このモジュールは更に、遺伝子発現の実際の測定値、すなわち観察値を予測値と比較することにより、各遺伝子についての一貫性スコアを計算する。第４ステップでは、モジュール４は、標的遺伝子発現レベルの不一致と、その特定の遺伝子の上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連を発見するために、すべてのテストサンプルについて得られた活性化及び不一致スコアを評価する。したがって、モジュール４は、その発現レベルが調節ネットワークから得られた予測値と有意に不一致である遺伝子を識別する。これらの遺伝子は、遺伝子内のコピー数異常又はその祖先の体細胞突然変異に起因して機能不全である可能性が高い。モジュール４は更に、子遺伝子発現レベルの不一致と潜在的に関連付けられる祖先遺伝子突然変異の有意性を評価するための統計を提供する。

モジュール１：経路ネットワークの組込み−調節ツリー構築
遺伝子転写は複雑な生物学的プロセスであり、生物学的経路データベースに注釈を付けるように、複数の相互作用するタンパク質及び複合体、並びにDNAメチル化の程度及び保持しているDNAセグメントのコピー数を通して異なるレベルで調節される。経路ネットワークは、細胞内相互作用及び遺伝子調節ネットワークをネットワーク形式で提示するために広く使用されている。ネットワークはノードとエッジの有向グラフを構築する。ノードは、遺伝子、タンパク質、RNA、miRNA、タンパク質複合体、シグナル受容体のよう多様な範囲のエンティティから、更にはアポトーシス、減数分裂、有糸分裂及び細胞増殖といった抽象的プロセスからも構成され得る。ネットワークエッジは、相互作用するノードのペアを決定し、各相互作用のタイプを指定する。いくつかの公に入手可能な経路ネットワークが、様々な種と組織タイプとの間の細胞内活動をモデル化するために開発されている。

本発明では、NCI-PID、Biocarta及びReactomeを含む様々な十分にキュレートされた経路源からの経路をまとめた包括的ネットワークを使用する。この「スーパー経路ネットワーク（super pathway network）」は、タンパク質又は対応する遺伝子、RNA、タンパク質複合体、遺伝子ファミリー、miRNA及び抽象概念（abstracts）を含む６つのノードタイプからなる。これらのノードは、i）陽性転写、ii）陰性転写、iii）陽性活性化、iv）陰性活性化、v）遺伝子ファミリーメンバーシップ及びvi）タンパク質複合体の成分にする、という６つの異なる方法で相互に作用する。通常、転写は、対応するタンパク質によって表される遺伝子に対してのみ終了し、一方、活性化はすべてのノードタイプに適用可能である。

遺伝子のmRNA発現レベルを、そのエピジェネティックパラメータ（DNAメチル化及びコピー数多型）並びにその調節ネットワークに関連付ける機能を学習するために、スーパー経路ネットワークデータベースから各遺伝子の調節ネットワークを抽出し、それを「ツリー」（図２）として表す。その後、調節ツリーを形成するすべてのノードの影響を集合的に捕捉する調節因子又は調節遺伝子と呼ばれる「調節祖先遺伝子」のリストを抽出する。調節因子の一部は、標的遺伝子の直接的な親であり、よって、その転写を直接調節するが、一方、他の調節因子は、タンパク質複合体及び直接調節因子の翻訳後修飾を介して間接的に標的遺伝子発現に影響する。

各遺伝子について調節ツリーを展開する際に、特定の標的遺伝子から出発し、上流ネットワークをリンクの反対方向にトラバース（traverse）してすべての上流ノードを収集し、その深さとともに調節遺伝子を捕捉する。その深さは、遺伝子転写調節の生物学に基づくいくつか修正、及び調節ネットワークに参加する他の遺伝子の発現を使用して標的遺伝子発現を予測することに関心を持つという事実とともに、周知の深さ優先探索（以下の疑似コードを参照されたい）のような深さ優先トラバースアルゴリズムを使用して、図２に図示されるように、ルートノードへのリンクの数として定義される。

最初に、事前定義された最大深さレベルに達すると、分岐をトラバースすることを終了する。この場合、深さは、訪問しているノードからルートノードまでのリンクの数として定義される。次に、遺伝子ノードで終了しない分岐をすべて取り除き、したがって、ツリーの葉は常に遺伝子である。不必要なネットワークの複雑性と不適切な相互作用の包含を避けるために、概念的抽象プロセスを表す抽象ノードを除くすべてのノードを通過する。遺伝子ノードに到達している間、「転写」タイプではないリンクを通過するだけである。なぜなら、「転写」リンクを介して遺伝子ノードで終わる上流調節ネットワークの一部は、この特定の遺伝子レベルを考慮することにより既に説明されているためである。このルールの唯一の例外はルートノードであり、この場合、以下のように正反対の手順を行う。

ルート近傍（root neighborhood）の最初の環においてルートノードから直接近傍（direct neighbors）まで通過することは、接続エッジが「転写」タイプである場合にのみ許可され、親を、ツリールート内に存在する遺伝子の発現レベルに影響を与えるものに限定する。また、葉から、機能学習プロセスで更に使用されるルートノードまでの距離も追跡する。最終的に、２つの分離した経路を介してノードに届く場合、最も短い経路が考慮される。モジュール１の選択プロセスのための擬似コードを以下に要約し、ネットワークから遺伝子PPP3CAについて抽出されたサンプル上流ツリーを図２に示す。

図２は、サンプル遺伝子PPP3CAの経路データベースから導出される調節相互作用を使用して生成された調節ツリーの例である。サブネットワークは、第３レベルまでの深さ１の祖先遺伝子を含む。形状はノード型を定義しており、遺伝子（楕円形）、タンパク質複合体（矩形）、遺伝子ファミリー（五角形）、抽象概念（ひし形）である。エッジは、その調節機能に従って色付けされており、正の活性化（黄色）、負の活性化（赤色）、陽性転写（緑色）、陰性転写（青色）、タンパク質複合体の成分（黒色）及び遺伝子ファミリーメンバー（グレー）である。追加の調節親（additional regulatory parents）と考えられるルートノードのエピジェネティック及びsCNA測定値（角丸矩形）は、緑色の矢印で接続される。調節因子は最大でレベル３（d_max=3）まで収集される。ルートノードPPP3CAの第１レベルの祖先（直接の親）は、遺伝子発現レベルを調節する「転写」エッジを介して接続されるように示される。例えば複合CAM/Ca++は活性化リンクを介してルートノードに接続されており、よって遺伝子発現レベルを調節しない。したがって、図２の左側の複合CAM/Ca++を介して接続している遺伝子はすべて、最終的な祖先リストから除外される。他の遺伝子を通過している間、非転写リンクのみが許容される。例えばMYBの上流サブネットワークは、PIAS3及びMAP3K7遺伝子のような非転写ノードに限定され、その影響は、MYB発現レベルを介してまだ捕捉されていない。遺伝子GATA3及びE2F1の影響は、遺伝子MYBの発現レベルによって暗示的に説明される。

例として、図３に、ルートノードの上流の最大７つまでのリンクをトラバースするときの祖先の数の経験的分布が、対数スケールで提示されている。多数の遺伝子が上流の単離されたオーファン遺伝子（upstream isolated orphan genes）である。遺伝子CDKN1Aについて、839の遺伝子のみが、23の遺伝子について１つのみの祖先から最大で1152の祖先までにわたる祖先を有する。祖先がゼロの遺伝子は、経路ネットワーク内には現れなかった。

モジュール２：各遺伝子についての非線形関数の学習
本発明の方法の第２ステップは、ルートノードに存在する遺伝子の発現レベルをその調節ネットワーク及びそれ自体のエピジェネティック情報（例えばDNAメチル化及びCNV）に関連付ける関数（function）を学習することである。関数を「学習する」ことは、標的遺伝子の発現に対する調節遺伝子の発現レベルの影響を定量化することを意味する。また、この方法では、標的遺伝子のモデルをトレーニングする。そのようなモデルでは、（以下のベイジアンモデル推定、特にβ_gを推定する方法で説明されるように）トレーニングデータで観察されるようなペアワイズの影響に基づいて、親遺伝子について異なる係数を割り当てる。複数のDNAメチル化プローブは、遺伝子のコード又は調節領域と重複する可能性があるため、本発明は、最小、最大及び加重平均値といったいくつかの代表的な統計値を含めることによってメチル化測定値を活用する。ここで、加重平均を計算する際には、更なる正確性のために１０プローブ未満の領域は除外する。したがって、遺伝子gが、

領域と重複し、各領域が、プローブ数

を有し、メチル化測定値

に対応する場合、次いで加重平均は次のように計算される：

ここで、I(.)は識別関数である。

コピー数多型を含めるために、本発明は、特定の遺伝子を有する領域に提供されるセグメント平均値を使用する。ほとんどの遺伝子は単一のCNVセグメントに入る。そうではなく、遺伝子が２つのセグメントの境界にある場合、単に双方のセグメント測定の平均値をとる。

各遺伝子について関数を学習するために、モジュール２は、n_gサンプルについて、その祖先のmRNA発現、体細胞コピー数変異及びDNAメチル化測定を使用し、以下の古典的回帰モデルを形成する：

ここで、y_gは、n_gサンプル全体にわたる遺伝子gについての発現レベルのn×1ベクトルである。

は、

（自己メチル化及びCNVデータ）と、

（祖先遺伝子の発現レベル）とを含む２つの部分から構成されるn×pデータ行列である。ここで、

である。項

は、長さn_gのすべて１つの列ベクトルであり、εはi.i.dゼロ平均単位分散ガウス要素（i.i.d zero-mean unit-variance Gaussian elements）を有するモデル雑音である。μ_gは遺伝子gの発現レベルの期待値である。

ここでの目的は、平均二乗誤差（MSE：Mean Squared Error）を最小化することにより最良の予測力を提供する最適モデルパラメータβ_i、i= 1, 2,.., pを見つけることである。１つには、高度に汚染された／不規則な癌細胞内のかなり攪乱した相互作用に起因するモデルクラッシュを回避するために、学習段階で正常サンプルを使用してもよい。しかしながら、これは、予測因子の数がサンプル数（n<O(p)）に対して大きいか又は同等であるときに予測力を弱くする可能性がある。ほとんどの研究では、プロファイルされる癌サンプルの数は、正常サンプルの数よりもかなり多い傾向がある。例えば乳癌のTCGAデータの場合、癌サンプルの数は正常サンプルの１０倍を超える。したがって、すべての癌サンプルを除外することはかなり非効率的である。一方、トレーニングセット内に癌サンプルを含めることは、上記のようなゲノム事象に起因して、一部のサンプルにおいて真の基礎となる生物学的機能から大きく逸脱する特定の遺伝子についてモデル性能を低下させる可能性がある。したがって、我々は、予測関数を学習するために、この特定の遺伝子及びその祖先の体細胞突然変異によって影響を受けていないすべての正常サンプルと癌サンプルの一部を含める。このアプローチは、各遺伝子について異なるトレーニングセットサイズにつながるが、モデル予測力においてかなりの改善を提供する。

最小２乗誤差（LSE：Least Squared Error）の解決策は、モデルパラメータβ_iに関して事前情報が利用可能でないときに、トレーニングセットについて二乗誤差を最小化する。

モデルパラメータに関する事前情報が存在するとき、LSE解決策は最適ではない。ここで、モデルの精度を高めるために使用することができる、モデルに関する事前知識がある。第１に、必ずしもすべての祖先遺伝子が所与の遺伝子の発現レベルの実質的な影響を有しているわけではない可能性が高い。したがって、かなりの数のモデルパラメータβ_iをゼロに向けて小さくすることができる。したがって、スパース性（sparsity）を課すことは、雑音オーバーフィットを回避することによってモデル一般化特性を強化する。スパース性の一部は、経路ネットワークを使用し、入力データとしてすべての遺伝子を使用する代わりに祖先遺伝子のみを含めることにより既に示されているが、祖先遺伝子の数が多くなると（数十及び数百の順で）、また更に高いレベルのスパース性が期待される。

スパース性を課す一般的な最適化ベースの解決策の１つは、モデルパラメータのノルムを正規化することである。ペナリゼーションを、係数ベクトルβ= [β₁,β₂,…,β_p]^TのノルムLp、p≧0に適用することができ、これはブリッジ回帰と呼ばれる。このアプローチの重要な特別なケースは、L、L₂、L₀のノルムペナリゼーションについてのそれぞれLasso、Ridge及びサブセット選択である。弾性ネット（elastic net）では、ペナルティ項は、L₁及びL₂のペナルティの線形結合である；

ここで、λ₁及びλ₂は収縮パラメータであり、スパース性及び一般化可能性（generalizability）を課す。凸最適化、Basis追跡（Basis pursuit）、LARS、座標降下、Dantzigセレクタ、直交マッチング追跡（Orthogonal matching pursuit）及び近似メッセージパッシング（approximate message passing）に基づく効果的なアルゴリズムを用いてこの問題を解決することができる。しかしながら、これらの方法の最も制限的な欠点は、回帰係数の点推定値を提供するだけであるということである。

対照的に、本発明は、後の整合性チェック分析で使用される事後分布（posterior distribution）を介してモデルパラメータに関するより詳細な情報を提供するベイジアンフレームワークを使用する。また、以下で説明されるように、スパース性に加えて他の事前知識の組込みも可能である。

歴史的には、遺伝子発現研究を解析する際に、生物学的測定の間の潜在的に非線形の関係は無視されていた。そのような非線形関係を捕捉するために、本発明のモジュール２は、平均値の周囲の感度を捕捉する中心シグモイド関数

及び
極端に高い又は低い値のみがモデルに寄与する場合を説明するソフト閾値処理関数

を使用する。f₂(x；c)は、一般的に使用される区分別（peace-wise）線形ソフト閾値処理関数 f₂(x；c) = sign(x)(|x|-c)₊のよりソフトなバージョンと考えることができる。これらの関数は、線形関数と対比して図４に示されている。我々は、要素別非線形拡張

のみを自己データ（例えばメチル化及びCNVデータ）に適用し、したがって、予測因子の数は、各遺伝子の祖先の数に比べてわずかに増加する。注目すべきことに、実際の基礎となる関数が線形であれば、非線形項の係数は、提案したモデルでゼロになる傾向があり、したがって、真の線形関係について非線形関数を学習している間は性能の低下は見られない。

経路ネットワークを上方に向かってトラバースすることにより各遺伝子の祖先セットを展開する上で重要な生物学的考慮事項は、ルートノードに対する葉ノードの距離の変動である。１つには、より近い祖先は、中間ノードの長い鎖を介して接続されるより遠いノードよりも、子孫の下流遺伝子発現レベルにより寄与するよう期待することができる。したがって、より近いノードは、回帰モデルにおいてより高い係数を提示する傾向がある。モジュール２は、この事実を、以下で説明されるベイジアンモデルにおいて、

によって捕捉されるように、ベイジアンフレームワークの深さペナリゼーション機構を通じて方法に利用する。

ここで、本発明は、ベイジアンフレームワークを使用して、その自己エピジェネティックデータ並びにその調節祖先遺伝子の発現レベルの非線形変換／投影を介して遺伝子発現レベルを予測する。ベイジアンフレームワークは、モデルパラメータの完全な事後分布を介して所望の統計値（例えば中央値、平均値、モーメント（moments）等．．．）を提供する。さらに、我々は、階層ベイジアンモデルを使用してモデルパラメータに関する事前知識を組み込む。結果として得られる事後分布は、経路における異常の機能的影響に関する重要な洞察を提供する。

本発明は、その発現が予測されている遺伝子からの祖先遺伝子の距離（すなわち、調節ネットワークにおける葉からルートへのリンクの数）に基づくペナリゼーションを伴うグローバル及びローカル収縮のアイディアを使用する。以下のモデルが構築される。なお、表記の便宜のために添え字gを省略している：

上記式は、リンク深さ情報をガンマ事前構築（gamma prior construction）に組み込むために、正常なガンマ事前構築を拡張する。この情報は、モデルパラメータの分散に含まれる係数kを介して利用される。したがって、βiの分散は、設定

を介して対応する祖先のリンク深さに反比例するように選択され、ここで、σ²はグローバル収縮を制御し、

は、ローカル収縮を示し、

は、リンク深さの影響を強化する。更なる柔軟性を提供するために、

についてガンマ事前分布の使用は、更なる柔軟性を提供する。ガンマを事前に使用すること（Using gamma prior）は、k_iについての閉形式（closed-form）事後分布を得るという利点を有し、よって、計算効率の良いGibbsサンプラーの利用を容易にする。したがって、分散の平均が深さパラメータに反比例するように

を使用する、すなわち、

である。定数ｃは、

を保証する正規化項であり、

を設定することにより得られる。したがって、我々は、k_i事前分布についての１つのフリーハイパーパラメータ

のみを有し、第２パラメータ

は、

から自動的に取得される。

であることに留意する。

を小さな値に設定することは、k_iについてより高い分散を提供し、したがって、あまり形成的ではないが、一方、大きな値の

は、ネットワークトポロジに関する高い確実性と、より短い経路を有するノード対が互いにより高い影響に関連付けられるという事実を反映する低い分散を提供する。この場合、ガンマ分布は、d_iの周囲に集中したガウス分布に近似する。我々は、基礎となる生物学的ネットワークの有意性を強調するために、

という比較的大きな値を選択する。

上記の階層モデルは、以下の完全同時分布をもたらす：

これは、結果として得られる確率が１に統合されることを保証するために、各パラメータについての完全条件付き事後分布が単に、その変数を含む項と正規化定数として機能する他の項との積であるという事実を利用して、以下の事後分布を直ちに提供する。この方法は項の補完（completion of terms）と呼ばれる：

ウッドベリー行列反転公式は、より安定した結果を得るためにn<pであるときにA^-1を計算するために使用され、p×pの正方行列反転をn×nに変換することによって計算を節約する。モデルパラメータβi、σの近似事後分布を得るために、バーンイン反復1000回と計算反復5000回でGibbsサンプラーを適用する。このプロセスは、すべてのサンプルs∈Sを使用してすべての遺伝子g∈Gに対して繰り返され、ここで、G及びSは、それぞれ遺伝子ID及びサンプルIDのセットである。

モジュール３：新たなサンプルの遺伝子レベル発現を予測し、すべての遺伝子の活性化及び一貫性レベルを報告
所与のサンプルについて標的遺伝子gの破壊を評価するために、活性化スコアA_g ^(new)及び不一致スコアC_g ^(new)を取得する。ここで最初のものは、その調節ネットワークと一貫性のある遺伝子発現レベルを示し、２つめは、（潜在的に体細胞突然変異に関連付けられる）遺伝子の調節解除（deregulation）を指示する遺伝子の期待値からの偏差を示す。

正常及び癌コホートの双方からのトレーニングサンプルを使用してモジュール２を実行することは、各機能が特定の遺伝子に対応する機能バンクの形で結果を提供する。この機能バンクを次いでモジュール３で使用して、テストサンプルを分析して潜在的な不一致を識別する。したがって、このモジュールは、すべての遺伝子について遺伝子発現レベル予測を行う。遺伝子ごとに、すべてのサンプルについて祖先遺伝子の発現レベル及び自己エピジェネティック情報を抽出する。次に、この遺伝子について学習された対応する機能を使用して、すべてのサンプルについてこの特定の遺伝子の発現レベルを予測する。予測プロセスは、この遺伝子の発現レベルについて条件付き事後分布を提供する。期待される遺伝子発現レベルを取得するために、最大事後確率（MAP：maximum-posteriori）法を使用する。

機能が学習された非単離標的遺伝子の一貫性スコアを計算するために、各々の新たなテストサンプルy^newについての任意の遺伝子のRNA発現の予測分布は、所与の入力x^new（自己エピジェネティック情報及び祖先の発現レベル）について条件付き事後分布からモデルパラメータを除外することにより得られることに留意する：

条件付き分布である第１項は閉形式で利用可能であるが、モデルパラメータの事後分布である第２項はそうではない。この分布を、以下の式で近似することができ、この場合、モデルパラメータ

の実現はギブスサンプリング法を用いて得られる。

上記の分布は、平均(Ψ(x^new)^Tβ⁽ⁱ⁾)と分散

の多数の等確率成分（equi-probable components）を有するガウス混合モデル（GMM：Gaussian mixture model）である。Gibbsサンプラーが収束する場合、β⁽ⁱ⁾は、共分散行列

を伴うβ_MAPの周囲に集中する。ここで、エンティティ

は、

と比べて小さい。したがって、Ψ(x^new)β⁽ⁱ⁾は、中心極限定理に従って、βi分布に関わらず、多数の予測因子について正規分布を近似する。計算と記憶を節約するために、予測分布の代用として次の正規分布を使用する：

ここで、

は、行列誘導ノルムである。この分布に基づいて、観察された値についてのz-スコア又は等価尤度（equivalent likelihood）を次のように計算する：

加えて、各遺伝子の基礎となる生物学的プロセスの複雑さと、未知の因子についての異なるレベルの遺伝のランダム性（inherit randomness）、自然の規則性（natural regularity）及び影響に起因して、学習した機能の予測力は、遺伝子ごとに大きく異なる可能性がある。したがって、正常サンプルについての各遺伝子の平均的な経験的予測性を、一貫性チェックのためのグラウンドレベルとみなす。よって、正常サンプルの平均的な不一致よりはるかに低い一貫性レベルを有する癌サンプルのみが、不一致サンプルとして報告される。以下の正規化された尤度が使用される：

ここで、n₀及びn₁は正常サンプル及び癌サンプルの数であり、αは、正常コホートと癌コホートに対して異なる強調をプッシュするための0と1の間の調整パラメータである。正常な癌をより強調し、より少ない数の正常なサンプルを補償するために、αについてより低い値を選択する。本発明では任意にα=1/10を設定する。これは、TCGA乳癌データセットのためのトレーニングセットにおける癌サンプルに対する正常サンプルの比にほぼ等しい。予測分布の分散がすべてのサンプルについて等しい場合、不等式（inequality）は等しいもの（equality）となる。上記のプロセスは、すべての遺伝子について並行して繰り返される。

一貫性スコアに加えて、各遺伝子の活性化スコアは、正規分布としてモデル化された遺伝子発現レベル分布を使用して得られる；

ここで、μ及びσは、外れ値を反復して除外した後の各遺伝子発現レベルについて学習された正規分布の平均及び標準偏差である。添え字gは、表記の便宜のために省略されている。同様の正規化が活性化スコアのために使用される。

上記で検討したように、このモジュールの用途は、トレーニングモデルを調節ネットワークの上部に対して使用して、標的遺伝子エピジェネティクス、並びに転写調節の役割を果たす遺伝子の発現レベルに基づいて、所与のサンプルについて所望の標的遺伝子発現レベルを予測することである。図５では、TCGA結腸癌データセットから導出される４２個の正常サンプルと４２個の腫瘍サンプルを含むテストサンプルにわたる遺伝子JUN発現レベルを予測するための実例が示されている。このモデルは、モジュール１及び２を使用して、５倍クロス検査（5-fold cross validation）を伴う３３８個の正常サンプル及び３６８個の癌サンプルを使用してトレーニングされる。遺伝子JUNは、モジュール１を使用して導出されるように、利用した経路ネットワークにおいてレベル２までの５１個の上流調節因子を有する。図５では、モジュール３内においてモジュール２で学習したモデルを用いて取得されるような正常サンプルと腫瘍サンプルの両方について、事後平均（posterior mean）の周囲の標準偏差とともに予測値が示されている。この図に示される信頼区間の提示は、予測値のみが得られ、予測の信頼性についての統計値が提供されないという点推定方法と比較して、遺伝子発現レベルを予測することにおける本発明の方法の利点である。第２の観察は、遺伝子JUNが正常サンプルにわたって厳密に調節されるということである。なぜなら、その調節因子の発現レベルに基づくその予測値は、癌サンプルと比較して、正常サンプルについてより正確であるからである。実際に、同様のレベルの偏差を有する１４個の腫瘍サンプルと比較して、予測値から３つの標準偏差を超えて逸脱するJUN発現レベルを経験するのは５つの正常サンプルのみである。

このモジュールで確立されるような遺伝子発現レベルの不一致と体細胞突然変異との間の関連を更に説明するために、図６は、調節ネットワークが利用可能であるすべての遺伝子にわたるBRCA及びCRCの両方についてのグローバル統計分析を提供する。これに関して、各遺伝子について、腫瘍サンプルを２つのサブセット、すなわち：ｉ）第１及び第２レベルの調整因子のうち関心のある遺伝子又はその一部についての遺伝子が突然変異しているもの；及びｉｉ）すべての調節因子が野生型（wild-type）であるもの、に分けられる。次に、突然変異サブセットと非突然変異サブセットの双方について絶対的な不一致レベルの平均をとる（図６Ａ、図６Ｃ）。２つのサブセットについての不一致スコアのヒストグラム（図６Ｂ及び図６Ｄ）は、双方の癌における突然変異サブセットの不一致スコアが非突然変異サブセットのものよりもかなり高いことを明らかにする。

図６Ａ及び図６Ｃにおいて、各ステムは特定の遺伝子に対応し、赤色のステムはその標的遺伝子又はその調節ネットワーク（レベル２まで）における突然変異のあるサンプルの平均絶対不一致（average absolute inconsistencies）であり、緑色のステムは、関心のある遺伝子及びその近い親（close parents）の遺伝子が野生型であるすべてのサンプルにわたる負の平均絶対一貫性スコア（negative of the average absolute consistency score）である。野生型調節遺伝子を有するサンプルについての緑色のステムは、提示の容易性のために垂直にひっくり返される。遺伝子は、野生型サンプルにおけるそれらの平均不一致レベルに基づいてソートされる。また、図６Ｂ及び図６Ｄは、平均不一致スコアについて得られるヒストグラムである。上段及び下段の列はそれぞれ乳癌及び結腸直腸癌についてである。結果は、調節ネットワーク内の標的遺伝子又はその近い親が体細胞突然変異を有するサンプル間にわたって、平均不一致のレベルがより高いことを示す。

モジュール４：体細胞突然変異と不一致との間の関連
遺伝子の発現レベルは、調節ネットワークにおける体細胞突然変異の存在によって予測値から逸脱する可能性があり、これにより調節機能の損失／増加をもたらすことがある。すなわち、調節遺伝子のいずれかにおける突然変異は、遺伝子発現の調節におけるその適切な役割に影響を与え、標的遺伝子発現に偏差を課すことがある。本方法のモジュール４は、下流標的遺伝子の不一致スコアに対する調節遺伝子の体細胞突然変異の影響を評価する方法論を提供する。したがって、このモジュールは、モジュール３によって提供される活性化及び一貫性スコアをとり、新たなテストサンプルごとに、有意に不一致な遺伝子を識別し、このような遺伝子が現在の遺伝子又は調節サブネットワーク内のCNV異常又は体細胞突然変異によって潜在的に引き起こされるかどうかを検査する。

まず、CNV異常事象によって引き起こされる不一致が識別される。不一致が遺伝子の過剰発現に起因し、遺伝子がコピー数重複（CNV>0.5）を経験する場合、CNV重複は、不一致の主原因として報告される。同様に、コピー数欠失（CNV<-0.5）が遺伝子の発現低下（down expression）と関連付けられる場合、CNV欠失は不一致のドライバー（driver）であると考えられる。

関連するコピー数異常を経験しない遺伝子については、下流遺伝子の転写に影響を与える遺伝子の上流調節ネットワークにおける突然変異から生じる可能性がある。調節遺伝子が下流の標的遺伝子に近づくほど、下流の遺伝子発現レベルの不一致に対するより大きな影響が期待される。したがって、モジュール４は、ルートノードgへのホップd_i,gを伴う突然変異した遺伝子iの影響が、値

でスケーリングされるように、グローバル深さペナリゼーションパラメータ0<α≦1を割り当てる。→1のように、深さの影響はあまり重要ではなくなる。α=1/2を結果セクションとして選択する。

調節ツリーにおける突然変異の影響を定量化するために、絶対不一致レベル及び深さペナリゼーション因子によってスケールされる癌サンプルの各々について、標的遺伝子又はその調節因子のいずれかに作用するすべての非サイレント突然変異をカウントする。一般に、遺伝子gの発現に対する遺伝子h突然変異の機能的影響はf_g(h)によって示され、以下のように計算される：

ここで、P_gは、遺伝子gの調節祖先遺伝子のセット（すなわち、対応する調節ツリーの葉）であり、M^(j)は、サンプルj内で突然変異した遺伝子のセットであり、

は、サンプルjにおける遺伝子gの不一致スコアであり、1(.)は指示関数である。分母の役割は、

を正規化することである。したがって、f_g(h)は、標的遺伝子gにおいて調節ネットワークに属するすべての遺伝子の突然変異の相対的影響を定量化する。

図７のフローチャートは、本方法におけるサンプルの不一致ごとの解釈を要約している。すべてのサンプルについてこの手順を繰り返し、遺伝子に割り当てられた体細胞突然変異の影響プロファイル（f_g(h)，∀g∈G，∀h∈P_g）に基づいて遺伝子をソートすることは、パッセンジャー事象（passenger events）をフィルタアウトし、その突然変異が下流の転写因子遺伝子に機能的に影響を与える最も影響力のある親遺伝子を決定する。したがって、本発明は、下流の遺伝子発現に影響を与える機能的突然変異の識別を可能にする。疾患の状況にわたる大部分の観察されたミスセンス突然変異の機能的影響がほとんど知られていないとすると、この進歩性は、臨床医及び／又は研究者が、所与の状況において最も可能性の高い機能的な疾患関連の突然変異に焦点を当てることを可能にし、したがって新規なバイオマーカー及び潜在的な治療目標の識別を可能にすることができる。

図８は、グラフ形式で示される、モジュール４で生成された結果の一例である。具体的には、図８Ａは、結腸癌で識別される遺伝子についてWnt経路標的遺伝子発現に対するAPC内の体細胞突然変異の相対的な影響を表示する。プロットされているのは、結腸癌サンプル内のAPCに影響を与える突然変異と標的遺伝子の活性化及び不一致の関連の有意性についての-log10（Pvalue）である。緑色で強調された遺伝子は有意に影響される（FDR≦15％）。図８Ｂでは、遺伝子発現の不一致に対する、PTENの上流の調節サブネットワークにおける体細胞突然変異の影響が表示されている。深さペナリゼーションパラメータはα=1/2に設定される。PTENの親における体細胞突然変異の組合せについて、その調節における変調効果（modulation effect）が示されている。この場合、遺伝子セット{PTEN、DYRK2、E4F1及びATF2}における突然変異は、PTENの発現低下との有意な関連を示す。したがって、これらの遺伝子は、PTENにおける体細胞突然変異の影響を変調する。したがって、DYRK2、E4F1及びATF2における突然変異は一緒にPTENの発現に影響を与え、したがって、これらの突然変異の組合せは、腫瘍におけるPTENのより正確な機能的状態を提供する。PTENの破壊がAKT経路の発癌性活性化をもたらすとすれば、これらの遺伝子の突然変異は治療の選択のための予後及び／又はバイオマーカーである。

実施例
本発明の方法の予測力を説明するために、その性能を、LASSO、RIDGE及び弾性ネット回帰を含む、いくつかの最適に近い最先端の点推定器（point-estimators）と比較する。

本発明の方法の正確さを実証するために、まず、有意な外れ値を反復して除外した後に、最尤法を介して各遺伝子発現レベルについてのガウス分布を学習する。各反復でサンプルについてガウス分布を学習することにより開始し、次いで、平均値の２番目の標準偏差近傍（second standard deviation neighborhood）にないサンプルを削除する。後続の反復では、アルゴリズムが収束し、更なる外れ値が存在しなくなるまで、残りのサンプルについて処理を繰り返す。サンプル遺伝子PTENについての経験的分布及び学習した正規分布が図９に提示される。比較目的のためにスチューデントTの分布も学習する。スチューデントTの分布は、外れ値に対するロバスト性の利点を有し、図９に示されるように外れ値の除外の後は正規分布に非常に近い。

次に、遺伝子発現レベルを、所定の閾値に基づいて３つの状態（ニュートラル、過剰発現及び過小発現）に分ける。閾値は、発現低下状態、ニュートラル状態及び過剰発現の確率がそれぞれ10%、80%及び10%になるように任意に設定される。モジュール３は、839個の非単離遺伝子のすべてについて患者固有の遺伝子発現予測を提供する。状態変化率は、すべての遺伝子及び患者にわたって状態変化事象を平均することを通して計算される。結果は、コホートごとに別々に計算される。サンプルi及び遺伝子gについての観察及び予測された発現状態がそれぞれ

である場合、状態変化率は次のように計算される：

テーブル１では、癌と高度に関連し、かつグローバル経路ネットワークにおいて上流調節遺伝子の有効なセットを有するいくつかの重要な遺伝子について予測誤差が計算される。本方法は、遺伝子発現レベルの完全な事後分布を提供するという追加の利点を有する最新のスパース性重視の回帰モデル（sparsity-imposing regression models）よりも優れていることが分かる。

別の重要な観察は、正常サンプルに対して癌サンプルの数が多いことに起因してモデルトレーニングに対する癌サンプルの寄与が高いという事実にもかかわらず、正常コホートがより良好な予測可能性を提示していることである。この観察はすべてのモデルに当てはまり、正常組織内の遺伝子発現の機能状態が上流調節ネットワークとより一致することを明らかにする。

癌サンプルと比べて正常サンプル内の標的遺伝子発現レベルの予測値と観察値との間の一貫性が高いという事実も図１０で観察される。図１０では、サンプル遺伝子MYB、GATA3、PTEN及びERBB2の観察値及び予測値が提示されている。ここで、正常サンプルにおける遺伝子発現レベルは、遺伝子自己エピジェネティックデータ及びその上流転写調節ネットワークから得られる予測とより一致する。この図は、遺伝子の発現レベルのみを分析する方法に関して、異なる源から生じ得る癌サンプルについての不一致分析の重要性を示し、経路攪乱（pathway perturbations）及び遺伝子調節不全に関する更なる情報を明らかにする。この不一致は、標的遺伝子におけるコピー数の重複及び欠失のような様々な源、並びに調節ネットワークの役割の正常な挙動を妨害し、結果として調節ネットワークのルートに存在する標的遺伝子の発現レベルに影響を与える調節ネットワーク内の突然変異に起因して生じる可能性がある。モデル係数に関して更なる洞察を得るために、２つの遺伝子ERBB2及びGATA3について得られたモデルパラメータをテーブル２及びテーブル３に提示する。各行は、異なる学習方法によって非線形ベイジアン法（within non-linear Bayesian method）について得られた、対応する係数値を提示する。事後分布についての標準偏差も最後の列の括弧内に提示されている。ERBB2の発現レベルは、提案した非線形ソフト閾値処理機能のモデルパラメータに見られるように、遺伝子座（locus）に影響を与えるコピー数の異常事象に大きく依存することが示される。この非線形性は、測定ノイズの可能性があるゼロ付近の小さな乱れに対するモデルの無知（ignorance）を反映する。したがって、SNRアレイから導出されるコピー数関連の対数比（logRatio）の値は、対数比を重複／ニュートラル／欠失状態へ離散化する必要なく、モデルで直接使用することができる。非線形関数の関連性は、興味深いことにすべての学習方法によってピックアップされる。図１１はこの関連性を検証し、観察されたRNAと予測されたRNA対CNVとの間の関係を遺伝子ERBB2について示している。図１１では、青色及び赤色の点は、モデルから得られる観察値及び予測値に対応する。黒い曲線は、テーブル２のモデルパラメータによって得られる非線形のRNA CNV関係である。

この図は、学習プロセスから得られる係数を有する非線形CNV項が、ERBB2についてのRNA発現レベルを、DNAメチル化及び祖先遺伝子発現レベルのような他の項に起因して多少の変動を伴って定義することを示している。実際、DNAメチル化の係数及び大部分の祖先は、LASSO及び弾性ネット法によって予測因子リストから明示的に除去され、また、本発明ではDNAメチル化について無視できる係数を割り当てるということも注目に値する。

一方、GATA3についてのRNA発現レベルは、DNAメチル化並びに上流調節ネットワークによってより影響を受ける。DNAメチル化係数の期待される負の符号は、両方の遺伝子について遺伝子発現レベルとDNAメチル化との間の逆の関係を示唆する。最後に、GATA3については、上流調節ネットワークが、この遺伝子の発現を調節するのに極めて重要な役割を果たし、乳癌におけるこの遺伝子の発現の変異の大部分が主として転写因子の活性に起因することを示唆している。テーブル２及びテーブル３に提供される２つの遺伝子ERBB2とGATA3のための方法によって推定される回帰係数は、遺伝子調節機能の高い異種性に起因して、回帰係数が遺伝子に関して有意に異なる可能性があることを明らかにする。

不一致の重要な原因は、標的遺伝子の上流調節ネットワークにおける突然変異に起因する。調節遺伝子の発現レベルの影響が既にこの方法によって捕捉されるという事実に留意すると、標的遺伝子の発現レベルの予測値が観察値と一致しない場合には、調節ネットワークがその調節役割を適切に果たしていないと推測する。この調節ネットワークの機能不全のほとんどは、調節ネットワークにおける体細胞突然変異から生じ、そのような機能不全は、その遺伝子又は産物のタンパク質が、その機能（複合体形成、遺伝子転写、タンパク質活性化等）を適切に行うことを妨げ、このことは下流標的遺伝子発現レベルに影響を与えることになる。説明的な例として、遺伝子PTENの調節不全に対する体細胞突然変異の機能的影響が図１２に示されており、PTEN発現における不一致が、TP53、PTEN、PIK3CA、MAP3K1及びMAP2K4の突然変異に高く関連付けられることを明らかにしている。PIK3CAがTP53より頻繁に（それぞれ387サンプル対333サンプル）突然変異されるとすると、TP53の突然変異対PIK3CAのより高い影響は、特に興味深い。我々は、luminalの乳癌に関連付けられることが以前に示されたMAP3K1及びMAP2K4の突然変異が、PTEN不活性化に影響を与え、したがって、乳癌の主要なサブタイプを引き起こす際のこれらの遺伝子間の興味深いつながりを提供することを観察する。我々はまた、PTENについての不一致スコアに対するタンパク質切断型（protein-truncating）及び他の非同義型の突然変異の相対的な影響を計算する。PTENにおけるタンパク質切断の突然変異は、ナンセンス介在性のPTEN mRNAの崩壊と一貫して、その調節解除に対してより大きな影響を有するが、このモデルは、PTENの調節遺伝子のいずれかに影響するときに、２種類の突然変異が同様の影響を有すると判断する。深さペナリゼーションパラメータはα=1/2に設定される。

Claims

コンピュータシステムによって、調節不全遺伝子を引き起こす患者固有の体細胞異常を識別する方法であって：
前記コンピュータシステムのプロセッサが、生物学的ネットワーク経路情報を取得することにより、各標的遺伝子について上流調節親遺伝子情報の第１データセットを決定するステップと；
前記プロセッサが、前記標的遺伝子の各々について、前記第１データセットから調節サブネットワークを決定するステップであって、前記調節サブネットワークは、前記標的遺伝子の発現レベルと前記標的遺伝子のゲノム及びエピジェネティック状態、並びにその遺伝子の上流転写調節因子の関係を備える、ステップと；
前記プロセッサが、ＲＮＡｓｅｑ発現データ、コピー数多型データ及びＤＮＡメチル化データのうちの少なくとも１つを含む、測定ベースのオミックスデータの第２データセットを決定するステップと；
前記プロセッサが、前記第１データセットと前記第２データセットを統合するステップと；
前記プロセッサが、前記統合された第１及び第２データセットから、前記標的遺伝子の各々について非線形関数を生成するステップであって、前記非線形関数は、前記調節サブネットワークに関連付けられる測定値に、前記遺伝子の発現レベルを関連付ける、ステップと；
前記プロセッサが、前記標的遺伝子の各々について、該標的遺伝子についての前記非線形関数を使用して、期待される発現レベルを計算するステップと；
前記プロセッサが、前記標的遺伝子について観察された遺伝子発現レベルに関連し、かつ前記標的遺伝子のうちの１つ以上についての前記決定された調節サブネットワーク内の１つ以上の親遺伝子のシーケンスを備える、患者固有の情報の第３データセットを決定するステップと；
前記プロセッサが、前記標的遺伝子の各々について、期待される遺伝子発現レベルと観察された患者固有の発現レベルとの間の患者固有の不一致スコアを計算するステップと；
前記プロセッサが、前記標的遺伝子の各々について、患者固有の活性化スコアを計算するステップであって、前記活性化スコアは、正規分布としてモデル化された遺伝子発現レベル分布を使用して得られる、ステップと；
前記プロセッサが、発現レベルが前記期待される発現レベルと有意に不一致な患者固有の標的遺伝子を識別するために、すべての患者サンプルについての活性化及び不一致スコアを評価するステップと；
前記プロセッサが、標的遺伝子の発現レベルにおける不一致と、その特定の標的遺伝子の上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連を識別するステップであって、突然変異が識別された特定の標的遺伝子の上流調節ネットワークにおける各親遺伝子について、前記計算された患者固有の不一致スコア、前記特定の標的遺伝子の上流調節ネットワーク内の遺伝子及び１つ以上の突然変異を含む１組の遺伝子に少なくとも部分的に基づいて、体細胞突然変異の機能的影響スコアを計算するステップを含む、ステップと；
前記プロセッサが、前記特定の標的遺伝子の前記上流調節ネットワーク内の２つ以上の親遺伝子についての前記計算された機能的影響スコアに基づいて、最も影響力のある親遺伝子を決定するステップであって、該最も影響力のある親遺伝子は、前記特定の標的遺伝子の前記上流調節ネットワーク内の他の親遺伝子と比較して、前記標的遺伝子の発現に影響を与えた可能性が最も高い、突然変異した親遺伝子を備える、ステップと；
前記プロセッサが、前記の有意な不一致を有する標的遺伝子を異常又は調節不全遺伝子として報告するステップであって、該報告は、前記の有意な不一致を有する標的遺伝子のうちの１つ以上について最も影響力のある標的遺伝子の識別を含む、ステップと；
を含む、方法。
前記非線形関数は、前記測定ベースのオミックスデータから得られる前記遺伝子のエピジェネティック情報及び前記遺伝子の調節サブネットワークの状態に基づいて決定される、
請求項１に記載の方法。
前記非線形関数は、サブネットワーク内の調節遺伝子の潜在的により強い影響を捕捉するグローバル深さペナルティメカニズムを用いて決定される、
請求項２に記載の方法。
前記患者固有の情報は、ＲＮＡ発現データ、ＣＮＶデータ、メチル化データ及び体細胞突然変異データのような癌サンプルデータを含む、
請求項１に記載の方法。
個々の患者サンプルにおける遺伝子発現レベルの有意な偏差及び不一致を識別するための統合された統一ネットワークを備えるシステムであって：
キュレートされた生物学的ネットワーク経路情報から取得される各標的遺伝子についての上流調節親遺伝子情報の第１データセットであって、前記経路情報を受け取るように構成されたプロセッサ上に配置され、前記標的遺伝子の発現レベルと前記標的遺伝子のゲノム及びエピジェネティック状態、並びにその遺伝子の上流転写調節因子の関係を備える前記第１データセットと；
前記標的遺伝子の発現レベルと前記標的遺伝子のゲノム及びエピジェネティック状態、並びにその上流転写調節因子の関係を捕捉する各固有の標的遺伝子についての調節ツリーであって、前記第１データセットから決定される調節ツリーと；
ＲＮＡｓｅｑ発現データ、コピー数多型データ及びＤＮＡメチル化データのうちの少なくとも１つを含む、測定ベースのオミックスデータの第２データセットであって、そのようなデータを受け取るように構成されるプロセッサ上に配置される第２データセットと；
各標的遺伝子についての非線形関数であって、該非線形関数のパラメータは、修正ベイジアン推論法を使用して決定される、非線形関数と；
前記標的遺伝子についての観察された発現レベルに関連し、かつ前記標的遺伝子のうちの１つ以上についての前記決定された調節サブネットワーク内の１つ以上の親遺伝子のシーケンスを備える、患者固有の情報の第３データセットであって、前記患者固有の情報は、新たな癌サンプルデータを含む、患者固有の情報の第３データセットと；
を備え、
前記標的遺伝子の発現レベルは前記非線形関数を利用して決定され、関連する患者固有の不一致スコアは、所与のサンプル内の前記標的遺伝子について予測された発現レベルと観察された発現レベルとの間で決定され、
活性化及び不一致スコアがすべてのテストサンプルについて決定され、これにより、前記標的遺伝子の発現レベルの不一致とその特定の遺伝子の上流調節ネットワークにおける体細胞突然変異との間の統計的に有意な関連は、（ｉ）突然変異が識別された特定の標的遺伝子の上流調節ネットワークにおける各親遺伝子について、計算された患者固有の不一致スコア、前記特定の標的遺伝子の上流調節ネットワーク内の遺伝子及び１つ以上の突然変異を含む１組の遺伝子に少なくとも部分的に基づいて、体細胞突然変異の機能的影響スコアを計算するステップと；（ｉｉ）前記特定の標的遺伝子の前記上流調節ネットワーク内の２つ以上の親遺伝子についての前記計算された機能的影響スコアに基づいて、最も影響力のある親遺伝子を決定するステップであって、該最も影響力のある親遺伝子は、前記特定の標的遺伝子の前記上流調節ネットワーク内の他の親遺伝子と比較して、前記標的遺伝子の発現に影響を与えた可能性が最も高い、突然変異した親遺伝子を備える、ステップと；を備えるプロセスによって、識別され、前記活性化スコアは、正規分布としてモデル化された遺伝子発現レベル分布を使用して得られる、システム。
前記非線形関数は、前記測定ベースのオミックスデータから得られる前記遺伝子のエピジェネティック情報及び前記遺伝子の調節サブネットワークの状態に基づいて決定される、
請求項５に記載のシステム。
前記修正ベイジアン推論法によって決定される前記非線形関数は、サブネットワーク内の調節遺伝子の潜在的により強い影響を捕捉するグローバル深さペナルティメカニズムを組み込む、
請求項６に記載のシステム。
前記患者固有の情報は、ＲＮＡ発現データ、ＣＮＶデータ、メチル化データ及び体細胞突然変異データのような癌サンプルデータを含む、
請求項５に記載のシステム。