(57)【要約】  エネルギー関数最適化システムは、分子の安定構造とその分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いるエネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標としてエネルギー関数の最適解を探索する多変数最適化指標最適解探索手段を有する。
2009090796.tif 000002
【特許請求の範囲】
【請求項１】
分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索手段を有するエネルギー関数最適化システム。
【請求項２】
前記エネルギー関数最適化システムが、  前記第１の指標を生成する第１の指標生成手段を有する、請求項１記載のエネルギー関数最適化システム。
【請求項３】
前記エネルギー関数最適化システムが、  前記第２の指標を生成する第２の指標生成手段を有する、請求項１または２記載のエネルギー関数最適化システム。
【請求項４】
前記第１の指標生成手段が、  以下の式１によって求められるｚスコアを求める、請求項２または３に記載のエネルギー関数最適化システム。【数１】
2009090796.tif 000010
（Ｅｅｘｐは前記分子の安定構造におけるエネルギーを表す。＜Ｅｃａｌｃ＞は前記分子の不安定構造におけるエネルギーの平均値を表す。σｃａｌｃは前記分子の不安定構造におけるエネルギーの標準偏差を表す。）
【請求項５】
前記第２の指標生成手段が、  前記複数種類の分子の結合エネルギーの実験値と、前記エネルギー関数によって求められる前記複数種類の分子の結合構造の計算値と、の相関係数を求める、請求項３または４に記載のエネルギー関数最適化システム。
【請求項６】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標の単調関数との積である、請求項１乃至５のいずれか一項に記載のエネルギー関数最適化システム。
【請求項７】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標のＮ乗（Ｎは任意）との積である、請求項１乃至６のいずれか一項に記載のエネルギー関数最適化システム。
【請求項８】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標の二乗との積である、請求項１乃至７のいずれか一項に記載のエネルギー関数最適化システム。
【請求項９】
前記多変数最適化指標最適解探索手段が、  前記エネルギー関数の最適解を与えうる全パラメータに対して評価を行い、前記エネルギー関数の最適解を求める、請求項１乃至８のいずれか一項に記載のエネルギー関数最適化システム。
【請求項１０】
前記多変数最適化指標最適解探索手段が、  最適化アルゴリズムを用いて前記エネルギー関数の近似解を求める、請求項１乃至８のいずれか一項に記載のエネルギー関数最適化システム。
【請求項１１】
前記分子及び前記複数種類の分子が、  タンパク質及び前記タンパク質に結合する化合物の少なくともいずれか一方である、請求項１乃至１０のいずれか一項に記載のエネルギー関数最適化システム。
【請求項１２】
コンピュータに、  分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索処理を実行させるためのエネルギー関数最適化プログラム。
【請求項１３】
前記エネルギー関数最適化プログラムが、  前記第１の指標を生成する第１の指標生成処理を実行させる、請求項１２記載のエネルギー関数最適化プログラム。
【請求項１４】
前記エネルギー関数最適化プログラムが、  前記第２の指標を生成する第２の指標生成処理を実行させる、請求項１２または１３記載のエネルギー関数最適化プログラム。
【請求項１５】
前記第１の指標生成処理が、  式１によって求められるｚスコアを求める、請求項１３または１４に記載のエネルギー関数最適化プログラム。【数２】
2009090796.tif 000011
（Ｅｅｘｐは前記分子の安定構造におけるエネルギーを表す。＜Ｅｃａｌｃ＞は前記分子の不安定構造におけるエネルギーの平均値を表す。σｃａｌｃは前記分子の不安定構造におけるエネルギーの標準偏差を表す。）
【請求項１６】
前記第２の指標生成処理が、  前記複数の分子の結合エネルギーの実験値と、前記エネルギー関数によって求められる前記複数の分子の結合構造の計算値と、の相関係数を求める、請求項１４または１５に記載のエネルギー関数最適化プログラム。
【請求項１７】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標の単調関数との積である、請求項１２乃至１６のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項１８】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標のＮ乗（Ｎは任意）との積である、請求項１２乃至１７のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項１９】
前記多変数最適化指標が、  前記第１の指標と前記第２の指標の二乗との積である、請求項１２乃至１８のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項２０】
前記多変数最適化指標最適解探索処理が、  前記エネルギー関数の最適解を与えうる全パラメータに対して評価を行い、前記エネルギー関数の最適解を求める、請求項１２乃至１９のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項２１】
前記多変数最適化指標最適解探索処理が、  最適化アルゴリズムを用いて前記エネルギー関数の近似解を求める、請求項１２乃至１９のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項２２】
前記分子及び前記複数種類の分子が、  タンパク質及び前記タンパク質に結合する化合物の少なくともいずれか一方である、請求項１２乃至２１のいずれか一項に記載のエネルギー関数最適化プログラム。
【請求項２３】
分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索ステップを有することを特徴とするエネルギー関数最適化方法。
【請求項２４】
前記エネルギー関数最適化方法が、  前記第１の指標を生成する第１の指標生成ステップを有する、請求項２３記載のエネルギー関数最適化方法。
【請求項２５】
前記エネルギー関数最適化方法が、  前記第２の指標を生成する第２の指標生成ステップを有する、請求項２３または２４記載のエネルギー関数最適化方法。
【請求項２６】
前記第１の指標生成ステップが、  式１によって求められるｚスコアを求める、請求項２４または２５に記載のエネルギー関数最適化方法。【数３】
2009090796.tif 000012
（Ｅｅｘｐは前記分子の安定構造におけるエネルギーを表す。＜Ｅｃａｌｃ＞は前記分子の不安定構造におけるエネルギーの平均値を表す。σｃａｌｃは前記分子の不安定構造におけるエネルギーの標準偏差を表す。）
【請求項２７】
前記第２の指標生成ステップが、前記複数種類の分子の結合エネルギーの実験値と、前記エネルギー関数によって求められる前記複数種類の分子の結合構造の計算値と、の相関係数を求める、請求項２５または２６に記載のエネルギー関数最適化方法。
【請求項２８】
前記多変数最適化指標最適解探索ステップが、  前記第１の指標と前記第２の指標の単調関数との積である、請求項２３乃至２７のいずれか一項に記載のエネルギー関数最適化方法。
【請求項２９】
前記多変数最適化指標が、前記第１の指標と前記第２の指標のＮ乗（Ｎは任意）との積である、請求項２３乃至２８のいずれか一項に記載のエネルギー関数最適化方法。
【請求項３０】
前記多変数最適化指標最適解探索ステップが、  前記第１の指標と前記第２の指標の二乗との積である、請求項２３乃至２９のいずれか一項に記載のエネルギー関数最適化方法。
【請求項３１】
前記多変数最適化指標最適解探索ステップが、  前記エネルギー関数の最適解を与えうる全パラメータに対して評価を行い、前記エネルギー関数の最適解を求める、請求項２３乃至３０のいずれか一項に記載のエネルギー関数最適化方法。
【請求項３２】
前記多変数最適化指標最適解探索ステップが、  最適化アルゴリズムを用いて前記エネルギー関数の近似解を求める、請求項２３乃至３０のいずれか一項に記載のエネルギー関数最適化方法。
【請求項３３】
前記分子及び前記複数種類の分子が、  タンパク質及び前記タンパク質に結合する化合物の少なくともいずれか一方である、請求項２３乃至３２のいずれか一項に記載のエネルギー関数最適化方法。
【請求項３４】
分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーション手段を有する分子シミュレーションシステム。
【請求項３５】
前記分子シミュレーションシステムが、  前記第１の指標を生成する第１の指標生成手段、前記第２の指標を生成する第２の指標生成手段、前記多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索手段、の少なくとも１つを有する、請求項３４記載の分子シミュレーションシステム。
【請求項３６】
コンピュータに、  分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーション処理を実行させる、分子シミュレーションプログラム。
【請求項３７】
前記分子シミュレーションプログラムが、  前記第１の指標を生成する第１の指標生成処理、前記第２の指標を生成する第２の指標生成処理、前記多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索処理、の少なくとも１つを実行させる、請求項３６記載の分子シミュレーションプログラム。
【請求項３８】
分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーションステップを有する分子シミュレーション方法。
【請求項３９】
前記分子シミュレーション方法が、  前記第１の指標を生成する第１の指標生成ステップ、前記第２の指標を生成する第２の指標生成ステップ、前記多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索ステップ、の少なくとも１つを有する、請求項３８記載の分子シミュレーション方法。
【発明の詳細な説明】
【技術分野】
【０００１】  本発明は、分子シミュレーションにおいて、特に分子シミュレーション用のエネルギー関数の最適化方法に関する。
【背景技術】
【０００２】  分子シミュレーション技術には、ドッキングシミュレーション、分子力学法、分子軌道法などある。これらの分子シミュレーション技術は、効率的な薬剤探索を目的としたｉｎ  ｓｉｌｉｃｏ創薬スクリーニングに利用されている。【０００３】  ドッキングシミュレーションや分子力学法ではエネルギー関数を用いてシミュレーションが行われる。このエネルギー関数にはパラメータが含まれる。分子シミュレーションを実行する際にこのパラメータが不適当であると、実際の分子の挙動を正確に表現することはできない。【０００４】  そこでドッキングシミュレーションや分子力学法による分子シミュレーションを行う際には、まず始めにエネルギー関数に含まれるパラメータを最適化する。最適化は、実験等で得られた既知の分子構造を、エネルギー関数によって再現できるように行われる。このような最適化方法は経験的手法と呼ばれる。【０００５】  以下、エネルギー関数に含まれるパラメータを最適化することを、エネルギー関数の最適化と呼ぶ。また、経験的手法を用いて最適化されるパラメータのことを経験的パラメータと呼ぶ。【０００６】  エネルギー関数の最適化では、分子シミュレーションの考え方に応じて適切な指標が設けられる。【０００７】  分子シミュレーションのためのエネルギー関数の最適化指標として、大きく次の２つが挙げられる。第１に、分子の安定構造とその分子の不安定構造とにおけるエネルギーを用いることである。第２に、複数種類の分子の安定構造におけるエネルギーを用いることである。【０００８】  第１の指標である、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いてエネルギー関数を最適化することについて説明する。例えば、分子Ａと分子Ｂが結合して形成される分子ＡＢのシミュレーションを行う場合、分子Ａ上における分子Ｂの配置を考慮する必要がある。分子Ａと分子Ｂが結合する場合、分子Ｂが分子Ａに結合できる箇所は複数ある。そのため、考えられる分子ＡＢの構造は複数種類ある。分子ＡＢの安定性は、分子Ｂが分子Ａのどこに結合するかで決定されるため、複数種類考えられる分子ＡＢの構造には安定なものもあれば不安定なものをある。第１の指標を用いてエネルギー関数の最適化を行う場合、分子ＡＢが最も安定な構造であるときのエネルギーと、その他の不安定な構造であるときのエネルギーを利用する。すなわち、第１の指標を用いてエネルギー関数の最適化を行うとは、同一の原子又は分子からなる分子における、安定構造と不安定構造のエネルギーに基づき最適化を行うことである。【０００９】  分子の安定構造におけるエネルギーとは、その分子の結合エネルギーが最も低い状態にあるときのエネルギーを指す。逆に、分子の不安定構造におけるエネルギーとは、分子の安定構造におけるエネルギーよりも、エネルギーが高いものを指す。【００１０】  第１の指標の例として、分子の安定構造におけるエネルギーとその分子の不安定構造におけるエネルギーとの平均二乗誤差（ＲＭＳＤ）を用いた例が、Antes,  I.  et  al.  J.  Chem.  Inf.  Model.  2005,  45,  1291-1302.に示される。また、ドッキングシミュレーションではないが、第１の指標の別の例として、Fujitsuka,  Y.  et  al.  Proteins.  2004,  54,  88-103.に示さ���るｚスコアがある。本文献で定義されているｚスコアは、タンパク質フォールディングシミュレーション用エネルギー関数の最適化に用いられており、フォールディング後の分子の安定構造におけるエネルギーに対して、その分子の不安定構造におけるエネルギーの関係性を考慮するものである。【００１１】  第２の指標である、複数種類の分子の安定構造におけるエネルギーを用いてエネルギー関数を最適化することについて説明する。複数種類の分子とは、例えば、任意の分子Ａと、分子Ａに結合する分子Ｂ、分子Ｃ、分子Dがある場合、それぞれの結合状態である分子ＡＢ、分子ＡＣ、分子ＡＤを指す。さらに、この場合、前述した分子ＡＢ，分子ＡＣ、分子ＡＤはそれぞれにおいて安定な構造を形成している。第２の指標を用いて最適化するとは、安定構造であるこれらの分子ＡＢ、分子ＡＣ、分子ＡＤの各々が持つエネルギーを利用するということである。第２の指標を用いて最適化するとは、言い換えれば、異なる分子において、それぞれの分子が有する安定な構造におけるエネルギーを用いてエネルギー関数の最適化を行うことである。【００１２】  なお、結合構造である分子ＡＢ、分子ＡＣ、分子ＡＤが分子であれば良いため、分子Ａ、分子Ｂ、分子Ｃ、分子Dは原子でも良い。【００１３】  第２の指標を用いた例が、Bohm,  H.  J.  J.  Comput.-Aided  Mol.  Des.  1998,  12,  309-323に示される。本文献では、結合状態にある分子の安定構造におけるエネルギーの実験値と、その安定構造を予測した計算値との相関係数を用いている。【００１４】  分子シミュレーションの目的と、シミュレーションに用いるエネルギー関数の最適化手段は一致させる必要がある。したがって、最適化指標は第１の指標と第２の指標のうち、分子シミュレーションで行う予測に対応する方を選ぶ。分子の安定構造とその分子の不安定構造を扱う場合は、第１の指標を用いるのが望ましい。一方、複数種類の分子の安定構造を扱う分子シミュレーションにおいては、エネルギー関数の最適化に第２の指標を用いるのが望ましい。【００１５】  しかしながら、実際には、分子シミュレーションで行う予測に合致しない最適化指標を用いて最適化されたエネルギー関数を利用して、分子シミュレーションが行われている場合もある。例えば、Bissantz,  C.  et  al.  J.  Med.  Chem.  2000,  43,  4759-4767.では、分子シミュレーションの予測で扱うのは第１の指標である、分子の安定構造とその分子の不安定構造におけるエネルギーであるのにも関わらず、シミュレーションで用いるＦｌｅｘＸのエネルギー関数は第２の指標を用いて最適化されている。
【発明の開示】
【００１６】  上述のように、分子シミュレーションにおいて、実際行う予測と合致しない方法で最適化されたエネルギー関数を用いることは論理的に正しくなく、予測精度は著しく悪くなる。
【００１７】  また、第１の指標と第２の指標のいずれかを用いて最適化されたエネルギー関数には、それぞれの指標で用いたエネルギーのみしか考慮していないという問題がある。第１の指標を用いたエネルギー関数の最適化では、分子の安定構造とその分子の不安定構造とにおけるエネルギーしか考慮されていない。一方、第２の指標を用いたエネルギー関数の最適化では、複数種類の分子の安定構造におけるエネルギーしか考慮されていない。そのため、第１の指標と第２の指標の両方を考慮すべき子シミュレーションにおいて、いずれか一方のみで最適化されたエネルギー関数を用いると分子シミュレーションの予測精度が低くなってしまう。
【００１８】  以上のことから、分子シミュレーションでは、第１の指標で用いるエネルギーと第２の指標で用いるエネルギーの両方を考慮できることが望ましい。したがって、分子シミュレーションで用いるエネルギー関数の最適化では、上記２つのエネルギーを最適化指標として扱えることが望ましい。両方のエネルギーを考慮して最適化したエネルギー関数を用いることによって、創薬スクリーニングのような、分子の安定構造とその分子の不安定構造におけるエネルギーと、複数種類の分子の安定構造におけるエネルギーとの両方を扱うシミュレーションにて、高い予測精度が実現できる。さらに、結合後の分子の安定構造と不安定構造を扱う予測と、複数種類の分子の安定構造を扱う予測との両方の分子シミュレーションを扱えるようになる。
【００１９】  本発明の目的は、分子の安定構造とその分子の不安定構造とにおけるエネルギーと、複数種類の分子の安定構造におけるエネルギーと、の両方を考慮してエネルギー関数の最適化を行う、エネルギー関数最適化システム及び方法、並びにプログラムを提供することにある。加えて、このエネルギー関数最適化方法によって最適化されたエネルギー関数を利用する分子シミュレーションシステム及び方法、並びにプログラムを提供することにある。
【００２０】  本発明の一態様によるエネルギー関数最適化システムは、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索手段を有する。
【００２１】  また、本発明の一態様によるエネルギー関数最適化プログラムは、コンピュータに、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索処理を実行させるものである。
【００２２】  また、本発明の一態様によるエネルギー関数最適化方法は、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索する多変数最適化指標最適解探索ステップを有する。
【００２３】  加えて、本発明の一態様による分子シミュレーションシステムは、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーション手段を有する。
【００２４】  また、本発明の一態様による分子シミュレーションプログラムは、コンピュータに、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーション処理を実行させるものである。
【００２５】  また、本発明の一態様による分子シミュレーション方法は、分子の安定構造と前記分子の不安定構造とにおけるエネルギーを用いるエネルギー関数の第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる前記エネルギー関数の第２の指標と、の多変数関数である多変数最適化指標を最適化指標として最適化された前記エネルギー関数を用いる分子シミュレーションステップを有する。
【図面の簡単な説明】
【００２６】
【図１】第１の実施形態の構成を示すブロック図である。
【図２】第１の実施形態の動作を示すフローチャートである。
【図３】第２の実施形態の構成を示すブロック図である。
【図４】第２の実施形態の動作を示すフローチャートである。
【図５】第３の実施形態の構成を示すブロック図である。
【図６】第３の実施形態の動作を示すフローチャートである。
【図７】第４の実施形態の構成を示すブロック図である。
【図８】第４の実施形態の動作を示すフローチャートである。
【図９】実施例１で使用したタンパク質と、各タンパク質と結合する化合物とタンパク質との結合構造における結合エネルギーの一覧である。
【図１０】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて計算された、第１の指標の値、第２の指標の値、及び多変数最適化指標の値を示す表である。
【図１１】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１２】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒ１を使用）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１３】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒ２を使用）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１４】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒ３を使用）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１５】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒ４を使用）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１６】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒ５を使用）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。
【図１７】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒｎを使用、ｎ＝１，２，３，）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。なお、ｘ軸は２５％、ｙ軸は５０％までしか示していない。
【図１８】本発明の実施形態に係る最適化方法によって最適化されたエネルギー関数（第２の指標としてｒｎを使用、ｎ＝２，４、５）、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、のそれぞれを用いて行った分子シミュレーションの結果である。なお、ｘ軸は２５％、ｙ軸は５０％までしか示していない。
【発明を実施するための最良の形態】
【００２７】  以下に、本発明の実施形態に係るエネルギー関数最適化システム及び方法、並びにプログラムとエネルギー関数最適化方法を利用した分子シミュレーションシステム及び方法、並びにプログラムについて、図１乃至図１８を用いて説明する。【００２８】  ＜第１の実施形態＞  第１の実施形態におけるエネルギー関数最適化システム及び方法について説明する。【００２９】  本実施形態におけるエネルギー関数最適化システムを図１に示す。本システムは、キーボード等の入力装置１と、プログラム制御により動作する処理装置２と、ディスプレイ装置や印刷装置などの出力装置３と、を含む。【００３０】  処理装置２は、多変数最適化指標最適解探索部２０１を含む。多変数最適化指標最適解探索部２０１は、エネルギー関数の最適解を探索する。その際、多変数最適化指標を最適化指標として用いる。多変数最適化指標最適解探索部２０１は、分子の安定構造とその分子の不安定構造とにおけるエネルギーを用いる第１の指標と、複数種類の分子の安定構造におけるエネルギーを用いる第２の指標との多変数関数である多変数最適化指標を生成する。そして、生成した多変数最適化指標を評価し、それが最適解であるか否かを判断する。最適解を探索し終えたら，多変数最適化指標最適解探索部２０１は最適解を与えるパラメータを含むエネルギー関数を出力装置３に出力する。最適解を与えるパラメータを含むエネルギー関数とは、最適化されたエネルギー関数のことである。【００３１】  本システムは、第１の指標と第２の指標の多変数関数である多変数最適化指標を最適化指標として前記エネルギー関数の最適解を探索し、エネルギー関数の最適化を行う。したがって、本システムは第１と第２の指標の両方を考慮してエネルギー関数の最適化を行っている。【００３２】  このため、本システムによって最適化されたエネルギー関数は、第１及び第２の指標のいずれか一方のみで最適化されたエネルギー関数よりも汎用性が高く、優れている。すなわち、分子シミュレーションを実行するにあたって、第１及び第２の指標のいずれか一方のみで最適化されたエネルギー関数を使用するよりも、本システムによって最適化されたエネルギー関数を使用する方が実際の分子の挙動を正確に表現でき、予測精度の高いシミュレーションを行うことができる。【００３３】  本システムによって最適化されたエネルギー関数を用いて行う分子シミュレーションが高い予測精度で実際の分子の挙動を表現できる理由について説明する。第１の指標の関数をｆ（ｐ）、第２の指標の関数をｇ（ｐ）、多変数関数最適化指標であるｆ（ｐ）とｇ（ｐ）の多変数関数をｈ（ｆ，ｇ）とする。ｐはエネルギー関数に含まれるパラメータを表す。ここで、ｆ（ｐ）の最適解を与えるパラメータ｛ｐ１｝と、ｇ（ｐ）の最適解を与えるパラメータ｛ｐ２｝は、それぞれｈ（ｈ，ｇ）空間のある断面における局所的最適解である。対して、本システムによって求められたｈ（ｆ，ｇ）の最適解を与えるパラメータ｛ｐ３｝は、ｈ（ｈ，ｇ）空間における大域的最適解である。それゆえ、本システムによって得られるエネルギー関数は、第１と第２の指標両方を最大限に考慮して最適化されている。よって、分子シミュレーション��おいて、本システムによって最適化されたエネルギー関数を使用する方が、第１及び第２の指標のいずれか一方のみで最適化されたエネルギー関数を使用するよりも、シミュレーションの予測精度が高くなる。【００３４】  次に、各構成について説明する。【００３５】  多変数最適化指標最適解探索部２０１で用いる第１と第２の指標は、第１の実施形態に係るシステムに含まれておらず、本システムの外から入力される。例えば、第１及び第２の指標を生成する外部システムを本システムとは別に用意して、その外部システムに第１及び第２の指標を生成させる。多変数最適化指標最適解探索部２０１は外部システムから第１と第２の指標を受け取った後、計算等を行ったりハッシュテーブルを用いたりして、多変数最適化指標を生成する。或いは、ユーザーが入力装置１を利用して第１と第２の指標を本システムに入力しても良い。【００３６】  分子の安定構造とその分子の不安定構造とにおけるエネルギーを用いる第１の指標とは、分子の安定構造とその分子の不安定構造とにおけるエネルギーを、エネルギー関数の最適化の際に考慮するという意味である。同様に、複数種類の分子の安定構造におけるエネルギーを用いる第２の指標とは、複数種類の分子の安定構造におけるエネルギーを、エネルギー関数の最適化の際に考慮するという意味である。【００３７】  第１と第２の指標は、既知の分子情報を利用して求められる。第１の実施形態に係るシステムには既知の分子情報は含まれておらず、既知の分子情報を本システムの外から取得する。既知の分子情報は、第１と第２の指標の生成に用いることができれば良く、その情報の取得方法は特に限定されない。例えば、第３の実施形態にて説明するように、システム内に既知の分子情報を記憶する記憶装置を設けても良い。或いは、ユーザーが入力装置１を利用して分子構造に関する情報を本システムに入力しても良い。【００３８】  既知の分子情報とは、分子構造及びエネルギー値またはエネルギー由来の物理量等を指す。【００３９】  本システムで使用するエネルギー関数は、分子シミュレーション用であって経験的手法によって最適化されるものであれば良い。例えば、ドッキングシミュレーション用のソフトウェアであるＦｌｅｘＸに含まれるエネルギー関数を用いることができる。【００４０】  第１の指標は、分子の安定構造とその分子の不安定構造とにおけるエネルギーを考慮して求められ、最適化指標として利用できるものであればその求め方は特に限定しない。例えば、式１によってｚスコアを計算し、求められたｚスコアを第１の指標として用いることができる。【００４１】
【数１】
2009090796.tif 000003
式１中に含まれるＥｅｘｐは分子の安定構造におけるエネルギーを表す。＜Ｅｃａｌｃ＞は分子の不安定構造におけるエネルギーの平均値を表す。σｃａｌｃは分子の不安定構造におけるエネルギーの標準偏差を表す。
【００４２】  ｚスコアを計算する場合、分子の不安定構造におけるエネルギーの平均値及び標準偏差を計算するために、分子の不安定構造は複数あるのが望ましい。【００４３】  分子の安定構造とその分子の不安定構造が含まれるデータセットが複数ある場合は、各々の分子について上記ｚスコアを計算し、それらｚスコアの平均を第１の指標として用いるのが望ましい。【００４４】  第１の指標として、分子の安定構造と不安定構造におけるエネルギーの差、ＲＭＳＤ等を用いることも可能である。【００４５】  第２の指標は、複数種類の分子の安定構造におけるエネルギーを考慮して求められ、最適化指標として利用できるものであればその求め方は特に限定しない。例えば、結合状態である分子の結合エネルギーの実験値と、エネルギー関数によって求められる分子の結合構造の計算値との相関係数を第２の指標として用いることができる。別の例として、第１の指標の例として挙げたｚスコアがある。【００４６】  多変数最適化指標最適解探索部２０１は、第１と第２の指標の多変数関数である多変数最適化指標を求めて評価するものであれば良く、その関数の形は特に限定しない。ユーザーが任意に一方の最適化指標に重みをつけたりすることもできる。多変数最適化指標の例として、第１と第２の指標の２変数関数である、第１の指標と第２の指標の二乗との積が挙げられる。【００４７】  多変数最適化指標最適解探索部２０１は、多変数最適化指標を用いて、本実施形態で用いるエネルギー関数の最適解を探索するものであれば良く、その探索方法は特に限定しない。例えば、エネルギー関数の最適解を与えうる全パラメータに対して計算を行うことによって、エネルギー関数の最適解を求めることができる。【００４８】  また、最適解候補を与えるパラメータの量が膨大で最適解探索のための時間を抑えたい場合は、最適化アルゴリズムを用いて近似解を求めて、その近似解を最適解としても良い。最適化アルゴリズムは、シミュレーティド・アニーリング手法、遺伝的アルゴリズム、タブ・サーチ手法、シミュレーティド・エボリューション、確率的進化手法など、一般的なものであれば良く、用いる最適化アルゴリズムは特に限定しない。【００４９】  前述した方法を用いて最適解を探索する場合、必ずしも全工程を終えてから最適解を求める必要はない。最適解を探索する途中で探索を中断し、それまでに求めた解の中から最も良い解を最適解として用いても良い。【００５０】  エネルギー関数に含まれるパラメータは複数の場合もある。その場合は、全てのパラメータを考慮してエネルギー関数の最適化を行う。【００５１】  次に、図１及び図２を参照して本システムの動作について説明する。図２は本システムの動作を示す。【００５２】  入力装置１から実行の命令が与えられると、処理装置２は、多変数最適化指標最適解探索部２０１によって、第１の指標と第２の指標を用いて多変数最適化指標の生成し評価する（ステップ５０１）。【００５３】  そして、ステップ５０１の動作は、多変数最適化指標最適解探索部２０１によって、最適化が終了するまで繰り返される（ステップ５０２）。最適化の終了はユーザーが選んだ最適解探索方法によって異なる。エネルギー関数の最適解を与えうる全パラメータに対して計算を行いエネルギー関数の最適解を探索する場合は、全てのパラメータに対して計算し終えたときに最適化が終了する。最適化アルゴリズムを用いてエネルギー関数の近似解を求める場合は、ユーザーが指定した最適化終了条件を満たしたときに最適化が終了する。【００５４】  最適化が終了した後、多変数最適化指標最適解探索部２０１が最適化されたエネルギー関数を出力装置３に出力する（ステップ５０３）。【００５５】  以上、説明したとおり、本システムは、第１の指標と第２の指標の多変数関数である多変数最適化指標を最適化指標として、前記エネルギー関数の最適解を探索し、エネルギー関数の最適化を行う。本システムは、タンパク質等の高分子の分子シミュレーションにおいて特に有効である。タンパク質等の高分子では、その立体構造が非常に複雑である。そのため、第１と第２の指標のいずれかのみを用いて最適化されたエネルギー関数では、実際の分子の挙動を正確に再現するのは難しい。予測精度の高い分子シミュレーションを行うためには、第１及び第２の指標両方を考慮してエネルギー関数を最適化するのが望ましい。本システムは第１及び第２の指標の両方を用いてエネルギー関数を最適化することによって、信頼性の高い分子シミュレーションを実現する。【００５６】  ＜第２の実施形態＞  第２の実施形態におけるエネルギー関数最適化システム及び方法であって、処理装置内に第１の指標生成部と第２の指標生成部を含むものについて説明する。本実施形態は、第１の実施形態の応用であるため、第１の実施形態と同様な点については説明を省略する。【００５７】  本実施形態と第１の実施形態の異なる点について図３、図４を使って説明する。図３は本実施形態におけるエネルギー関数最適化システムを示す。図４は本システムの動作を示す。【００５８】  処理装置２は、第１の指標生成部２０２と、第２の指標生成部２０３と、多変数最適化指標最適解探索部２０１と、を含む。第１の指標生成部２０２はエネルギー関数の第１の指標を生成し、求められた第１の指標を多変数最適化指標最適解探索部２０１に出力する。第２の指標生成部２０３はエネルギー関数の第２の指標を生成し、求められた第２の指標を多変数最適化指標最適解探索部２０１に出力する。多変数最適化指標最適解探索部２０１は第１の指標と第２の指標を受け取り、多変数最適化指標を求め、評価する。なお、ここで言う第１の指標及び第２の指標とは、第１の実施形態で説明したものと同じである。【００５９】  本システムの処理装置２は、第１の指標生成部２０２と第２の指標生成部２０３とを含む。多変数最適化指標最適解探索部２０１は、第１の指標の値と第２の指標の値とを用いて多変数最適化指標を求めて評価する。第１の指標生成部２０２と第２の指標生成部２０３を本システムに含むことによって、第１の指標の値と第２の指標の値を本システム内で準備することができる。そのため、多変数最適化指標の評価及び多変数最適化指標最適解の探索を円滑に行うことができる。【００６０】  ただし、第２の実施形態のシステムは第１の実施形態のシステムと同様に既知の分子情報は含んでいないため、第１の指標生成部２０２と第２の指標生成部２０��は、第１及び第２の指標を生成するにあたって、既知の分子情報を本システムの外から取得する。或いは、ユーザーが入力装置１を利用して本システムに既知の分子情報を入力しても良い。【００６１】  次に、図３及び図４を参照して本システムの動作について説明する。【００６２】  入力装置１から実行の命令が与えられると、処理装置２は、第１の指標生成部２０２によって、一つのパラメータに対して第１の指標を生成する（ステップ５００Ａ）。そして、生成された値を多変数最適化指標最適解探索部２０１に出力する。【００６３】  第２の指標生成部２０３は、ステップ５００Ａで用いたパラメータと同じものに対して第２の指標を生成する（ステップ５００Ｂ）。そして、生成された値を多変数最適化指標最適解探索部２０１に出力する。【００６４】  次に、多変数最適化指標最適解探索部２０１が、ステップ５００Ａ、ステップ５００Ｂで得られた第１及び第２の指標の値を用いて多変数最適化指標を求め、評価する（ステップ５０１）。【００６５】  図４のステップ５００Ａからステップ５０１までの動作は、多変数最適化指標最適解探索部２０１によって、最適化が終了するまで繰り返される（ステップ５０２）。【００６６】  そして、最適化が終了した後、多変数最適化指標最適解探索部２０１が最適化されたエネルギー関数を出力装置３に出力する（ステップ５０３）。【００６７】  なお、上述の動作説明においては第１の指標を生成（ステップ５００Ａ）した後に第２の指標を生成（ステップ５００Ｂ）したが、生成する順序はどちらが先でも構わない。第２の指標を生成した後に第１の指標を生成しても良い。或いは、第１と第２の指標の生成を並列に実行しても良い。【００６８】  以上、説明したとおり、本システムの多変数最適化指標最適解探索部２０１は第１の指標生成部２０２と第２の指標生成部２０３とを含む。これによって、第１の指標の値と第２の指標の値を本システム内で準備することができ、多変数最適化指標の評価及び多変数最適化指標最適解の探索を円滑に行うことができる。【００６９】  ＜第３の実施形態＞  第３の実施形態におけるエネルギー関数最適化システム及び方法であって、記憶装置を含むものについて説明する。本実施形態は、第１及び第２の実施形態の応用であるため、第１及び第２の実施形態と同様な点については説明を省略する。【００７０】  本実施形態と第１及び第２の実施形態の異なる点について図５、図６を使って説明する。図５は本実施形態におけるエネルギー関数最適化システムを示す。図６は本システムの動作を示す。【００７１】  本システムは記憶装置４を含む。記憶装置４は、トレーニングデータ記憶部４０１と探索過程発生データ記憶部４０２を含む。トレーニングデータ記憶部４０１は、第１の指標生成部２０２及び第２の指標生成部２０３でトレーニングデータとして用いる既知の分子情報を記憶している。第１及び第２の指標生成部は各指標を生成する際にトレーニングデータ記憶部４０１から既知の分子情報を読み込む。探索過程発生データ記憶部４０２は、多変数最適化指標最適解探索部２０１が多変数最適化指標の最適解を探索する際に生じる一時的な情報を記憶する。多変数最適化指標最適解探索部２０１��必要に応じて探索過程発生データ記憶部４０２に記憶されている情報を読み取る。【００７２】  本システムは、トレーニングデータ記憶部４０１を含む。トレーニングデータ記憶部４０１に記憶されるトレーニングデータとは、経験的手法によりエネルギー関数を最適化する際に使用する情報である。システムにトレーニングデータ記憶部４０１を含むことによって、トレーニングデータをシステム外に準備したり、ユーザーが入力したりする必要がなくなる。したがって、本システムでは、第１の指標生成部２０２及び第２の指標生成部２０３による第１及び第２の指標の生成を円滑に行うことができる。【００７３】  本システムは、探索過程発生データ記憶部４０２を含む。探索過程発生データ記憶部４０２は、多変数最適化指標最適解探索部２０１で行う多変数最適化指標の最適解の探索時における一時的な情報を保存する。例えば最適化アルゴリズムを用いて近似解を求める場合、探索中において、あるパラメータを用いて生成された解が最適解であるか否かが判断できない。そこで、そのようなパラメータは最適解候補を与えるものとして記憶しておく必要がある。そのようなパラメータの情報量が少量であれば処理装置２内に記憶しても良いが、莫大になる場合には、探索過程発生データ記憶部４０２に保存するのが望ましい。【００７４】  上述の方法で近似解を求める場合、最適化が終了した後に探索過程発生データ記憶部４０２に記憶されているデータの中から最も良い解を与えるパラメータを抽出すれば良い。【００７５】  トレーニングデータ記憶部４０１にはトレーニングデータとして既知の分子情報が記憶される。これは、第１の指標及び第２の指標を求めるために利用される情報である。すなわち、分子の安定構造とその不安定構造とにおけるエネルギーや、複数種類の分子の安定構造におけるエネルギーの情報である。例えば、複数種類の分子の結合エネルギーの実験値やエネルギー関数によって求められる複数種類の分子の結合構造の計算値、或いはそれらを求めるための情報や、結合後の分子の安定構造におけるエネルギーや結合後の分子の不安定構造におけるエネルギー、或いはそれらを求めるための情報が挙げられる。【００７６】  探索過程発生データ記憶部４０２は本システムに含まれなくても良い。例えば、エネルギー関数の最適解を与えうる全パラメータに対して計算を行ってエネルギー関数の最適解を探索する際に、あるパラメータに対して求められた最適解が、その直前のパラメータに対して求められた最適解よりも良かった場合にのみ、記録してある最適解を与えるパラメータの情報を更新するようにする。全てのパラメータに対して前述のように計算することによって、探索過程において発生する一時的な情報を最小限に抑えることが可能となる。このような場合においては、探索過程発生データ記憶部４０２を本システムから除くことができる。【００７７】  次に、図５及び図６を参照して本実施形態の動作について説明する。【００７８】  入力装置１から実行の命令が与えられると、処理装置２は第１の指標生成部２０２により、トレーニングデータ記憶部４０１からトレーニングデータを読み込み（ステップ５００a）、エネルギー関数を用いて一つのパラメータ��対して第１の指標を生成する（ステップ５００Ａ）。そして、第１の指標生成部２０２は、生成した第１の指標を多変数最適化指標最適解探索部２０１に出力する。【００７９】  第２の指標生成部２０３は、トレーニングデータ記憶部４０１からトレーニングデータを読み込み（ステップ５００ｂ）、エネルギー関数を用いてステップ５００Ａと同じパラメータに対して第２の指標を生成する（ステップ５００Ｂ）。そして、第２の指標生成部２０３は、生成した第２の指標を多変数最適化指標最適解探索部２０１に出力する。【００８０】  ステップ５０１ａからステップ５０１までの動作は、多変数最適化指標最適解探索部２０１によって、最適化が終了するまで繰り返される（ステップ５０２）。その際、多変数最適化指標最適解探索部２０１は必要に応じて探索中に生じる一時的な情報を探索過程発生データ記憶部４０２に記憶する。【００８１】  以上、説明したとおり、本システムは記憶装置４を含む。記憶装置４は、トレーニングデータ記憶部４０１と探索過程発生データ記憶部４０２とを含む。本システムはトレーニングデータ記憶部４０１をシステム内に含むことによって、トレーニングデータをシステム外に準備したり、ユーザーが入力したりする必要がなくなる。したがって、本システムでは、第１の指標生成部２０２及び第２の指標生成部２０３による第１及び第２の指標の生成を円滑に行うことができる。また、本システムは探索過程発生データ記憶部４０２をシステム内に含むことによって、最適化アルゴリズムを用いて近似解を求める際に発生する大量の一時的な情報を保存することが可能となる。【００８２】  ＜第４の実施形態＞  第４の実施形態では、上述したエネルギー関数最適化方法によって最適化されたエネルギー関数を利用する分子シミュレーションシステム及び方法について説明する。本実施形態は、第１、第２、及び第３の実施形態の応用であるため、第１、第２、及び第３の実施形態と同様な点については説明を省略する。【００８３】  本実施形態と第１、第２、及び第３の実施形態の異なる点について図７、図８を使って説明する。図７は本実施形態に係るエネルギー関数最適化方法によって最適化されたエネルギー関数を利用する分子シミュレーションシステムを示す。図８は本システムの動作を示す。【００８４】  本システムの処理装置２は分子シミュレーション部２０４を含む。分子シミュレーション部２０４は、エネルギー関数最適化方法によって最適化されたエネルギー関数４０３を読み込む。そして、分子シミュレーション部２０４は、その最適化されたエネルギー関数を用いて分子シミュレーションを実行する。分子シミュレーション部２０４は分子シミュレーションを実行後、その結果を出力装置３に出力する。【００８５】  本システムは、エネルギー関数最適化方法によって最適化されたエネルギー関数を利用する分子シミュレーション部２０４を含む。これによって、本システムは、第１の指標のみで最適化されたエネルギー関数または第２の指標のみで最適化されたエネルギー関数のいずれか一方のみを用いる分子シミュレーションシステムよりも、実際の分子の挙動を正確に表現でき、予測精度の高いシミュレーションを行うことができる。【０���８６】  第４の実施形態に係るシステムにはエネルギー関数を最適化する機能が含まれていないため、最適化されたエネルギー関数４０３は、本システムの外で用意する。例えば、本システムとは別のシステムを用意し、その別のシステムにエネルギー関数の最適化を行わせても良い。【００８７】  最適化されたエネルギー関数４０３は、本システムで用いる分子シミュレーションのためのエネルギー関数であって、上述した方法によって最適化されたエネルギー関数であれば何でも良い。【００８８】  処理装置２には、第１の指標生成部、第２の指標生成部、多変数最適化指標最適解探索部、トレーニングデータ記憶部、探索過程発生データ記憶部が含まれていても良い。これらが全て本システム内に含まれていれば、エネルギー関数の最適化から分子シミュレーションまでの一連の動作を本システムのみで実行できる。このような場合においては、最適化されたエネルギー関数４０３を本システムの外で用意する必要はない。【００８９】  分子シミュレーション部２０４で用いる分子シミュレーション技術は、経験的手法を用いて最適化されるエネルギー関数を利用するものであれば特にその方法は限定しない。例として、ドッキングシミュレーション、分子力学法などがある。【００９０】  次に、図７及び図８を参照して本実施形態の動作について説明する。【００９１】  入力装置１から実行の命令が与えられると処理装置２は分子シミュレーション部２０４により、最適化されたエネルギー関数４０３を読み込む（ステップ５０６）。そして、その最適化されたエネルギー関数を用いて分子シミュレーションを実行する（ステップ５０７）。分子シミュレーション部２０４は分子シミュレーションの実行後、その結果を出力装置３に出力する（ステップ５０８）。【００９２】  以上、説明したとおり、本システムは、上述したエネルギー関数最適化方法によって最適化されたエネルギー関数を利用する分子シミュレーション部を含んでいる。これによって、第１の指標のみで最適化されたエネルギー関数及び第２の指標のみで最適化されたエネルギー関数のいずれか一方のみを用いる分子シミュレーションシステムよりも、予測精度が高いシミュレーションを行うことができる。【００９３】  ＜実施例１＞  次に、本発明の実施例を、図９，１０を参照して説明する。【００９４】  本実施例では、エネルギー関数の最適化を、本発明の実施形態に係る最適化方法、第１の指標のみを用いる方法、第２の指標のみを用いる方法の３通りで行い、その結果を比較した。【００９５】  本実施例で用いる分子シミュレーションシステムは、第３の実施形態における分子シミュレーションシステムと同様の構成を成している。【００９６】  記憶装置４に含まれるトレーニングセット記憶部４１には、エネルギー関数の最適化に用いる３８種類のタンパク質の情報が記憶されている。それぞれのタンパク質には、そのタンパク質と結合することが既知の化合物（リガンド）が存在する。タンパク質とリガンドが結合することによって形成される複合体は、タンパク質‐リガンド複合体と呼ばれる。それぞれのタンパク質の情報には、そのタンパク質‐リガンド複合体のＸ線結晶構造（結合構造）、結合構造生成シミュレーシ���ンによって人工的に生成された非結晶構造、及び結合エネルギーを求めるための実験データが含まれる。タンパク質‐リガンド複合体の非結晶構造の数は各タンパク質によって異なり、その数は６から９２である。【００９７】  分子シミュレーション用エネルギー関数は、ＦｌｅｘＸに含まれるエネルギー関数の一部を用いた。ＦｌｅｘＸのエネルギー関数は以下のように表される（Rarey,  M.ら  J.  Mol.  Biol.  1996,  261,  470-89.）。【００９８】
【数２】
2009090796.tif 000004
ΔＧｂｉｎｄはタンパク質‐リガンド複合体の結合状態を表す結合構造スコア、Ｆｉは位置に依存する関数を、ΔＧｉはエネルギーパラメータを、Σは相互作用に関わる全ての原子ペアの和を表す。ｍａｔｃｈは水素結合と、金属コンタクトと、芳香族間の相互作用から成る項である。ｌｉｐｏは疎水性相互作用を表す項、ａｍｂｉｇは極性原子と非極性原子の相互作用を表す項、ｃｌａｓｈは原子の衝突に対するペナルティ項、ｒｏｔは化合物がタンパク質と結合することによって失うエントロピー項を表す。ｎｒｏｔは化合物の回転可能単結合数である。ΔＧ０はこのエネルギー関数特有の定数項である。
【００９９】  本実施例で用いたエネルギー関数を以下に示す。【０１００】
【数３】
2009090796.tif 000005
式３は、式２からｃｌａｓｈとｒｏｔの項を除いたものである。これらの項は、実際の分子間の相互作用エネルギーに直接関わる項ではなく、分子シミュレーションの予測精度を上げるために付加された項である。本実施例では、式２からｃｌａｓｈとｒｏｔの項を除いた、相互作用エネルギーに関わる項のみで形成されるエネルギー関数を用いた。
【０１０１】  本実施例では、エネルギー関数の最適なパラメータを生成するために以下のスコアリング関数を定義した。すなわち、エネルギー関数の最適化では、このスコアリング関数の最適化を行った。【０１０２】
【数４】
2009090796.tif 000006
式４中の（ａ、ｂ、ｃ）は（ΔＧｍａｔｃｈ、ΔＧｌｉｐｏ、ΔＧａｍｂｉｇ）のスケーリングパラメータである。
【０１０３】  最適化の際には、各スケーリングパラメータを−３から３まで０．１毎に区切り、全ての組み合わせについて計算した。この場合、（ａ、ｂ、ｃ）の全ての組み合わせは６１３通りある。【０１０４】  第１の指標と第２の指標の生成では、タンパク質‐リガンド複合体の結合エネルギー、実験結合構造スコア、非結合構造スコアを用いた。【０１０５】  タンパク質‐リガンド複合体の結合エネルギーとして、結合エネルギーを求めるための実験データの平衡定数から求められる−ｌｏｇＫｄを用いた。【０１０６】  実験結合構造スコアとして、タンパク質‐リガンド複合体のＸ線結晶構造を用いてエネルギー関数を解くことによって得られる値を用いた。【０１０７】  非結合構造スコアとして、結合構造生成シミュレーションによって人工的に生成された非結晶構造を用いてエネルギー関数を解くことによって得られる値を用いた。【０１０８】  エネルギー関数の最適化に用いる３８個のタンパク質を図９に示す。タンパク質の情報はＰｒｏｔｅｉｎ  Ｄａｔａ  Ｂａｎｋ  （http://www.wwpdb.org/）に登録されているものを利用した。図９中に記載されたＰＤＢとはＰｒｏｔｅｉｎ  Ｄａｔａ  ＢａｎｋのＩＤを表す。図９中に記載された‐ｌｏｇＫｄは結合エネルギーを表す。【０１０９】  第１の指標は、タンパク質‐リガンド複合体の平均ｚスコア(ｚａｖｅ）とした。ｚａｖｅを求めるために、まずは式１に従ってタンパク質‐リガンド複合体におけるｚスコアを計算した。その際、Ｅｅｘｐとして実験結合構造スコアを用いた。＜Ｅｃａｌｃ＞として非結合構造スコアの平均値を用いた。σｃａｌｃとして非結合構造スコアの標準偏差を用いた。そして、それぞれのタンパク質‐リガンド複合体において得られたｚスコアの平均を求め、ｚａｖｅを求めた。【０１１０】  この場合、実験によって得られるタンパク質‐リガンド複合体の構造を正確に再現するためには、ｚａｖｅは負に大きい方が良い。したがって、エネルギー関数の最適化はｚａｖｅが負に大きくなるように実行された。【０１１１】  第２の指標として、複数種類の分子の結合エネルギーの実験値と、エネルギー関数によって求められる複数種類の分子の結合構造の計算値との相関係数を用いた。複数種類の分子の結合エネルギーの実験値として、タンパク質‐リガンド複合体の結合エネルギーを用いた。エネルギー関数によって求められる複数種類の分子の結合構造の計算値として、実験結合構造スコアを用いた。結合エネルギーと実験結合構造スコアのデータ列をそれぞれｘ＝｛ｘｉ｝、ｙ＝｛ｙｉ｝、（ｉ＝１，２，…，３８）とすると、相関係数ｒは次の式で表される。【０１１２】
【数５】
2009090796.tif 000007
ｘａｖｅとｙａｖｅは各データ列の相加平均を表す。
【０１１３】  この場合、実験によって得られるタンパク質‐リガンド複合体の構造を正確に再現するためには、相関係数は高い方が良い。したがって、エネルギー関数の最適化は相関係数が高くなるように実行された。【０１１４】  多変数最適化指標（ｏｐｔ）は、次の式で表されるものを用いた。【０１１５】
【数６】
2009090796.tif 000008
上述したｒとｚａｖｅの特徴から、実験によって得られるタンパク質‐リガンド複合体の構造を正確に再現するためには、ｏｐｔは負に大きい方が良い。したがって、エネルギー関数の最適化は、ｏｐｔを負に最大化させるパラメータを探索するように実行された。
【０１１６】  本実施例では、本発明の実施形態に係る最適化指標としてｏｐｔを用いてエネルギー関数の最適化を行った。そして、比較例１として第１の指標であるｚａｖｅのみを用いてエネルギー関数の最適化を行った。加えて、比較例２として第２の指標であるｒのみを用いてエネルギー関数の最適化を行った。【０１１７】  説明の便宜上、ｏｐｔを用いて最適化されたエネルギー関数を、本発明の実施形態に係るエネルギー関数と呼ぶ。第１の指標のみを用いて最適化されたエネルギー関数を第１のエネルギー関数と呼ぶ。第２の指標のみを用いて最適化されたエネルギー関数を第２のエネルギー関数と呼ぶ。【０１１８】  それぞれの最適化の結果について、図１０を用いて説明する。図１０は、それぞれの最適化指標を用いて最適化されたエネルギー関数を用いて、ｚａｖｅ、ｒ、ｏｐｔを求めた結果を示す。図中の第１の指標の値はｚａｖｅを表す。図中の第２の指標の値はrを表す。図中の多変数最適化指標の値はｏｐｔを表す。図中の「本実施形態」は、本発明の実施形態に係るエネルギー関数を表す。比較例１は、第１のエネルギー関数を表す。比較例２は第２のエネルギー関数を表す。【０１１９】  比較例１では、他のものと比べて、第１の指標の値が最も良かった。しかし、この場合、最適化の際に第２の指標については一切考慮されていないため、第２の指標の値は悪かった。一方、比較例２では第２の指標の値が最も良かった。しかし、この場合、最適化の際に第１の指標については一切考慮されていないため、第１の指標の値は悪かった。【０１２０】  本発明の実施形態に係るエネルギー関数では、比較例１と比較例２に比べて、第１の指標の値及び第２の指標の値共に良い結果を得られた。本発明の実施形態に係るエネルギー関数の第２の指標の値においては、比較例２における第２の指標の値と同等であった。多変数最適化指標の値は本発明の実施形態に係るエネルギー関数によって求めたものが最も良かった。【０１２１】  第１または第２の指標のみを用いて最適化した場合、他の指標は全く考慮されていない。対して、本発明の実施形態では第１の指標と第２の指標の多変数関数を最適化指標として用いたため、両指標を考慮して最適化が行われている。そのため、本発明の実施形態に係るエネルギー関数では両指標において良い結果を得ることができた。【０１２２】  したがって、分子シミュレーションにおいて本発明に係るエネルギー関数を用いれば、第１の指標と第２の指標との両方を考慮してシミュレーションを実行することができる。【０１２３】  以上の結果から、本発明の実施形態によって、エネルギー関数の第１の指標である、分子の安定構造とその分子の不安定構造とにおけるエネルギーと、エネルギー関数の第２の指標である、複数種類の分子の安定構造におけるエネルギーと、の両方を考慮してエネルギー関数の最適化を行うことが可能となる。本発明の実施形態によって最適化されたエネルギー関数は前述した二つの要素を考慮して最適化されているため、それぞれを単独で用いて最適化されたエネルギー関数よりも、分子シミュレーションに適している。すなわち、本システムによって最適化されたエネルギー関数を用いることによって、特に、両方のエネル���ーを考慮する必要のある分子シミュレーションにおいて、高い予測精度が実現できる。【０１２４】  ＜実施例２＞  本実施例では、実施例１にて最適化された本発明の実施形態に係るエネルギー関数の精度を評価するために以下のドッキングシミュレーションのテストを行った。【０１２５】  テストに用いる標的タンパク質としてＥｓｔｒｏｇｅｎ  Ｒｅｃｅｐｔｏｒ（ＥＲ）を使用した。また、ＥＲにドッキングさせる化合物として、結合することが分かっている化合物（既知リガンド）３２個、結合しないことが分かっている化合物１３４０個の計１３７３個を用いた。【０１２６】  化合物の多数の配座を発生させるプログラムとして、市販のドッキングシミュレーションソフトウェアＦｌｅｘＳＩＳを用いた。【０１２７】  評価方法を以下に説明する。【０１２８】  はじめにＦｌｅｘＳＩＳを用いて全ての化合物に対して配座を発生させた。【０１２９】  次に、実施例１において最適化された本発明の実施形態に係るエネルギー関数を用いて、ある化合物の全配座に対してＥＲとの結合エネルギーを計算した。そして、結合エネルギーが最も低い配座を選抜し、その結合エネルギーをＥＲと化合物の結合エネルギーとした。用意した１３７３個の全化合物に対して同様の計算をし、全化合物に対して結合エネルギーを決定した。その後、結合エネルギーが低い順に化合物を並べ、化合物に対して昇順に順位付けを行った。【０１３０】  結合時の分子の構造においては結合エネルギーが低い方が安定であると考えられる。したがって、上述のように化合物に対して順位付けを行った場合、上位の化合物ほどＥＲとの結合状態が安定であるといえる。すなわち、上位の化合物ほど、シミュレーションによる予測上ではＥＲと結合しやすいと考えられる。【０１３１】  さらに、第１のエネルギー関数、第２のエネルギー関数を用いて、同様の計算及び動作を行った。【０１３２】  結果を図１１に示す。図１１はエンリッチメントプロットである。ｘ軸は、順位付けされた全化合物の上位からの割合（％）を表す。ｙ軸は、全既知リガンドに対する、全化合物の上位ｘ％に含まれる既知リガンドの割合（％）を表す。図中の「本実施形態」は本発明の実施形態に係るエネルギー関数の結果を表す。図中の比較例１は第１のエネルギー関数の結果を表す。図中の比較例２は第２のエネルギー関数の結果を表す。【０１３３】  本実施例の場合、ｘの値が低いところでｙの値が高くなるような結果を示すエネルギー関数が最も優れていると評価できる。【０１３４】  本発明の実施形態に係るエネルギー関数、比較例１、比較例２のそれぞれの結果を比較したところ、図１１に示すとおり、本発明の実施形態に係るエネルギー関数の結果が最も良かった。すなわち、ＥＲと化合物の結合状態におけるエネルギーを予測して、ＥＲと結合しやすい化合物の順位付けをそれぞれのエネルギー関数を用いて行ったところ、最も多くの既知リガンドが本発明の実施形態に係るエネルギー関数による順位付けの上位に含まれていた。したがって、本発明の実施形態に係るエネルギー関数を用いて行った分子シミュレーションの予測精度が最も高かった。【０１３５】  以上、説明したとおり、第１の指標のみを用いて最適化されたエネルギー関数、第２の指標のみを用いて最適化されたエネルギー関数、本発明の実施形態に係る最適化方法を用いて最適化されたエネルギー関数のそれぞれを用いて分子シミュレーションを行ったところ、本発明の実施形態に係る最適化方法を用いて最適化されたエネルギー関数を用いた分子シミュレーションの予測精度が最も高かった。本発明の実施形態に係る最適化方法を用いて最適化されたエネルギー関数を分子シミュレーションに用いることによって、実際の分子の挙動をより正確に表現することが可能となる。【０１３６】  ＜実施例３＞  実施例３では、下記式によって求められるｏｐｔを用いて最適化されたエネルギー関数の精度を評価した。下記式に含まれるｒの指数ｎは任意の数である。本実施例では指数ｎを変化させることによるエネルギー関数の性能の変化を調べた。なお、ｎ＝２の場合における下記式は実施例１で用いたｏｐｔを求める式（式６）と同じものとなる。【０１３７】
【数７】
2009090796.tif 000009
初めにエネルギー関数の最適化を行った。エネルギー関数は実施例１で用いたものと同じものを用いた。エネルギー関数の最適化は実施例１と同様に行った。ただし、ｏｐｔとして、ｎ＝１、２，３，４，５の場合における上記式を用いた。最適化のために用いるタンパク質データとして、実施例１で用いた３８個のタンパク質データからｌａｂｅ、３ｆｘ２、４ｃｌａを除いた計３５個のデータを用いた。
【０１３８】  次に、最適化したエネルギー関数の精度評価を実施例２と同様の方法にて行い、比較例１，２と本発明の実施形態によって最適化されたエネルギー関数の結果とを比較した。【０１３９】  精度評価の結果を図１２〜図１８に示す。図１２〜図１６はｎ＝１，２、３，４，５の場合の結果を表すエンリッチプロットメントである。図１７はｎ＝１，２，３のそれぞれの場合における順位付けされた全化合物の上位２５％の結果を表すエンリッチプロットメントである。従って、図１７は図１２、１３，１４に示す順位付けされた全化合物の上位２５％を抽出した図である。図１８はｎ＝２，４、５のそれぞれの場合における順位付けされた全化合物の上位２５％の結果を表すエンリッチプロットメントである。従って、図１８は図１３、１５、１６に示す順位付けされた全化合物の上位２５％を抽出した図である。ｘ軸は順位付けされた全化合物の上位からの割合（％）を表す。ｙ軸は、全既知リガンドに対する、全化合物の上位ｘ％含まれる既知リガンドの割合（％）を表す。図中のｚｒ１は、本実施形態のエネルギー関数であってｎ＝１の場合におけるｏｐｔを用いて最適化されたものの結果である。ｚｒ２〜ｚｒ５もｚｒ１と同様に、本実施形態のエネルギー関数であって、ｎ＝２〜５の場合におけるｏｐｔを用いて最適化されたものの結果である。比較例１は第１のエネルギー関数の結果、比較例２は第２のエネルギー関数の結果を表す。【０１４０】  図１７、１８では、ｘ軸は２５％、ｙ軸は５０％までしか示していない。これは、薬剤開発の初期工程である化合物スクリーニング等における実用性を考慮すると、上位２０％程度までの結果が性能評価の対象として相応しいからである。【０１４１】  本発明の実施形態に係るエネルギー関数であってｎ＝１〜５のもの、比較例１、比較例２のそれぞれの結果を比較したところ、図１２〜図１６に示すとおり、本発明の実施形態に係るエネルギー関数であってｎ＝１〜５のものの結果が最も良かった。すなわち、ＥＲと化合物の結合状態におけるエネルギー関数を予測して、ＥＲと結合しやすい化合物の順位付けをそれぞれのエネルギー関数を用いて行ったところ、最も多くの既知リガンドが本発明の実施形態に係るエネルギー関数による順位付けの上位に含まれていた。【０１４２】  本実施例では、ｏｐｔの計算式（式７）に含まれるｒｎの指数ｎを１〜５として、エネルギー関数の最適化を行った。そして、ｎ＝１〜５の全ての場合におけるエネルギー関数にて、比較例よりも優れた結果を得ることができた。したがって、本実施形態のエネルギー関数の第２の指標として、ｒｎのような単調関数を用いることは有用であると考えられる。【０１４３】  なお、本実施例におけるｒｎはnが大きくなるとｒｎは小さくなることから、ｒｎは単調減少関数であると言える。【０１４４】  加えて、ｎ＝１〜５の場合におけるそれぞれの結果を比較したところ、図１７、１８に示すとおり、ｎ＝２の場合における結果が最も良かった。【０１４５】  なお、本実施例では最適化指標としてｎ＝１〜５におけるｏｐｔ（式７）を用いたが、ｎの値はこれらに限定されるものではない。例えば、ｎ＝６やｎ＝７におけるｏｐｔを用いても良い。【０１４６】  また��最適なｏｐｔは、エネルギー関数や用いるデータ等によって変わる。従って、ユーザは適宜、ｎの値や第１の指標、第２の指標等を変更して、最適なｏｐｔを設定するのが望ましい。【０１４７】  実施例３では、本発明の実施形態に係るエネルギー関数の最適化方法にて、単調減少関数ｒｎを第２の指標として用いた。そして、それによって最適化されたエネルギー関数の性能と比較例の性能とを比べて、本発明の実施形態に係る方法に従って最適化されたエネルギー関数の方が優れていることを確認した。従って、本実施形態のエネルギー関数の最適化方法で第２の指標として単調減少関数を用いることによって、高い性能を有するエネルギー関数が得られるようになる。【０１４８】  以上、実施形態および実施例を参照して本発明を説明したが、本発明は実施形態および実施例に限定されるものではない。クレームに定義された本発明の構成や詳細には、発明のスコープ内で当業者が理解し得る様々な変更をすることができる。【０１４９】  この出願は、２００８年１月１６日に出願された日本出願特願２００８−００６４２６と、２００８年４月２２日に出願された日本出願特願２００８−１１０８７１とを基礎として優先権の利益を主張するものであり、その開示の全てを引用によってここに取り込む。
【図１】
2009090796.tif 000013
【図２】
2009090796.tif 000014
【図３】
2009090796.tif 000015
【図４】
2009090796.tif 000016
【図５】
2009090796.tif 000017
【図６】
2009090796.tif 000018
【図７】
2009090796.tif 000019
【図８】
2009090796.tif 000020
【図９】
2009090796.tif 000021
【図１０】
2009090796.tif 000022
【図１１】
2009090796.tif 000023
【図１２】
2009090796.tif 000024
【図１３】
2009090796.tif 000025
【図１４】
2009090796.tif 000026
【図１５】
2009090796.tif 000027
【図１６】
2009090796.tif 000028
【図１７】
2009090796.tif 000029
【図１８】
2009090796.tif 000030
【国際調査報告】
2009090796.tif 000031
2009090796.tif 000032
2009090796.tif 000033
2009090796.tif 000034