(57)【要約】【課題】  遺伝子発現の時系列研究から遺伝子間のネットワーク関係を推定するためのノンパラメトリック回帰によるベイジアンモデルを使用する技術を提供する。【解決手段】  本発明の実施形態は、遺伝子ネットワーク含む複雑な生物学的情報を分析するための新しい推定法の応用を含む。いくつかの実施形態では、生物体におけるいくつかの遺伝子に対して経時変化データが同時に得られる。新しい方法は、ベイジアン推定法の修正と、発現遺伝子間の原因及び結果の関係を判断し、いくつかの実施形態では、調節遺伝子の上流エフェクターを判断するためのこれらの方法の応用とを含む。ベイジアン法の付加的な修正には、経時変化データを用いて発現遺伝子間の因果関係を推定することが含まれる。他の実施形態は、ブートストラップ法の使用と、発現遺伝子間のより正確なネットワーク情報をもたらすためのエッジ効果の判断とを含む。遺伝子ネットワークに関する情報は、メモリ装置に保存することができ、かつ出力装置に送信することができ、又は遠隔位置に送信することができる。
2006506746.tif 000002
【特許請求の範囲】
【請求項１】
遺伝子ネットワークを構築する方法であって、
（ａ）生物体の一組の遺伝子に対して、該遺伝子の組の各遺伝子の発現の経時変化に基づく発現結果を含み、該遺伝子の相互の平均的影響及び各時間点の変動性の尺度を定量化する定量的経時変化データライブラリを準備する段階と、
（ｂ）前記ライブラリから遺伝子発現マトリックスを作成する段階と、
（ｃ）前記遺伝子間のネットワーク関係を生成する段階と、
（ｄ）遺伝子の１つ又はそれ以上の群が該遺伝子の群の他のものと異なって発現するか判断する段階と、
を含むことを特徴とする方法。
【請求項２】
（ｅ）ＢＮＲＣdynamic判断基準を最小にする段階を含むベイジアン計算モデルを準備する段階、
を更に含むことを特徴とする請求項１に記載の方法。
【請求項３】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、多項式基底、フーリエ級数、ウェーブレット基底、回帰スプライン基底、及びＢ−スプラインから成る群から選択された非線形曲線当て嵌め法を用いる段階を含むことを特徴とする請求項２に記載の方法。
【請求項４】
前記データライブラリは、遺伝子発現を変化させる経時変化研究を用いて作成されることを特徴とする請求項１に記載の方法。
【請求項５】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、逆当て嵌めアルゴリズムを用いてベイジアンモードを選択する段階を更に含むことを特徴とする請求項２に記載の方法。
【請求項６】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、Ａｋａｉｋｅの情報判断基準を用いる段階を更に含むことを特徴とする請求項２に記載の方法。
【請求項７】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、最大尤度推定を用いる段階を更に含むことを特徴とする請求項２に記載の方法。
【請求項８】
前記遺伝子は、細胞サイクルに関連することを特徴とする請求項１に記載の方法。
【請求項９】
前記変動性の尺度は、分散であることを特徴とする請求項２に記載の方法。
【請求項１０】
前記非線形曲線当て嵌め法は、ノンパラメトリック法であることを特徴とする請求項３に記載の方法。
【請求項１１】
ＢＮＲＣdynamic判断基準を最小にするための前記ノンパラメトリック法は、不均一誤差分散を用いる段階を含むことを特徴とする請求項１０に記載の方法。
【請求項１２】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、
（１）（ｉ、ｊ）番目の要素がグラフ遺伝子i→遺伝子jのＢＮＲＣjdynamicスコアであるスコアマトリックスを作る段階と、
（２）最小ＢＮＲＣdynamicをもたらす追加、除去、及び逆転の１つ又はそれ以上を実施する段階と、
（３）前記ＢＮＲＣdynamicがそれ以上減少しなくなるまで段階２を繰り返す段階と、
を更に含むことを特徴とする請求項１１に記載の方法。
【請求項１３】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、ＢＮＲＣ(j)dynamicを最小にするためにヒルクライムアルゴリズムを適用する段階を更に含むことを特徴とする請求項１１に記載の方法。
【請求項１４】
エッジの強度は、ブートストラップ法を用いて判断されることを特徴とする請求項１１に記載の方法。
【請求項１５】
前記ブートストラップ法は、
（１）オリジナル遺伝子ライブラリ発現データから、置換を用いて何回か無作為にサンプリングすることにより、ブートストラップ遺伝子発現マトリックスを準備する段階と、
（２）遺伝子i及び遺伝子jに対する遺伝子ネットワークを推定する段階と、
（３）段階（１）及び（２）をＴ回繰り返し、それによってＴ個の遺伝子ネットワークを生成する段階と、
（４）遺伝子iと遺伝子jの間のブートストラップエッジ強度を（ｔ1＋ｔ2）／Ｔとして計算する段階と、
を含むことを特徴とする請求項１４に記載の方法。
【請求項１６】
遺伝子ネットワークを計算する方法であって、
（ａ）生物体の複数の遺伝子に対する経時変化遺伝子発現データの生データライブラリを準備する段階と、
（ｂ）前記生データライブラリから背景信号強度を差し引く段階と、
（ｃ）前記複数の遺伝子の各々に対して遺伝子発現の相対変化を計算する段階と、
（ｄ）スチューデントのｔ検定を用いて、遺伝子発現の前記相対変化の統計的有意性を分析する段階と、
（ｅ）遺伝子発現の前記変化を線形スプライン関数に当て嵌める段階と、
を含むことを特徴とする方法。
【請求項１７】
発現レベルが主としてノイズにより決定されるほど十分に低い遺伝子を考慮から除外する段階を更に含むことを特徴とする請求項１６に記載の方法。
【請求項１８】
前記グループ分けする段階は、前記遺伝子を１つ又はそれ以上の同等の組にグループ分けする段階を含むことを特徴とする請求項１に記載の方法。
【請求項１９】
遺伝子ネットワーク関係を推定する方法であって、
（１）（ｉ、ｊ）番目の要素がグラフ遺伝子i→遺伝子jのＢＮＲＣスコアであるｐ×ｐマトリックスを作り、
（２）小さなＢＮＲＣスコアを与える遺伝子jの親遺伝子の候補の組を選択し、
（３）前記親遺伝子の計算順序を選択し、
（４）以下の段階を繰り返し、
（４．１）遺伝子jに対して、親遺伝子を加えるか又は親遺伝子を削除する段階、
（４．２）ＢＮＲＣdynamicスコアを再計算する段階、及び
（４．３）適切な収束判断基準が満たされるまで段階４．１及び４．２を繰り返す段階、
（５）段階（３）の前記親遺伝子の計算順序を入れ変え、
（６）段階（４）を繰り返し、
（７）ＢＮＲＣdynamicが最小になるまで段階（５）及び（６）を繰り返す、
段階を含むことを特徴とする方法。
【請求項２０】
経時変化遺伝子発現データから遺伝子のネットワークを包含するシステムの遺伝子ネットワークモデルを構築する方法であって、
ＢＮＲＣdynamic判断基準を最小にする段階を含むベイジアン計算モデルを用いる段階、
を含むことを特徴とする方法。
【請求項２１】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、多項式基底、フーリエ級数、ウェーブレット基底、回帰スプライン基底、及びＢ−スプラインから成る群から選択された非線形曲線当て嵌め法を用いる段階を含むことを特徴とする請求項２０に記載の方法。
【請求項２２】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、逆当て嵌めアルゴリズムを用いてベイジアンモードを選択する段階を含むことを特徴とする請求項２０に記載の方法。
【請求項２３】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、Ａｋａｉｋｅの情報判断基準を用いる段階を含むことを特徴とする請求項２０に記載の方法。
【請求項２４】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、最大尤度推定を用いる段階を含むことを特徴とする請求項２０に記載の方法。
【請求項２５】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、ノンパラメトリック法である非線形曲線当て嵌め法を用いる段階を含むことを特徴とする請求項２０に記載の方法。
【請求項２６】
前記ノンパラメトリック法は、不均一誤差分散を用いる段階を含むことを特徴とする請求項２５に記載の方法。
【請求項２７】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、
（１）（ｉ、ｊ）番目の要素がグラフ遺伝子i→遺伝子jのＢＮＲＣjdynamicスコアであるスコアマトリックスを作る段階と、
（２）最小ＢＮＲＣdynamicをもたらす追加、除去、及び逆転の１つ又はそれ以上を実施する段階と、
（３）前記ＢＮＲＣdynamicがそれ以上減少しなくなるまで段階２を繰り返す段階と、
を更に含むことを特徴とする請求項２６に記載の方法。
【請求項２８】
前記ＢＮＲＣdynamic判断基準を最小にする段階は、ＢＮＲＣ(j)dynamicを最小にするためにヒルクライムアルゴリズムを適用する段階を更に含むことを特徴とする請求項２６に記載の方法。
【請求項２９】
エッジの強度は、ブートストラップ法を用いて判断されることを特徴とする請求項２６に記載の方法。
【請求項３０】
前記ブートストラップ法は、
（１）オリジナル遺伝子ライブラリ発現データから置換を用いて何回か無作為にサンプリングすることにより、ブートストラップ遺伝子発現マトリックスを準備する段階と、
（２）遺伝子i及び遺伝子jに対する遺伝子ネットワークを推定する段階と、
（３）段階（１）及び（２）をＴ回繰り返し、それによってＴ個の遺伝子ネットワークを生成する段階と、
（４）遺伝子iと遺伝子jの間のブートストラップエッジ強度を（ｔ1＋ｔ2）／Ｔとして計算する段階と、
を含むことを特徴とする請求項２９に記載の方法。
【請求項３１】
請求項２０に記載の方法により構築された遺伝子ネットワークモデルを含むデータファイル。
【請求項３２】
コンピュータ可読形態であることを特徴とする請求項３１に記載のデータファイル。
【請求項３３】
遠隔位置からアクセス可能であることを特徴とする請求項３１に記載のデータファイル。
【請求項３４】
インターネットウェブの位置からアクセス可能であることを特徴とする請求項３１に記載のデータファイル。
【請求項３５】
遺伝子ネットワークを包含するシステムにおいてターゲット遺伝子を識別する方法であって、
（ａ）ベイジアン計算モデルを用いて第１及び第２の遺伝子ネットワークモデルを構築する段階、
を含み、
前記ベイジアン計算モデルは、ＢＮＲＣdynamic判断基準を最小にする段階を含み、前記第１の遺伝子ネットワークモデルは、第１の遺伝子発現プロフィールを分析することによって得られ、前記第２の遺伝子ネットワークモデルは、第２の遺伝子発現プロフィールを分析することによって得られるものであり、
前記第１の遺伝子発現プロフィールは、第１の時間点でシステムから得られ、前記第２の遺伝子発現プロフィールは、該第１の時間点の後の第２の時間点で該システムから得られるものであり、
（ｂ）前記ベイジアン計算モデルを用いて前記第１及び第２の遺伝子ネットワークモデルを分析する段階、
を更に含み、
遺伝子発現の経時変化は定量化され、親遺伝子は、ターゲット遺伝子として識別される、
ことを特徴とする方法。
【請求項３６】
前記ターゲット遺伝子は、親遺伝子であることを特徴とする請求項３５に記載の方法。
【請求項３７】
前記ターゲット遺伝子は、親遺伝子から下流の遺伝子であることを特徴とする請求項３５に記載の方法。
【請求項３８】
請求項３５に記載の方法によって得られる１つ又はそれ以上のターゲット遺伝子のアイデンティティを包含するデータファイル。
【請求項３９】
コンピュータ可読形態であることを特徴とする請求項３８に記載のデータファイル。
【請求項４０】
遠隔位置からアクセス可能であることを特徴とする請求項３８に記載のデータファイル。
【請求項４１】
インターネットウェブの位置からアクセス可能であることを特徴とする請求項３８に記載のデータファイル。
【請求項４２】
サービスを提供する方法であって、
（１）一群の遺伝子の経時変化発現データを含むデータセットをパーティから受け取る段階と、
（２）ＢＮＲＣdynamic判断基準を最小にすることにより、前記群の遺伝子間のネットワーク関係を判断する段階と、
を含むことを特徴とする方法。
【請求項４３】
前記データセットを受け取る段階は、前記遺伝子の少なくとも１つのアイデンティティを受け取る段階を含むことを特徴とする請求項４２に記載の方法。
【請求項４４】
サービスを提供する方法であって、
パーティからエージェントを受け取る段階と、
請求項３５に記載の方法により構築した遺伝子ネットワークモデルを用いて、前記パーティに対するターゲット遺伝子を識別する段階と、
を含むことを特徴とする方法。
【発明の詳細な説明】
【技術分野】
【０００１】
関連出願
本出願は、「３５Ｕ．Ｓ．Ｃ．§１１９（ｅ）」の下で、２００２年１１月１９日出願の米国特許仮出願出願番号第６０／４２７，４４８号に対する優先権を請求するものである。この出願は、本明細書においてその全内容が引用により組み込まれている。
本発明は、遺伝子発現の時系列研究から遺伝子間のネットワーク関係を推定するためのノンパラメトリック回帰によるベイジアンモデルの使用に関するものである。特に、本発明は、最適ネットワーク関係を推定するために判断基準ＢＮＲＣdynamicを最小にする段階を伴う方法に関する。
【背景技術】
【０００２】
生命科学、医学、薬剤の発見及び開発、及び薬品工業における現在の研究及び開発の最も重要な態様の１つは、大量の生データを解釈し、このようなデータに基づいて結論を導出する方法及び装置を開発する必要性である。生命情報科学は、システム生物学の理解に大きく貢献しており、生体系の構成要素間の複雑な関係の理解を更に深めるのに有望である。特に、発現遺伝子をすばやく検出し、遺伝子の発現を定量化するための新しい方法の出現と共に、生命情報科学は、特定遺伝子が有機体の生物学で果たす正確な役割を確実に知ることなく潜在的な治療ターゲットを予測するのに用いることができる。
【０００３】
遺伝子系のシミュレーションは、システム生物学の中心的問題である。シミュレーションは、生物学的知識に基づくことができるために、ネットワーク推定法は、予め分からない関係を予測又は推定することにより生物学的シミュレーションをサポートすることができる。
特に、マイクロアレイ技術の開発により、様々な生物体からの多くの遺伝子の発現の研究が可能になった。生物体からのいくつかの遺伝子から非常に大量の生データを得ることができ、遺伝子発現は、突然変異、疾病、又は薬物の何れかで介入することによって研究することができる。特定の遺伝子の発現が特定の疾病で又は特定の介入に応答して増大することを見出すことにより、その遺伝子が疾病過程又は薬物反応に直接関わっていると考えることができる。しかし、生物学的な生物体では、多くの遺伝子が特定の介入により影響を及ぼされる可能性があるという点で、遺伝子がこのような介入のいずれによっても独立に調節されることはほとんどない。多くの異なる遺伝子がこのように影響される場合があるために、このような研究では、遺伝子間の原因及び結果の関係を理解することは非常に困難である。従って、遺伝子が生物学的現象の中心であり、遺伝子の発現が研究中の生物学的過程に対して末梢的であるところの遺伝子間の原因及び結果の関係を判断する方法の開発に多大な努力が為されている。このような末梢的遺伝子の発現は、生物学的又は病態生理学的状態のマーカとして有用である場合があるが、このような遺伝子が生理学的又は病態生理学的状態での中心的存在でない場合には、このような遺伝子に基づいて薬物を開発するように努めても報われないであろう。逆に、過程に中心的であると識別された遺伝子に対しては、薬物又は他の介入の開発は、遺伝子の発現の変化に付随する状態に対する治療の開発に決定的である場合がある。
【０００４】
マイクロアレイ遺伝子発現データから遺伝子ネットワークを推定するためのベイジアンネットワーク解析の開発は相当な注目を集めており、多くの成功した研究が報告されている（Ｆｒｉｅｄｍａｎ他［１３］、Ｉｍｏｔｏ他［１４］、Ｐｅ’ｅｒ他［１８］、及び本発明人の研究［本明細書においてその全内容が引用により組み込まれている米国特許出願出願番号第１０／２５９，７２３号］）。
しかし、従来のベイジアンネットワークモデルの欠点は、循環ネットワークを構築することができず、一方、いくつかの実際の遺伝子調節機構は循環構成要素を有することである。近年、循環調節構成要素を有する遺伝子ネットワークを構築するために、動的ベイジアンネットワークモデル（Ｂｉｌｍｅｓ他［３］、Ｆｒｉｅｄｍａｎ他［１２］、Ｓｏｍｅｒｅｎ他［１９］）が提案されている。動的ベイジアンネットワークは、時系列データに基づいており、通常は、データがいくつかの部類に区分化することができるものである。従って、動的ネットワークモデルは、区分化過程の閾値の設定に依存する可能性があり、残念なことに、区分化により情報が失われる可能性がある。Ｉｍｏｔｏ他［１４、１５］は、区分化を避けるためのソリューションとして、かつ遺伝子間の非線形関係を捕捉するためにベイジアンネットワーク及びノンパラメトリック回帰に基づくネットワーク推定法を提案している。しかし、ベイジアンネットワーク及びノンパラメトリック回帰モデル［１４、１５］は、依然として循環調節構成要素を有するネットワークを適切に解くことはできないであろう。
【０００５】
【特許文献１】米国特許仮出願出願番号第６０／４２７，４４８号
【特許文献２】米国特許出願出願番号第１０／２５９，７２３号
【非特許文献１】Ａｋａｉｋｅ、Ｊ．：「情報理論と最大尤度原理の拡張」、掲載：Ｐｅｔｒｏｖ、Ｂ．Ｎ．、Ｃｓａｋｉ、Ｆ．（編）：第２回情報理論国際シンポジウム、Ａｋａｄｅｍｉａｉ  Ｋｉａｄｏ、ブダペスト、２６７〜２８１頁（１９７３）
【非特許文献２】Ｂｅｒｇｅｒ、Ｊ．Ｏ．：「統計的確定理論とベイジアン分析」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ニューヨーク（１９８５）
【非特許文献３】Ｂｉｌｍｅｓ、Ｊ．Ａ．：「動的ベイジアンマルチネット」、第１６回人工知能の不確定性に関する会議紀要、３８〜４５頁（２０００）
【非特許文献４】Ｂｕｒｎｈａｍ、Ｋ．Ｐ．、Ａｎｄｅｒｓｏｎ、Ｄ．Ｒ．：「モデル選択と推定、実用情報−理論的手法」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ニューヨーク（１９８８）
【非特許文献５】Ｃｈｅｎ、ＴＩ．、Ｈｅ、Ｈ．Ｌ．、Ｃｈｕｒｃｈ、Ｇ．Ｍ．：「微分方程式を用いる遺伝子発現のモデリング」、バイオコンピューティングに関する太平洋シンポジウム紀要、第４巻、２９〜４０（１９９９）
【非特許文献６】Ｄａｖｉｓｏｎ，Ａ．Ｃ．：「近似予測尤度」、Ｂｉｏｍｅｔｒｉｋａ、第７３巻、３２３〜３３２（１９８６）
【非特許文献７】ＤｅＢｏｏｒ，Ｃ．：「スプラインの実用ガイド」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ベルリン（１９７８）
【非特許文献８】Ｄｅ  Ｈｏｏｎ，Ｍ．Ｊ．Ｌ．、Ｉｍｏｔｏ，Ｓ．、Ｋｏｂａｙａｓｈｉ，Ｋ．、Ｏｇａｓａｗａｒａ，ＮＨ．、Ｍｉｙａｎｏ，Ｓ．：「微分方程式を用いる「Ｂａｃｉｌｌｕｓ  ｓｕｂｔｉｌｉｓ」の時間順遺伝子発現データからの遺伝子調節ネットワークの推定」、バイオコンピューティングに関する太平洋シンポジウム紀要、第８巻、２００３、近刊
【非特許文献９】ＤｅＲｉｓｉ，Ｊ．、Ｌｙｅｒ、Ｖ．Ｒ．、Ｂｒｏｗｎ、Ｐ．Ｏ．：「ゲノム規模の遺伝子発現の代謝及び遺伝子制御の探求」、Ｓｃｉｅｎｃｅ、第２７８巻、６８０〜６８６頁（１９９７）
【非特許文献１０】Ｄｉｅｒｃｋｘ、Ｐ．：「スプラインを用いる曲線及び表面当て嵌め」、Ｏｘｆｏｒｄ（１９９３）
【非特許文献１１】Ｅｉｌｅｒ、Ｐ．Ｈ．Ｃ．、Ｍａｉｘ、Ｂ．：「Ｂスプラインを用いる柔軟な平滑化及び不利益（討論）」、Ｓｔａｔｉｓｔｉｃａｌ  Ｓｃｉｅｎｃｅ、第１１巻、８９〜１２１頁（１９９６）
【非特許文献１２】Ｆｒｉｅｄｍａｎ，Ｎ．、Ｍｕｒｐｈｙ，Ｋ．、Ｒｕｓｓｅｌｌ，Ｓ．：「動的確率ネットワークの構造の学習」、人工知能における不確定性に関する会議紀要、１３９〜１４７頁（１９９８）
【非特許文献１３】Ｆｉｒｅｄｍａｎ，Ｎ．、Ｌｉｎｉａｌ，Ｍｌ、Ｎａｃｈｍａｎ，Ｉ．、Ｐｅ’ｅｒ，Ｄ．：「発現データの分析のためのベイジアンネットワークの使用」、Ｊ．Ｃｏｍｐ．Ｂｉｏｌ．、第７巻、６０１〜６２０（２０００）
【非特許文献１４】Ｉｍｏｔｏ，Ｓ．、Ｇｏｔｏ，Ｔ．、Ｍｉｙａｎｏ，Ｓｌ．：「ベイジアンネットワーク及びノンパラメトリック回帰を用いる遺伝子間の遺伝子ネットワーク及び機能的構造の推定」、バイオコンピューティングに関する太平洋シンポジウム紀要、第７巻、１７５〜１８６（２００２）
【非特許文献１５】Ｉｍｏｔｏ，Ｓ．、Ｋｉｍ，Ｓ．、Ｇｏｔｏ，Ｔ．、Ａｂｕｒａｔａｎｉ，Ｓ．、Ｔａｓｈｉｒｏ，Ｋ．、Ｋｕｈａｒａ，Ｓ．、Ｍｊｉｙａｎｏ，Ｓ．：「遺伝子ネットワークの非線形モデリングのためのベイジアンネットワーク及びノンパラメトリック異分散回帰」、ＩＥＥＥ  Ｃｏｍｐｕｔｅｒ  Ｓｏｃｉｅｔｙ  Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ  Ｃｏｎｆｅｒｅｎｃｅ紀要、２１９〜２２７頁（２００２）
【非特許文献１６】Ｋｏｎｉｓｈｉ，Ｓ．：「統計学的モデルの評価及び情報判断基準」、掲載：Ｇｈｏｓｈ，Ｓ．（編）「多変量分析、実験計画及び調査サンプリング」、Ｍａｒｃｅｌ  Ｄｅｋｋｅｒ、ニューヨーク、３６９〜３９９頁（１９９９）
【非特許文献１７】Ｋｏｎｉｓｈｉ，Ｓ．、Ｋｉｔａｇａｗａ，Ｇ．：「モデル選択における汎用情報判断基準」、Ｂｉｏｍｅｔｒｉｋａ、第８３巻、８７５〜８９０頁（１９９６）
【非特許文献１８】Ｐｅ’ｅｒ，Ｄ．、Ｒｅｇｅｖ，Ａ．、Ｅｌｉｄａｎ，Ｇ．、Ｆｒｉｅｄｍａｎ，Ｎ．：「摂動発現プロフィールからのサブネットワークの推定」、Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ、第１７巻、２１５〜２２４頁（ＩＳＢＭ  ２００１）
【非特許文献１９】Ｓｏｍｅｒｅｎ，Ｅ．Ｖ．、Ｗｅｓｓｅｌｓ，Ｌ．、Ｒｅｉｎｄｅｒｓ，Ｍ．：「実験データからの遺伝子ネットワークの線形モデリング」、Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ、第１８巻、３５５〜３６６頁（ＩＳＢＭ  ２００２）
【非特許文献２０】Ｓｐｅｌｈｎａｎ，Ｐ．Ｔ．、Ｓｈｅｒｌｏｃｋ，Ｇ．、Ｚｈａｎｇ，Ｍ．Ｑ．、Ｉｙｅｒ，Ｖ．Ｒ．、Ａｎｄｅｒｓ，Ｋ．、Ｅｉｓｅｎ，Ｍ．Ｂ．、Ｂｒｏｗｎ，Ｐ．Ｏ．、Ｂｏｔｓｔｅｉｎ，Ｄ．、Ｆｕｔｃｈｅｒ，Ｂ．：「マイクロアレイハイブリッド形成による酵母菌「Ｓａｃｃｈａｒｏｍｙｃｅｓ  ｃｅｒｅｖｉｓｉａｅ」の細胞サイクル調節遺伝子の包括的同定」、Ｍｏｌｅｃｕｌａｒ  Ｂｉｏｌｏｇｙ  ｏｆ  ｔｈｅ  Ｃｅｌｌ、第９巻、３２７３〜３２９７（１９９８）
【非特許文献２１】Ｔｉｎｅｒｅｙ，Ｌ．、Ｋａｄａｎｅ，Ｊ．Ｂ．：「事後モーメント及び辺縁密度に対する正確な近似」、Ｊ．Ａｍｅｒ．Ｓｔａｔｉｓｔ．Ａｓｓｏｃ．、第８１巻、８２〜８６（１９８６）
【発明の開示】
【０００６】
いくつかの実施形態では、本発明は、ノンパラメトリック回帰によるベイジアンネットワークモデルでの時系列発現データの使用を含むものである。時系列発現データを用いて、循環調節構成要素を識別することができる。他の実施形態では、時間遅延情報をベイジアン／ノンパラメトリック回帰モデルに組込むことができ、これによって、次に、遺伝子間の非線形関係さえも抽出することができる。これらの実施形態のいくつかでは、常微分方程式モデルを代替として用いることができる。本発明人はまた、ベイジアン統計的見地から最適ネットワークを選択するための新しい判断基準を開発した。このような判断基準は、ノイズのあるデータに基づくネットワーク構造を最適化することができる。
本発明は、その特定的な実施形態を参照して説明するものである。本発明の付加的な態様は、実施例及び図面に見出される。
【発明を実施するための最良の形態】
【０００７】
一般的に、動的ベイジアンネットワークモデルは、遺伝子発現を判断するための任意の適切な方法を用いて得ることができる。いくつかの実施形態では、マイクロアレイ実験は、アレイに適用した単一のサンプルから多くの遺伝子を考察することができ、遺伝子発現の相対差を判断することが容易になるために望ましいものである。信号を反映する真の遺伝子発現から背景信号を差し引くことにより、及び／又は遺伝子発現（例えば、ｃｙ３／ｃｙ５）を測定するのに用いられるラベルの固有差を補正することにより、マイクロアレイ法の精度を向上させることが望ましい場合がある。
【０００８】
ベイジアンネットワークのフレームワークを用いて、遺伝子をランダム変数とみなし、その同時確率を条件付確率の積に分解する。例えば、ランダムベクトルの一連の観測値を有する場合には、所定の観測値を得る確率は、条件付確率密度に依存する可能性があると表すことができる。いくつかの実施形態では、変数間の関係を捕捉するためにノンパラメトリック回帰モデルを用いることができる。様々な図形的ツールを用いて、関係を解明することができる。例えば、遺伝子関係のグラフを形成するのに、多項式、フーリエ級数、回帰スプライン基底、Ｂ−スプライン基底、及びウェーブレット基底などを用いることができる。ネットワーク関係を解明するためのいくつかの方法は、本明細書においてその全内容が引用により組み込まれる米国特許出願出願番号第１０／２５９，７２３号に開示されている。適正なグラフを選択する際の困難な点の１つは、システムの分散及びノイズを適正に評価することである。
【０００９】
本発明のいくつかの実施形態では、ネットワークは、時系列研究からのデータを用いたノンパラメトリック回帰によるベイジアン推定を用いて構築することができる。多くの遺伝子ネットワークでは、介入により、他の遺伝子に変化が観測される前に特定の遺伝子の発現が変化することになる。時間的に遅く起こる介入の後のいくつかの遺伝子の発現は、発現が早い遺伝子が原因として関係する場合があると推定することができる。時系列情報は、「早い」つまり複数の遺伝子と「遅い」つまり遺伝子とを識別するのに有用である。遅い遺伝子の発現が変化しても、遅い遺伝子の発現よりも時間的に早く変化する早い遺伝子の発現の変化する原因となる可能性はあまりない。この仮定は、全ての場合に適用することができるわけではないが、早い遺伝子は、「下流」遺伝子である傾向が強い遅い遺伝子よりもネットワークの「上流」にある傾向が強いという可能性が高い。従って、遺伝子発現の時間関係は、ベイジアン推定及びノンパラメトリック回帰を修正して信頼性の高いネットワーク解を提供するのに有用であるとすることができる。
【００１０】
本発明の態様では、ベイジアンネットワーク及びノンパラメトリック回帰モデルを動的ベイジアンネットワークモデルにまで拡張し、これを用いて、時系列遺伝子発現データを有する場合の循環関係を構築することができる。遺伝子発現の変化の間の時間遅延に関する情報は、容易にモデルに含めることができ、このモデルは、遺伝子間の非線形関係さえも容易に抽出することができる。
いくつかの実施形態では、循環調節構成要素を有する遺伝子ネットワークを構築するために、常微分方程式モデル（Ｃｈｅｎ他［５］、ｄｅ  Ｈｏｏｎ他［８］）を用いることができる。しかし、このモデルは、線形システムに基づいており、複雑な現象を捕捉するには適切でない可能性がある。本発明人は、ベイジアン統計学的見地から最適ネットワークを選択するための新しい判断基準を導いた［２］。この判断基準により、ネットワーク構造を最適化することができ、これは、ノイズを伴うデータで表される遺伝子の相互作用を最も良く表すものである。新しい判断基準は、本明細書ではＢＮＲＣdynamicと呼ばれる。
【００１１】
ＢＮＲＣdynamicは、図１に示すような一次マルコフ関係を用いて評価することができる。このような関係では、上流遺伝子Ｘiは、１つ又はそれ以上の下流遺伝子Ｘ2に影響し（右矢印）、これがＸ3（図示せず）に影響し、Ｘnに影響するのが見られるまで連続して影響するように描かれている。Ｘ1がその「上流」遺伝子でない場合には、Ｘ1は、ネットワーク内の「親」遺伝子と呼ばれる。親遺伝子の影響を受けない遺伝子は、「ターゲット」遺伝子と呼ばれる。この関連の「ターゲット遺伝子」の使用は、潜在的薬物によるような介入のためのターゲットである遺伝子と混同しないように注意されたい。実際に、親遺伝子は、治療的介入のためのターゲットである場合がある。この体系の下では、Ｘnへの影響は、Ｘ1、Ｘ2等への影響が誘発されるまで見ることができない。図１は、並列又はフィードバックシステムが存在しない「直列」の原因／結果の関係を示し、一方、多くの遺伝子系では直列効果が存在し、かつ、２つ又はそれ以上の遺伝子が上流遺伝子により影響される可能性があり及び／又はそれら自身が下流遺伝子に影響を及ぼすことができる「並列」効果が存在することに注意されたい。更に、遺伝子Ｘaが別の遺伝子Ｘbに影響を及ぼすことができ、これがＸcに影響を及ぼすことができ、これ自身がＸa（又はＸb）に影響を及ぼすことができる循環効果（フィードバック）が存在することも可能である。更に、このようなフィードバックは、ＸcがＸaを刺激する正、又はＸcがＸaを阻害する「負」の何れかとすることができる。直列、並列、正のフィードバック及び負のフィードバックの両方の関係が存在する状況では、更に複雑になる可能性がある。
【００１２】
一般的に、時間点の間の関係は任意とすることができるが、いくつかの場合では、研究中の遺伝子及びその発現力学の生物学的影響の知識に基づいて予め選択した時間点を用いると有利である可能性がある。一次条件の下では、同時確率は、以下の実施例１の式（１）に示すように分解することができる。次に、実施例１の式（２）を用いて条件付確率を条件付確率の積に分解することができる。式（１）及び（２）が成り立つことが可能であり、確率測定の代わりに密度関数を用いることができる。従って、例えば、実施例１に説明する密度を用いて動的ベイジアンネットワークを表し、実施例１の式（３）による遺伝子及びその親遺伝子のローカルネットワーク構造に到達することができる。
【００１３】
例えば、実施例２に説明するように、ノンパラメトリック回帰による動的ベイジアンモデルを適用することができる。実験データが収集された状態で、ネットワークの解は、統計的モデル選択の問題と考えることができる。いくつかの実施形態では、ベイジアン手法を用いてこの問題を解き、動的ベイジアンネットワーク及びノンパラメトリック回帰法の適合度を評価するための判断基準を導出することができる。事前分配を仮定すれば、辺縁尤度及び事後確率は、実施例２の式（４）により判断することができる。引き続く遺伝子ネットワークの構築には、式（４）に説明するような高次元積分の計算が含まれる。いくつかの実施形態では、例えば、積分のためのラプラス法を用いて積分の近似値を判断することができる。従って、実施例２の式（５）に示すような判断基準ＢＮＲＣdynamicを解くことができる。
【００１４】
ＢＮＲＣdynamicを実験システムに適用するために、例えば、遺伝子系に影響を及ぼした後のいくつかの時間点でｃＤＮＡマイクロアレイデータを実験的に得ることができる。曲線を平滑化するために、スプライン関数、例えば実施例３に示すようなＢスプラインを用いることができる。ＢＮＲＣdynamicは、実施例３の式（６）により分解することができる。最適なネットワーク関係は、ＢＮＲＣdynamicを最小にした時に得られるものである。
【００１５】
ノンパラメトリック回帰による動的ベイジアンネットワークモデル及び判断基準ＢＮＲＣdynamicを用いて、ネットワーク学習過程を定式化することができる。しかし、どの遺伝子が親遺伝子であり、どれがターゲット遺伝子であるかを判断する段階は、全ての可能性のある遺伝子の組合せ及び関係を考慮する場合は時間を要する可能性がある。必要な解析の数を低減するために、候補親遺伝子を選択することができる。続いて、貪欲ヒルクライムアルゴリズムを用いることができる。ＢＮＲＣdynamicを計算し、次に、付加的な親遺伝子を追加するか又は削除し、実施例３の段階２により、ＢＮＲＣdynamicを再評価する。適切な収束が見つかるまでこの過程を繰り返す。次に、計算の順序を並べ替えてＢＮＲＣdynamicを再評価する。最適なネットワークでは、ＢＮＲＣdynamicが最も小さくなる。
【００１６】
上述の方法の特定の解説を図２ａ及び２ｂの実施例４に示す。この方法の効率は、「Ｓａｃｃｈａｒｏｍｙｃｅｓ  ｃｅｒｅｖｉｓｉａｅ」からの遺伝子発現データの分析を通じて示される。図２ａは、細胞サイクルの調節に関わる「Ｓ．ｃｅｒｅｖｉｓｉａｅ」遺伝子の群を示すものである。遺伝子は、関連する全代謝経路に基づく群で示されており、サイクリン依存プロテインキナーゼ遺伝子（ＹＢＲ１６０ｗ）に注目している。親／ターゲット遺伝子ネットワーク関係は、図２ａに基づいて未知であることに注意されたい。これと対照的に、本発明の方法を用いて、これらの遺伝子のネットワーク関係を評価することができ、これは図２ｂに示されている。
【００１７】
別の実施例を図３ａ〜３ｃに示す。図３ａは、代謝経路に関わる遺伝子を示すものである。図３ａは、遺伝子ネットワーク関係を示してはいない。図３ｂは、ノンパラメトリック回帰を有するが、ＢＮＲＣdynamicを考慮していないベイジアンネットワーク分析を用いて得られるネットワークの解を示すものである。図３ｃは、ＢＮＲＣdynamicを最小にすることにより得られるネットワーク解を示すものである。図３ｂでは、ネットワーク関係が比較的単純であり、図３ｂに示すものと比べて偽の正の関係（ｘ）が遥かに少ないことに注意されたい。
【００１８】
ネットワークでの遺伝子群の間の境界は、当業技術で公知の方法、例えば、ブートストラップ法を用いて判断することができる。このような方法には、以下の段階を用いて辺縁の強度を判断する段階が含まれる。
（１）オリジナル遺伝子ライブラリ発現データから置換を用いて何回か無作為にサンプリングすることにより、ブートストラップ遺伝子発現マトリックスを準備する段階、
（２）遺伝子i及び遺伝子jのための遺伝子ネットワークを推定する段階、
（３）段階（１）及び（２）をＴ回繰り返すことにより、Ｔの遺伝子ネットワークを生成する段階、及び
（４）（ｔ1＋ｔ2）／Ｔとして遺伝子ｉと遺伝子ｊの間のブートストラップ辺縁強度を計算する段階。
【００１９】
ベイジアン及びブーリアンネットワークのような他のネットワーク推定法と比較した新しい方法の利点には、（１）時間情報を容易に組み込むことができること、（２）マイクロアレイデータは、区分化のような余分なデータ予備処理を行うことなく連続データとして分析することができること、及び（３）見出される偽の正の関係が少ないことが含まれる。本発明の実施形態により、非線形関係でさえも検出してモデリングすることができる。本発明の方法は、遺伝子ネットワークの分析及び重要な遺伝子の遺伝子発現を制御する特定の遺伝子をターゲットにする新しい医薬品の開発に有用である。すなわち、本発明の方法により、薬物ターゲットを識別するのに必要な時間を低減することができ、従って、新しい治療法を開発するのに必要な時間を低減することができる。
本発明の方法の他の態様は、以下の実施例で説明する。
【００２０】
実施例
以下に示す実施例は、本発明の特定の実施形態を表すものである。当業者は、過度の実験を行うことなく本発明の他の態様を作り出すことができる。全てのこのような実施形態は、本発明の一部と見なされるものである。
【実施例１】
【００２１】
ベイジアンネットワーク及びノンパラメトリック回帰
ｎ及びｐが各々マイクロアレイ及び遺伝子の数であるｎ×ｐマイクロアレイ遺伝子発現データマトリックスＸがあるとする。通常、遺伝子の数ｐは、マイクロアレイの数ｎよりも遥かに大きい。ベイジアンネットワークに基づく遺伝子ネットワークの推定では、遺伝子は、ランダム変数と見なされる。密度又は確率関数で説明される統計的モデルを用いることにより遺伝子ネットワークをモデリングする場合には、統計的モデルは、ｐランダム変数を含む必要がある。しかし、ここではｎサンプルしかなく、かつｎは通常はｐよりも遥かに小さい。このような場合には、モデルの推定は、モデルが多くのパラメータを有し、サンプルの数がパラメータを推定するのに十分でないために非常に困難であるか又は不可能である。このようなモデリングにおいて、ベイジアンネットワークモデルが提唱されている。
【００２２】
動的ベイジアンネットワークにおいては、時系列データを考えて、Ｘのｉ番目の列ベクトルｘiが時間ｉでのｐ遺伝子の状態に対応するものとする。時間依存性に関しては、図１に説明した一次マルコフ関係を考える。この条件では、同時確率は、次のように分解することができる。
Ｐ（Ｘ11，・・・，ＸnP）＝Ｐ（Ｘ1）Ｐ（Ｘ2｜Ｘ1）×・・・×Ｐ（Ｘn｜Ｘn-1）      （１）
ここで、Ｘi＝（Ｘi1，・・・，Ｘip）は、時間ｉでのｐ遺伝子のランダム変数ベクトルである。また、条件付確率Ｐ（ｘi｜Ｘi-1）は、以下の形の条件付確率の積に分解することができる。
Ｐ（Ｘi｜Ｘi-1）＝Ｐ（Ｘil｜Ｐi-l,l）×・・・×Ｐ（Ｘip｜Ｐi-1,p）      （２）
ここで、Ｐi-1,jは、時間ｉ−１でのｊ番目の遺伝子の親遺伝子の状態ベクトルである。式（１）及び（２）は、確率測度の代わりに密度関数を用いる場合にも成り立つものである。従って、次に、動的ベイジアンネットワークは、以下のように密度を用いて表すことができる。
【００２３】
【数１】
2006506746.tif 000003
【００２４】
ここで（２）から分解が得られる。
ｆi（ｘi｜ｘi-1）＝ｇ1（ｘi1｜Ｐi-1,1）×・・・×ｇp（ｘip｜Ｐi-1,p）
ここで、
【００２５】
【数２】
2006506746.tif 000004
【００２６】
は、親遺伝子のｑj次元観測ベクトルである。
【００２７】
ｘijとｐi-1,jの間の関係をモデリングするために、次のようなノンパラメトリック付加回帰モデルを用いる。
【００２８】
【数３】
2006506746.tif 000005
【００２９】
ここで、εijは、独立にかつ通常は平均０及び分散
【００３０】
【数４】
2006506746.tif 000006
【００３１】
に依存する。ここで、ｍjk（・）はＲからＲまでの平滑関数であり、基底関数の線形結合を用いて表すことができる。
【００３２】
【数５】
2006506746.tif 000007
【００３３】
ここで、
【００３４】
【数６】
2006506746.tif 000008
【００３５】
は未知の係数パラメータであり、
【００３６】
【数７】
2006506746.tif 000009
【００３７】
は、基底関数の規定の組である。次に、以下の形の動的ベイジアンネットワーク及びノンパラメトリック回帰モデルを定義する。
【００３８】
【数８】
2006506746.tif 000010
【００３９】
ここで、
【００４０】
【数９】
2006506746.tif 000011
【００４１】
である。ｊ番目の遺伝子に親遺伝子がない場合は、μ（ｐi-1,j）は、定数μjになる。
更に、ｆ1（ｘ1）＝ｇ1（ｘ11）×・・・×ｇ1（ｘ1p）と仮定すると、結合密度ｆ（ｘ11，・・・，ｘnpiθG）は、以下のように書き直すことができる。
【００４２】
【数１０】
2006506746.tif 000012
【００４３】
ここで、Ｐoj＝θである。従って、ｇj（ｘij｜ｐi-1j；θj）は、ｊ番目の遺伝子及びその親遺伝子のローカル構造を表すものである。
【実施例２】
【００４４】
ネットワークを選択するための判断基準の導出
上述の節で導入された動的ベイジアンネットワーク及びノンパラメトリック回帰モデルは、ネットワーク構造が固定される時に構築され、適切な手順により推定することができる。しかし、遺伝子ネットワークは、一般的に未知であり、データに基づき最適なネットワークを推定すべきである。この問題は、統計的モデルの選択の問題としてみることができる（例えば、Ａｋａｉｋｅ［１］、Ｋｏｎｉｓｈｉ及びＫｉｔａｇａｗａ［１７］、Ｂｕｒｎｈａｍ及びＡｎｄｅｒｓｏｎ［４］、Ｋｏｎｉｓｈｉ［１６］を参照）。本発明人は、この問題をベイジアン統計学的手法で解決し、動的ベイジアンネットワーク及びノンパラメトリック回帰モデルが良好であるか否かを評価するための判断基準を導くものである。
【００４５】
π（θG｜λ）を動的ベイジアンネットワーク及びノンパラメトリック回帰モデルのパラメータθGの事前分布とし、ｌｏｇπ（θG｜λ）＝Ｏ（ｎ）とする。辺縁尤度は、以下のように表すことができる。
【００４６】
【数１１】
2006506746.tif 000013
【００４７】
従って、データが与えられると、ネットワークＧの事後確率は、以下のようになる。
【００４８】
【数１２】
2006506746.tif 000014
【００４９】
ここで、πprior（Ｇ）は、ネットワークＧの事前確率である。（４）の分母はモデル評価には関係がない。従って、ネットワークの評価は、分子のマグニチュードに依存する。従って、以下を最大にするものとして最適ネットワークを選択することができる。
【００５０】
【数１３】
2006506746.tif 000015
【００５１】
ネットワーク選択判断基準を構築するための本質的な点は、高次元積分を計算する方法であることは明らかである。Ｉｎｍｏｔｏ他［１４，１５］は、積分のためのラプラス近似法を用いており（Ｔｉｎｅｒｅｙ及びＫａｄａｎｅ［２１］、Ｄａｖｉｓｏｎ［６］も参照されたい）、この技術は、動的ベイジアンネットワークモデル及びノンパラメトリック回帰に直接適用することができる。従って、以下の形のＢＮＲＣdynamicという名称の判断基準を有することになる。
【００５２】
【数１４】
2006506746.tif 000016
【００５３】
ここで、ｒはθGの次元であり、また、
【００５４】
【数１５】
2006506746.tif 000017
【００５５】
であり、
【００５６】
【数１６】
2006506746.tif 000018
【００５７】
は、ｌλ（θG｜Ｘn）のモードである。最適なグラフは、判断基準ＢＮＲＣdynamic（５）が最小になるように選択される。
【実施例３】
【００５８】
遺伝子ネットワークの推定
この節では、ｃＤＮＡマイクロアレイ時系列遺伝子発現データから遺伝子ネットワークを推定するための具体的な戦略を示す。
３．１．ノンパラメトリック回帰
第２節に説明した平滑関数ｍjk（・）を構築するために基底関数を用いる。本稿では、基底関数としてＢスプライン（ｄｅ  Ｂｏｏｒ［７］）を用いる。Ｄｅ  Ｂｏｏｒのアルゴリズム（ｄｅ  Ｂｏｏｒ［７］、第１０章１３０頁（３））は、任意の次数のＢスプラインを計算するために有用な方法である。本発明人は、等距離ノットを有する２０個のＢスプラインを用いる（Ｂスプラインの詳細は、Ｄｉｅｒｃｋｘ［１０］、Ｅｉｌｅｒ及びＭａｒｘ［１１］も参照されたい）。
【００５９】
３．２．モデル内のパラメータに対する事前分布
パラメータθGに関する事前分布に関して、パラメータベクトルθjは互いに独立であるとすると、事前分布は、以下のように分解することができる。
【００６０】
【数１７】
2006506746.tif 000019
【００６１】
事前分布πj（θj｜λj）は、以下のように因数分解される。
【００６２】
【数１８】
2006506746.tif 000020
【００６３】
ここで、λjkは超パラメータである。本発明人は、γjkに関する事前分布として以下のような特異Ｍjk変量正規分布を用いる。
【００６４】
【数１９】
2006506746.tif 000021
【００６５】
ここで、Ｋjkは、以下を満たすＭjk×Ｍjk対称正半有限マトリックスである。
【００６６】
【数２０】
2006506746.tif 000022
【００６７】
このθGに関する事前確率の設定は、ｌｍｏｔｏ他［１４，１５］と同じであり、詳細は、その論文に説明されている。
【００６８】
３．３．提案判断基準
第４．２節の事前分布を用いることにより、ＢＮＲＣdynamicは、次のように分解することができる。
【００６９】
【数２１】
2006506746.tif 000023
【００７０】
ここで、
【００７１】
【数２２】
2006506746.tif 000024
【００７２】
は、ｊ番目の遺伝子のローカル判断基準スコアであり、以下のように定義される。
【００７３】
【数２３】
2006506746.tif 000025
【００７４】
ここで、γjはθjの次元であり、
【００７５】
【数２４】
2006506746.tif 000026
【００７６】
であり、
【００７７】
【数２５】
2006506746.tif 000027
【００７８】
は、
【００７９】
【数２６】
2006506746.tif 000028
【００８０】
のモードである。ここで、πprior（Ｌj）は、以下を満足する事前確率である。
【００８１】
【数２７】
2006506746.tif 000029
【００８２】
ローカル構造の事前確率πprior（Ｌj）は、πprior（Ｌj）＝ｅｘｐ｛−（ｊ番目の遺伝子の親遺伝子の数）｝と設定される。
【００８３】
動的ベイジアンネットワーク及びノンパラメトリック回帰モデルをこの提案判断基準ＢＮＲＣｄ１／ｎａｍｉｃと共に用いることにより、ネットワーク学習過程を次のように定式化することができる。ネットワーク構造の最適化は、ターゲット遺伝子を調節する親遺伝子の選択と同等であることは（３）及び（６）から明らかである。しかし、親遺伝子として全ての可能性のある遺伝子を考えると、それは時間を要する仕事である。従って、本発明人は、候補親遺伝子を選択することにより学習空間を削減する。この段階の後に、良好なネットワークを見つけるために貪欲ヒルクライムアルゴリズムを用いる。本発明人のアルゴリズムは、次のように表すことができる。
段階１：前処理段階
（ｉ，ｊ）番目の要素がグラフ「遺伝子i→遺伝子j」のＢＮＲＣスコアであるｐ×ｐマトリックスを作り、ＢＮＲＣスコアが小さくなる遺伝子の親遺伝子の候補の組を定義する。親遺伝子１０の候補の組の要素の数を設定する。
段階２：学習段階
貪欲ヒルクライムアルゴリズムに対しては、空のネットワークから開始し、以下の段階を繰り返す。
段階２−１：遺伝子iに対して、親遺伝子を「加える」又は親遺伝子を「削除する」という２つの手順のうちＢＮＲＣ動的スコアが小さくなる方を実行する。
段階２−２：遺伝子の規定の計算順序に対して段階２−１を適切な収束判断基準を満たすまで繰り返す。
段階２−３：良好な解を見つけるために計算順序を入れ変え、段階２−１及び段階２−２を繰り返す。
段階２−４：ＢＮＲＣdynamicスコアが最も小さくなる最適ネットワークを選択する。
【実施例４】
【００８４】
計算実験
Ｓｐｅｌｌｍａｎ他［２０］により收集された「Ｓａｃｃｈａｒｏｍｙｃｅｓ  ｃｅｒｅｖｉｓｉａｅ」細胞サイクル遺伝子発現データの分析を通じて本発明の一実施形態を明らかにする。このデータは、２つの短い時系列（２つの時間点：ｃｌｎ３、ｃｌｂ２）及び４つの中程度の時系列（１８、２４、１７、及び１４個の時間点：ａｌｐｈａ、ｃｄｃ１５、ｃｄｃ２８、及びｅｌｕ）を包含するものである。遺伝子ネットワークの推定では、４つの中程度の時系列を用いた。４つの時系列を組み合わせるために、ノンパラメトリック回帰モデルを当て嵌める時に各時系列に対してターゲット遺伝子の最初の観測値及び親遺伝子の最後の観測値を無視した。
【００８５】
最初に、ＫＥＧＧデータベース［２２］に集められた細胞サイクル経路に着目した。ターゲットネットワークは、ＣＤＣ２８（ＹＢＲ１６Ｏｗ：サイクリン依存プロテインキナーゼ）の周辺にある。このネットワークは４５遺伝子を含み、ＫＥＧＧに登録された経路が図２（ａ）に示され、推定ネットワークが図２（ｂ）に示されている。点線の円内のエッジは正確なエッジと見なすことができる。従って、正確な関係のいくつかがモデリングされた。エッジに隣接する円は正確な推定を表すものである。三角形は、正しい方向の逆又は飛越しを表すものである。「ｘ」符号は不正確な関係を表すものである。
【００８６】
本発明の方法を明らかにするのに用いる第２の例は、ＤｅＲｉｓｉ他［９］により報告されている代謝経路である。このネットワークには５７遺伝子が含まれ、ターゲット経路は図３（ａ）に示されている。
このデータにベイジアンネットワーク及びノンパラメトリック回帰モデル［１４，１５］を適用し、得られるネットワークを図３（ｂ）に示す。図３（ｃ）のネットワークは、動的ベイジアンネットワーク及びノンパラメトリック回帰モデルにより得たものである。ｃＤＮＡマイクロアレイデータから代謝経路を推定するのは困難である。しかし、本発明のモデルでは、遺伝子間の関係が正確に検出された。ベイジアンネットワーク及びノンパラメトリック回帰に比較して、「ｘ」符号で図３（ｃ）に示す本方法の偽の正の関係の数は、図３（ｂ）に示す数よりも遥かに少なかった。
ここに列挙する全ての参考文献は、本明細書においてその全内容が組み込まれている。
【００８７】
参考文献
１．Ａｋａｉｋｅ、Ｊ．：「情報理論と最大尤度原理の拡張」、掲載：Ｐｅｔｒｏｖ、Ｂ．Ｎ．、Ｃｓａｋｉ、Ｆ．（編）：第２回情報理論国際シンポジウム、Ａｋａｄｅｍｉａｉ  Ｋｉａｄｏ、ブダペスト、２６７〜２８１頁（１９７３）
２．Ｂｅｒｇｅｒ、Ｊ．Ｏ．：「統計的確定理論とベイジアン分析」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ニューヨーク（１９８５）
３．Ｂｉｌｍｅｓ、Ｊ．Ａ．：「動的ベイジアンマルチネット」、第１６回人工知能の不確定性に関する会議紀要、３８〜４５頁（２０００）
４．Ｂｕｒｎｈａｍ、Ｋ．Ｐ．、Ａｎｄｅｒｓｏｎ、Ｄ．Ｒ．：「モデル選択と推定、実用情報−理論的手法」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ニューヨーク（１９８８）
５．Ｃｈｅｎ、ＴＩ．、Ｈｅ、Ｈ．Ｌ．、Ｃｈｕｒｃｈ、Ｇ．Ｍ．：「微分方程式を用いる遺伝子発現のモデリング」、バイオコンピューティングに関する太平洋シンポジウム紀要、第４巻、２９〜４０（１９９９）
６．Ｄａｖｉｓｏｎ，Ａ．Ｃ．：「近似予測尤度」、Ｂｉｏｍｅｔｒｉｋａ、第７３巻、３２３〜３３２（１９８６）
７．ＤｅＢｏｏｒ，Ｃ．：「スプラインの実用ガイド」、Ｓｐｒｉｎｇｅｒ−Ｖｅｒｌａｇ、ベルリン（１９７８）
８．Ｄｅ  Ｈｏｏｎ，Ｍ．Ｊ．Ｌ．、Ｉｍｏｔｏ，Ｓ．、Ｋｏｂａｙａｓｈｉ，Ｋ．、Ｏｇａｓａｗａｒａ，ＮＨ．、Ｍｉｙａｎｏ，Ｓ．：「微分方程式を用いる「Ｂａｃｉｌｌｕｓ  ｓｕｂｔｉｌｉｓ」の時間順遺伝子発現データからの遺伝子調節ネットワークの推定」、バイオコンピューティングに関する太平洋シンポジウム紀要、第８巻、２００３、近刊
９．ＤｅＲｉｓｉ，Ｊ．、Ｌｙｅｒ、Ｖ．Ｒ．、Ｂｒｏｗｎ、Ｐ．Ｏ．：「ゲノム規模の遺伝子発現の代謝及び遺伝子制御の探求」、Ｓｃｉｅｎｃｅ、第２７８巻、６８０〜６８６頁（１９９７）
１０．Ｄｉｅｒｃｋｘ、Ｐ．：「スプラインを用いる曲線及び表面当て嵌め」、Ｏｘｆｏｒｄ（１９９３）
１１．Ｅｉｌｅｒ、Ｐ．Ｈ．Ｃ．、Ｍａｉｘ、Ｂ．：「Ｂスプラインを用いる柔軟な平滑化及び不利益（討論）」、Ｓｔａｔｉｓｔｉｃａｌ  Ｓｃｉｅｎｃｅ、第１１巻、８９〜１２１頁（１９９６）
１２．Ｆｒｉｅｄｍａｎ，Ｎ．、Ｍｕｒｐｈｙ，Ｋ．、Ｒｕｓｓｅｌｌ，Ｓ．：「動的確率ネットワークの構造の学習」、人工知能における不確定性に関する会議紀要、１３９〜１４７頁（１９９８）
１３．Ｆｉｒｅｄｍａｎ，Ｎ．、Ｌｉｎｉａｌ，Ｍｌ、Ｎａｃｈｍａｎ，Ｉ．、Ｐｅ’ｅｒ，Ｄ．：「発現データの分析のためのベイジアンネットワークの使用」、Ｊ．Ｃｏｍｐ．Ｂｉｏｌ．、第７巻、６０１〜６２０（２０００）
１４．Ｉｍｏｔｏ，Ｓ．、Ｇｏｔｏ，Ｔ．、Ｍｉｙａｎｏ，Ｓｌ．：「ベイジアンネットワーク及びノンパラメトリック回帰を用いる遺伝子間の遺伝子ネットワーク及び機能的構造の推定」、バイオコンピューティングに関する太平洋シンポジウム紀要、第７巻、１７５〜１８６（２００２）
１５．Ｉｍｏｔｏ，Ｓ．、Ｋｉｍ，Ｓ．、Ｇｏｔｏ，Ｔ．、Ａｂｕｒａｔａｎｉ，Ｓ．、Ｔａｓｈｉｒｏ，Ｋ．、Ｋｕｈａｒａ，Ｓ．、Ｍｊｉｙａｎｏ，Ｓ．：「遺伝子ネットワークの非線形モデリングのためのベイジアンネットワーク及びノンパラメトリック異分散回帰」、ＩＥＥＥ  Ｃｏｍｐｕｔｅｒ  Ｓｏｃｉｅｔｙ  Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ  Ｃｏｎｆｅｒｅｎｃｅ紀要、２１９〜２２７頁（２００２）
１６．Ｋｏｎｉｓｈｉ，Ｓ．：「統計学的モデルの評価及び情報判断基準」、掲載：Ｇｈｏｓｈ，Ｓ．（編）「多変量分析、実験計画及び調査サンプリング」、Ｍａｒｃｅｌ  Ｄｅｋｋｅｒ、ニューヨーク、３６９〜３９９頁（１９９９）
１７．Ｋｏｎｉｓｈｉ，Ｓ．、Ｋｉｔａｇａｗａ，Ｇ．：「モデル選択における汎用情報判断基準」、Ｂｉｏｍｅｔｒｉｋａ、第８３巻、８７５〜８９０頁（１９９６）
１８．Ｐｅ’ｅｒ，Ｄ．、Ｒｅｇｅｖ，Ａ．、Ｅｌｉｄａｎ，Ｇ．、Ｆｒｉｅｄｍａｎ，Ｎ．：「摂動発現プロフィールからのサブネットワークの推定」、Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ、第１７巻、２１５〜２２４頁（ＩＳＢＭ  ２００１）
１９．Ｓｏｍｅｒｅｎ，Ｅ．Ｖ．、Ｗｅｓｓｅｌｓ，Ｌ．、Ｒｅｉｎｄｅｒｓ，Ｍ．：「実験データからの遺伝子ネットワークの線形モデリング」、Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ、第１８巻、３５５〜３６６頁（ＩＳＢＭ  ２００２）
２０．Ｓｐｅｌｈｎａｎ，Ｐ．Ｔ．、Ｓｈｅｒｌｏｃｋ，Ｇ．、Ｚｈａｎｇ，Ｍ．Ｑ．、Ｉｙｅｒ，Ｖ．Ｒ．、Ａｎｄｅｒｓ，Ｋ．、Ｅｉｓｅｎ，Ｍ．Ｂ．、Ｂｒｏｗｎ，Ｐ．Ｏ．、Ｂｏｔｓｔｅｉｎ，Ｄ．、Ｆｕｔｃｈｅｒ，Ｂ．：「マイクロアレイハイブリッド形成による酵母菌「Ｓａｃｃｈａｒｏｍｙｃｅｓ  ｃｅｒｅｖｉｓｉａｅ」の細胞サイクル調節遺伝子の包括的同定」、Ｍｏｌｅｃｕｌａｒ  Ｂｉｏｌｏｇｙ  ｏｆ  ｔｈｅ  Ｃｅｌｌ、第９巻、３２７３〜３２９７（１９９８）
２１．Ｔｉｎｅｒｅｙ，Ｌ．、Ｋａｄａｎｅ，Ｊ．Ｂ．：「事後モーメント及び辺縁密度に対する正確な近似」、Ｊ．Ａｍｅｒ．Ｓｔａｔｉｓｔ．Ａｓｓｏｃ．、第８１巻、８２〜８６（１９８６）
【図面の簡単な説明】
【００８８】
【図１】遺伝子発現における時間力学を示す概略図である。
【図２ａ】ＫＥＧＧに集められた酵母菌の細胞サイクル調節に関わる遺伝子のネットワーク関係を示し、サイクリン依存プロテインキナーゼ経路の遺伝子を示す図である。
【図２ｂ】ＫＥＧＧに集められた酵母菌の細胞サイクル調節に関わる遺伝子のネットワーク関係を示し、サイクリン依存プロテインキナーゼの調節に関わる図２ａに説明の遺伝子間のネットワーク関係を示す図である。
【図３ａ】代謝経路に関わる酵母菌遺伝子のネットワーク関係の図を示し、代謝経路に関わるいくつかの遺伝子を示す図である。
【図３ｂ】代謝経路に関わる酵母菌遺伝子のネットワーク関係の図を示し、ベイジアン／ノンパラメトリック回帰モデルからの導出された図３ａに説明の遺伝子間のネットワーク関係を示す図である。
【図３ｃ】代謝経路に関わる酵母菌遺伝子のネットワーク関係の図を示し、動的ベイジアン／ノンパラメトリック回帰モデルから導出された図３ａに説明の遺伝子間のネットワーク関係を示す図である。
【符号の説明】
【００８９】
ｘ  偽の正の関係
【図１】
2006506746.tif 000030
2006506746.tif 000031
2006506746.tif 000032
2006506746.tif 000033
2006506746.tif 000034
2006506746.tif 000035
【国際調査報告】
2006506746.tif 000036