(57)【特許請求の範囲】
【請求項１】
基準状態と標的状態との間の自由エネルギー差を計算するためにコンピュータで実施される方法において、前記基準状態及び前記標的状態はそれぞれ、原子ＰＡＢの共通の組を含み、前記基準状態は原子ＰＡの組を更に含み、前記標的状態は原子ＰＢの組を更に含み、前記組ＰＡは、前記基準状態のみに存在し、且つ前記標的状態には存在せず、及び前記組ＰＢは、前記標的状態のみに存在し、且つ前記基準状態には存在せず、少なくとも２つの原子Ａａ及びＡｂが存在し、Ａａ及びＡｂは、（１）前記基準状態では互いに原子価結合せず、且つ前記標的状態では原子価結合するか、又は（２）基準状態では互いに原子価結合し、且つ前記標的状態では互いに原子価結合しないかのいずれかであり、前記方法は、（ａ）ＰＡ、ＰＢ、及びＰＡＢ内の全ての前記原子について、前記原子間の結合接続及び前記原子の相対空間構成を含むトポロジを提供する手順と、
（ｂ）前記基準状態と前記標的状態との間の変換路に沿った１つ又は複数の中間系状態を特定することであって、前記変換路は、各系状態の原子間相互作用から生じるエネルギーを変調する結合パラメータλによって定義され、前記結合パラメータλは、複数の成分であって、それぞれが［０，１］に属する値を有する複数の成分を含み、且つ異なるタイプの相互作用エネルギーを変調する、前記基準状態と前記標的状態との間の変換路に沿った１つ又は複数の中間系状態を特定する手順と、
（ｃ）前記基準状態、前記標的状態、及び中間系状態のマイクロ状態のアンサンブルを取得するために、少なくとも１つのコンピュータプロセッサを使用して、分子シミュレーションを実行する手順であって、前記系状態のそれぞれに対して分子シミュレーションを実行する手順は、前記原子Ａａ及びＡｂ間の結合伸縮相互作用エネルギーを計算することを含み、前記結合伸縮相互作用エネルギーは、弱結合ポテンシャルによって定義され、
前記弱結合ポテンシャルは、前記結合パラメータλの結合伸縮成分λｓｂｓの関数であり、且つ［０，１］内のλｓｂｓの全ての値及びＡａとＡｂとの間の距離ｒの全ての値について、いかなる特異領域も含まず、前記弱結合ポテンシャルは以下の条件：λｓｂｓが（０，１）内にある場合、ＡａとＡｂとの間の前記距離が無限に近づくとき、前記弱結合ポテンシャルは平坦であり、Ａａ及びＡｂが前記基準状態又は前記標的状態のいずれかで原子価結合されないとき、前記弱結合ポテンシャルは、ＡａとＡｂとの間の全ての距離について平坦且つゼロであり、並びにＡａ及びＡｂが前記標的状態又は前記基準状態のいずれかで原子価結合されるとき、前記弱結合ポテンシャルは調和ポテンシャルに戻ることを更に満たす、分子シミュレーションを実行する手順と、
（ｄ）少なくとも１つのコンピュータプロセッサを使用して、前記標的状態、前記基準状態、及び前記中間系状態において得られるマイクロ状態の前記アンサンブルの解析により、前記基準状態と前記標的状態との間の前記自由エネルギー差を計算する手順と
を含んでなる、方法。
【請求項２】
λｓｂｓに対する前記弱結合ポテンシャルの一次導関数は、０〜１のλｓｂｓの全ての値及びＡａとＡｂとの間の前記距離ｒの全ての値について連続であり、且つ結合される、請求項１に記載の方法。
【請求項３】
ＡａとＡｂとの間の前記距離ｒに対する前記弱結合ポテンシャルの一次導関数及び二次導関数の両方が、λｓｂｓの全ての値について連続であり、且つ結合され、及びｒが無限に近づくとき、ゼロに近づく、請求項１に記載の方法。
【請求項４】
前記弱結合ポテンシャルは、（ｒ−ｒ０）２の関数であり、ここで、ｒ０はＡａとＡｂとの間の平衡距離である、請求項１に記載の方法。
【請求項５】
前記弱結合ポテンシャルは、
【数１】
0006298167.tif 000077
によって表され、式中、ｋは定数であり、並びに関数ｆ、ｇ、及びαはそれぞれ連続関数であり、且つ以下の条件：
ｆ（λｓｂｓ＝０）＝０、
ｆ（λｓｂｓ＝１）＝１、
ｇ（λｓｂｓ＝０）＝１、
ｇ（λｓｂｓ＝１）＝０、
ａ（ｋ，λｓｂｓ＜１）＞０
を満たす、請求項４に記載の方法。
【請求項６】
ｆ（λｓｂｓ）＝λｓｂｓ、ｇ（λｓｂｓ）＝１−λｓｂｓ、及びα（ｋ，λｓｂｓ）＝ｃｏｎｓｔである、請求項５に記載の方法。
【請求項７】
前記系状態のそれぞれに対して分子シミュレーションを実行する手順は、
Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において原子価結合しない場合、前記中間系状態のそれぞれでのＡａとＡｂとの間の前記結合伸縮相互作用エネルギーを計算するために、前記結合伸縮成分の第１のスケジュールλｓｂｓＡと対応する弱結合ポテンシャルとを使用することであって、λｓｂｓＡは、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、使用することと、
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において原子価結合する場合、ＡａとＡｂとの間の結合伸縮相互作用を計算するために、前記中間系状態のそれぞれ及びλｓｂｓＢに対応する弱結合ポテンシャルについて、前記結合伸縮成分の第２のスケジュールλｓｂｓＢと対応する弱結合ポテンシャルとを使用することであって、λｓｂｓＢは、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では０から１まで変化する、使用することと
を含んでなる、請求項１に記載の方法。
【請求項８】
前記系状態のそれぞれに対して分子シミュレーションを実行する手順は、
（ａ）適切なパラメータを力場の結合角相互作用に使用して、（ｉ）Ａａ及び別の原子Ａｃによって形成される結合と、（ｉｉ）前記基準状態から前記標的状態への変換により破断されているか、又は形成されているＡａとＡｂとの結合との間の結合角相互作用を計算することと、
（ｂ）Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において互いに原子価結合しない場合、（ａ）において得られる前記結合角相互作用を結合角結合パラメータλｂａＡで乗算することであって、λｂａＡは、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において互いに原子価結合する場合、（ａ）において得られる前記結合角相互作用を結合角結合パラメータλｂａＢで乗算することであって、λｂａＢは、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では０から１まで変化する、乗算することと、
（ｃ）（ｂ）において得られる前記結合角相互作用を前記対応する系状態のシミュレーションステップの総エネルギーに含めることと
を更に含んでなる、請求項７に記載の方法。
【請求項９】
前記系状態のそれぞれに対して分子シミュレーションを実行する手順は、
（ａ）力場の二面角相互作用に適切なパラメータを使用して、４つの接続原子群｛Ａｉ，Ａｊ，Ａｋ，Ａｌ｝の二面角相互作用を計算することであって、前記群はＡａ及びＡｂの両方を含む、計算することと、
（ｂ）Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において互いに原子価結合しない場合、（ａ）において得られる前記二面角相互作用を二面角結合パラメータλｂｄＡで乗算することであって、λｂｄＡは、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において互いに原子価結合する場合、（ａ）において得られる前記二面角相互作用を二面角結合パラメータλｂｄＢで乗算することであって、λｂｄＢは、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では０から１まで変化する、乗算することと、
（ｃ）（ｂ）において得られる前記二面角相互作用を前記対応する系状態の前記シミュレーションステップの前記総エネルギーに含めることと
を更に含んでなる、請求項８に記載の方法。
【請求項１０】
Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において互いに原子価結合しない場合、結合角相互作用結合パラメータλｂａＡ及び結合二面角相互作用結合パラメータλｂｄＡはそれぞれ、λｓｂｓＡが予め定義される閾値よりも小さい場合、０であるように選択され、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において互いに原子価結合する場合、前記結合角相互作用結合パラメータλｂａＢ及び前記結合二面角相互作用結合パラメータλｂｄＢはそれぞれ、λｓｂｓＢが予め定義される閾値よりも小さい場合、０であるように選択される、請求項９に記載の方法。
【請求項１１】
前記状態の全てに対して分子シミュレーションを実行する手順は、
（ａ）力場の静電相互作用及びファンデルワールスに適切なパラメータを使用して、２つの原子Ａｉ及びＡｊ間の非結合静電相互作用及びファンデルワールス相互作用を計算することであって、対（Ａｉ，Ａｊ）の非結合排斥ステータスは、前記基準状態から前記標的状態への変換による影響を受ける、計算することと、
（ｂ）Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において原子価結合せず、及びＡｉとＡｊとの間の前記非結合静電相互作用が、前記基準状態では排除されるが、前記標的状態では排除されない場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電相互作用及びファンデルワールス相互作用を結合パラメータλｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘでそれぞれ乗算することであって、λｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘの両方が、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では０から１まで変化する、乗算すること、
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において原子価結合し、及びＡｉとＡｊとの間の前記非結合静電相互作用が、前記基準状態では排除されないが、前記標的状態では排除される場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電相互作用及びファンデルワールス相互作用を結合パラメータλｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘでそれぞれ乗算することであって、λｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘの両方が、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算することと、
（ｃ）（ｂ）において得られる前記非結合静電相互作用及びファンデルワールス相互作用を前記対応する系状態の前記シミュレーションステップの前記総エネルギーに含めることと
を更に含んでなる、請求項９に記載の方法。
【請求項１２】
前記状態の全てに対して分子シミュレーションを実行する手順は、
Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において原子価結合せず、及びＡｉとＡｊとの間の前記非結合静電相互作用が、前記基準状態では排除されるが、前記標的状態では排除されない場合、λｖｄｗＡｅｘが中間系状態で１よりも小さいとき、λｅｌｅｃＡｅｘが前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従って、λｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘのそれぞれを変更し、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において原子価結合し、及びＡｉとＡｊとの間の前記非結合静電相互作用が、前記基準状態では排除されないが、前記標的状態では排除される場合、λｖｄｗＡｅｘが中間系状態で１よりも小さいとき、λｅｌｅｃＡｅｘが前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従って、λｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘのそれぞれを変更すること
を更に含む、請求項１１に記載の方法。
【請求項１３】
前記状態の全てに対して分子シミュレーションを実行する手順は、
（ａ）力場の静電１−４対相互作用及びファンデルワールス１−４対相互作用に適切なパラメータを使用して、別の２つの介在する原子と一緒に、前記基準状態又は前記標的状態のいずれかで結合二面角相互作用を形成する２つの原子Ａｉ及びＡｊ間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を計算することと、
（ｂ）Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において原子価結合しない場合、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれるが、前記標的状態では含まれない場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＡ１４及びλｖｄｗＡ１４でそれぞれ乗算することであって、λｅｌｅｃＡ１４及びλｖｄｗＡ１４の両方が、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算すること、及び
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれないが、前記標的状態では含まれる場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＢ１４及びλｖｄｗＢ１４でそれぞれ乗算することであって、λｅｌｅｃＢ１４及びλｖｄｗＢ１４の両方が、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において原子価結合する場合、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれないが、前記標的状態では含まれる場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＢ１４及びλｖｄｗＢ１４でそれぞれ乗算することであって、λｅｌｅｃＢ１４及びλｖｄｗＢ１４の両方が、前記基準状態では０であり、前記標的状態では１であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では０から１まで変化する、乗算すること、及び
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれるが、前記標的状態では含まれない場合、（ａ）において得られるＡｉとＡｊとの間の前記非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を前記結合パラメータλｅｌｅｃＡ１４及びλｖｄｗＡ１４でそれぞれ乗算することであって、λｅｌｅｃＡ１４及びλｖｄｗＡ１４の両方が、前記基準状態では１であり、前記標的状態では０であり、且つ前記基準状態から前記標的状態への変換に沿った各中間系状態では１から０まで変化する、乗算することと、
（ｃ）（ｂ）において得られる前記計算静電１−４対相互作用及びファンデルワールス１−４対相互作用を前記対応する系状態の前記シミュレーションステップの前記総エネルギーに含めることと
を更に含んでなる、請求項１１に記載の方法。
【請求項１４】
前記状態の全てに対して分子シミュレーションを実行する手順は、
Ａａ及びＡｂが、前記基準状態において互いに原子価結合し、且つ前記標的状態において原子価結合しない場合、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれるが、前記標的状態では含まれない場合、λｖｄｗＡ１４が中間系状態で１よりも小さいとき、λｅｌｅｃＡ１４が前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従ってλｅｌｅｃＡ１４及びλｖｄｗＡ１４のそれぞれを変更し、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれないが、前記標的状態では含まれる場合、λｖｄｗＢ１４が中間系状態で１よりも小さいとき、λｅｌｅｃＢ１４が前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従ってλｅｌｅｃＢ１４及びλｖｄｗＢ１４のそれぞれを変更し、並びに
Ａａ及びＡｂが、前記基準状態において互いに原子価結合せず、且つ前記標的状態において原子価結合する場合、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれないが、前記標的状態では含まれる場合、λｖｄｗＢ１４が中間系状態で１よりも小さいとき、λｅｌｅｃＢ１４が前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従ってλｅｌｅｃＢ１４及びλｖｄｗＢ１４のそれぞれを変更し、
ＡｉとＡｊとの間の前記非結合１−４対相互作用が、前記基準状態では含まれるが、前記標的状態では含まれない場合、λｖｄｗＡ１４が中間系状態で１よりも小さいとき、λｅｌｅｃＡ１４が前記中間系状態で０であるように、前記基準状態から前記標的状態までの変換に沿った前記中間系状態のそれぞれのスケジュールに従ってλｅｌｅｃＡ１４及びλｖｄｗＡ１４のそれぞれを変更すること
を更に含んでなる、請求項１３に記載の方法。
【請求項１５】
前記ファンデルワールス相互作用の前記計算の少なくとも１つは、ソフトコアＬＪ相互作用ポテンシャルを使用することを含んでなる、請求項１１〜１４のいずれか一項に記載の方法。
【請求項１６】
前記基準状態及び前記標的状態のうちの少なくとも一方は、環構造であって、前記原子Ａａ及びＡｂが互いに結合し、且つ前記環構造の一部を形成する、環構造を有する分子を含んでなる、請求項１に記載の方法。
【請求項１７】
前記分子シミュレーションは、分子動的シミュレーション及びモンテカルロシミュレーションのうちの少なくとも一方を含んでなる、請求項１に記載の方法。
【請求項１８】
前記基準状態と前記標的状態との間の前記自由エネルギー差を計算する手順は、結合パラメータλの前記変更に関連付けられた仕事の特定及び解析により、前記標的状態、前記基準状態、及び前記中間系状態において得られるマイクロ状態の前記アンサンブルの解析を実行することを含んでなる、請求項１に記載の方法。
【請求項１９】
前記基準状態と前記標的状態との間の前記自由エネルギー差を計算する手順は、結合パラメータλが、マイクロ状態の選択されたアンサンブルで瞬時に変更される際、前記標的状態、前記基準状態、及び前記中間系状態において得られる前記マイクロ状態の前記アンサンブルの熱力学的特性の差の解析により、前記標的状態、前記基準状態、及び前記中間系状態において得られるマイクロ状態の前記アンサンブルの解析を実行することを含んでなる、請求項１に記載の方法。
【請求項２０】
前記アンサンブルは、ＮＶＴアンサンブル、ＮＰＴアンサンブル、ＮＶＥアンサンブル、及びμＶＴアンサンブルから選択される、請求項１９に記載の方法。
【請求項２１】
前記熱力学的特性の差の前記解析を実行することは、ＢＡＲ、ＭＢＡＲ、ＷＨＡＭ、ツヴァンツィヒ平均推定量から選択される推定量を適用することを含む、請求項１９に記載の方法。
【請求項２２】
前記熱力学的特性の差の前記解析を実行することは、自由エネルギー摂動ファミリ推定量の１つを適用することを含む、請求項１９に記載の方法。
【請求項２３】
前記基準状態と前記標的状態との間の前記自由エネルギー差を計算する手順は、前記結合パラメータλに関して前記標的状態、前記基準状態、及び前記中間系状態で得られるマイクロ状態のアンサンブルの熱力学的特性の導関数の熱力学的積分解析を実行することを含んでなる、請求項１に記載の方法。
【請求項２４】
前記アンサンブルは、ＮＶＴアンサンブル、ＮＰＴアンサンブル、ＮＶＥアンサンブル、及びμＶＴアンサンブルから選択される、請求項２３に記載の方法。
【請求項２５】
ＰＡ及びＰＢのうちの少なくとも一方はヌルである、請求項１に記載の方法。
【請求項２６】
命令を記憶する非一時的コンピュータ可読媒体において、前記命令は、１つ又は複数のプロセッサによって実行されると、基準状態と標的状態との間の自由エネルギー差の計算方法を実行し、前記基準状態及び前記標的状態はそれぞれ、原子ＰＡＢの共通の組を含み、前記基準状態は原子ＰＡの組を更に含み、前記標的状態は原子ＰＢの組を更に含み、前記組ＰＡは、前記基準状態のみに存在し、且つ前記標的状態には存在せず、及び前記組ＰＢは、前記標的状態のみに存在し、且つ前記基準状態には存在せず、少なくとも２つの原子Ａａ及びＡｂが存在し、Ａａ及びＡｂは、（１）前記基準状態では互いに原子価結合せず、且つ前記標的状態では原子価結合するか、又は（２）基準状態では互いに原子価結合し、且つ前記標的状態では互いに原子価結合しないかのいずれかであり、前記方法は、
（ａ）ＰＡ、ＰＢ、及びＰＡＢ内の全ての前記原子について、前記原子間の結合接続及び前記原子の相対空間構成を含むトポロジを提供する手順と、
（ｂ）前記基準状態と前記標的状態との間の変換路に沿った１つ又は複数の中間系状態を特定する手順であって、前記変換路は、各系状態の原子間相互作用から生じるエネルギーを変調する結合パラメータλによって定義され、前記結合パラメータλは、複数の成分であって、それぞれが［０，１］に属する値を有する複数の成分を含み、且つ異なるタイプの相互作用エネルギーを変調する、前記基準状態と前記標的状態との間の変換路に沿った１つ又は複数の中間系状態を特定する手順と、
（ｃ）前記基準状態、前記標的状態、及び中間系状態のマイクロ状態のアンサンブルを取得するために、少なくとも１つコンピュータプロセッサを使用して、分子シミュレーションを実行する手順であって、前記系状態のそれぞれに対して分子シミュレーションを実行することは、前記原子Ａａ及びＡｂ間の結合伸縮相互作用エネルギーを計算することを含み、前記結合伸縮相互作用エネルギーは、弱結合ポテンシャルによって定義され、
前記弱結合ポテンシャルは、前記結合パラメータλの結合伸縮成分λｓｂｓの関数であり、且つ［０，１］内のλｓｂｓの全ての値及びＡａとＡｂとの間の距離ｒの全ての値について、いかなる特異領域も含まず、前記弱結合ポテンシャルは以下の条件：λｓｂｓが（０，１）内にある場合、ＡａとＡｂとの間の前記距離が無限に近づくとき、前記弱結合ポテンシャルは平坦であり、Ａａ及びＡｂが前記基準状態又は前記標的状態のいずれかで原子価結合されないとき、前記弱結合ポテンシャルは、ＡａとＡｂとの間の全ての距離について平坦且つゼロであり、並びにＡａ及びＡｂが前記標的状態又は前記基準状態のいずれかで原子価結合されるとき、前記弱結合ポテンシャルは調和ポテンシャルに戻ることを更に満たす、分子シミュレーションを実行する手順と、
（ｄ）少なくとも１つのコンピュータプロセッサを使用して、前記標的状態、前記基準状態、及び前記中間系状態において得られるマイクロ状態の前記アンサンブルの解析により、前記基準状態と前記標的状態との間の前記自由エネルギー差を計算する手順と
を含んでなる、非一時的コンピュータ可読媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は変更された結合伸縮ポテンシャルを使用して自由エネルギー差を計算する方法及びシステムに関する。
【背景技術】
【０００２】
自由エネルギーは、化学系及び生物系の特徴付けにおいて重要な役割を果たす基本的な分子特性である。タンパク質−リガンド結合等の多くの化学的過程及び生物学的過程の自由エネルギー挙動の理解は、合理的な薬剤設計（生体分子標的に結合する小分子の設計を含む）等の取り組みにおいて極めて重要なものであり得る。
【０００３】
コンピュータモデリング及びシミュレーションは、多くの場合、自由エネルギーの研究で使用される。大半の場合、シミュレーションから正確な絶対自由エネルギーを評価することは、仮に可能であったとしても極めて困難である。したがって、異なるが類似するリガンド（例えば、同種リガンド）の実測親和性に基づくリガンドの相対結合親和性等の特定の系への洞察を提供するために、２つの明確に線引きされる熱力学的状態間の自由エネルギー差、すなわち相対自由エネルギーが、研究系として使用されることが多い。
【０００４】
相対自由エネルギーの計算では、２つの熱力学的状態は、基準系状態及び標的系状態と呼ぶことができ、これらはそれぞれ、最初の分子等の分子系の初期状態及び１つ又は複数の変換（配座変更、トポロジ変更、又はある原子若しくは化学基が別の原子若しくは化学基で置換されること（すなわち、突然変異）等）が行われた後の分子の終了状態を表すことができる。そのような変換が、常に現実的な物理的変換を表すわけではなく、非物理的、すなわち「錬金術的」変換を含み得ることに留意する。自由エネルギー摂動（ＦＥＰ）、熱力学的積分（ＴＩ）、及びアンブレラサンプリング等の、自由エネルギー差を計算する様々な枠組みが開発されてきた。
【０００５】
自由エネルギー摂動の枠組み内では、２つの系状態ａ及びｂ間の自由エネルギー差ΔＦａ→ｂは、
【０００６】
【数１】
0006298167.tif 000002
によって表すことができ、式中、β−１＝ｋＢＴであり、ここで、ｋＢはボルツマン定数であり、Ｔは温度である。Ｈａ（ｘ，ｐｘ）及びＨｂ（ｘ，ｐｘ）はそれぞれ、状態ａ及びｂのハミルトニアン特性である。＜．．．＞ａは、初期の基準状態ａを表す構成にわたるアンサンブル平均を示す。
【０００７】
自由エネルギー摂動の実際用途では、２つの熱力学的状態間の変換は通常、ａとｂとを結ぶ明確な経路に沿った非物理的な中間状態間の一連の変換によって達成される。この経路は多くの場合、結合パラメータλと呼ばれることが多い汎用範囲パラメータによって特徴付けられ、結合パラメータλは、基準状態から標的状態まで０から１まで変化し、
Ｈ（λ）＝（１−λ）Ｈａ＋λＨｂ      （２）
により、２つの状態のハミルトニアンを関連付け、式中、Ｈ（λ）は、２つの状態間の系のλ結合又は混成ハミルトニアンである（λが０及び１の終了値を取る場合の２つの状態
を含む）。したがって、ａとｂとの間の自由エネルギー差ΔＦａ→ｂは、
【０００８】
【数２】
0006298167.tif 000003
になり、式中、Ｎは、基準（初期）状態と標的（最終）状態との間の隣り合う状態間の「窓」数を表し、λｉは、初期状態、中間状態、及び最終状態での結合パラメータの値である。
【０００９】
基準系状態ａと標的系状態ｂとの間の自由エネルギー差は、熱力学的積分法を使用して計算することもでき、この場合、自由エネルギー差は、以下の式：
【００１０】
【数３】
0006298167.tif 000004
を使用して計算され、式中、λは、基準状態から標的状態までにおいて０から１まで変化する結合パラメータであり、Ｈ（λ）は、２つの状態間の系のλ結合又は混成ハミルトニアンであり（λが０及び１の終了値を取る場合の２つの状態を含む）、
【００１１】
【数４】
0006298167.tif 000005
は、結合パラメータλに対する結合ハミルトニアンの一次導関数である。
【００１２】
ＴＩの実際用途では、基準系状態と標的系状態との間の変換は、ａとｂとを結ぶ明確な経路に沿った一連の変換によって達成され、アンサンブル平均
【００１３】
【数５】
0006298167.tif 000006
は、基準系状態、中間非物理的状態、及び標的系状態を含め、サンプリングされる全ての状態で計算される。次に、基準系状態と標的系状態との間の自由エネルギー差が、
【００１４】
【数６】
0006298167.tif 000007
の値に基づいて、上記積分の数値積分によって近似され、式中、λｉは、初期状態、中間状態、及び最終状態での結合パラメータの値である。
【発明の概要】
【発明が解決しようとする課題】
【００１５】
従来の方法では、分子の環を開いて線形構造にするか、又は分子の線形構造を閉じて環にする自由エネルギーの計算は困難であり得る。これらの場合、２つの環原子を結ぶ原子価結合をオン及びオフに切り換える可逆仕事を計算する必要がある。可能な一手法は、環全体をなくし、対応する線形構造を最初から成長させることであり得るが、これは計算的に非常に非効率的であり、環が非常に大きいか、又は環が他の環と融合する場合、不可能なことがある。
【課題を解決するための手段】
【００１６】
本発明の幾つかの実施形態では、基準状態と標的状態との間の自由エネルギー差を計算するコンピュータ実施方法が提供される。基準状態及び標的状態はそれぞれ、原子ＰＡＢの共通の組を含む。基準状態は原子ＰＡの組を更に含み、標的状態は原子ＰＢの組を更に含む。組ＰＡは、基準状態のみに存在し、且つ標的状態には存在せず、及び組ＰＢは、標的状態のみに存在し、且つ基準状態には存在しない。少なくとも２つの原子Ａａ及びＡｂが存在し、Ａａ及びＡｂは、（１）基準状態では互いに原子価結合せず、標的状態では原子価結合するか、又は（２）基準状態では互いに原子価結合し、且つ標的状態では互いに原子価結合しないかのいずれかである。本方法は、
（ａ）ＰＡ、ＰＢ、及びＰＡＢ内の全ての原子について、原子間の結合接続及び原子の相対空間構成を含むトポロジを提供することと、
（ｂ）基準状態と標的状態との間の変換路に沿った１つ又は複数の中間系状態を特定することであって、変換路は、各系状態の原子間相互作用から生じるエネルギーを変調する結合パラメータλによって定義され、結合パラメータλは、複数の成分であって、それぞれが［０，１］に属する値を有する複数の成分を含み、且つ異なるタイプの相互作用エネルギーを変調する、特定することと、
（ｃ）基準状態、標的状態、及び中間状態のマイクロ状態のアンサンブルを取得するために、少なくとも１つのコンピュータプロセッサを使用して、モンテカルロシミュレーション又は分子動力学シミュレーション等の分子シミュレーションを実行することであって、系状態のそれぞれに対して分子シミュレーションを実行することは、原子Ａａ及びＡｂ間の結合伸縮相互作用エネルギーを計算することを含み、結合伸縮相互作用エネルギーは、弱結合ポテンシャルによって定義される、実行することと、
（ｄ）少なくとも１つのコンピュータプロセッサを使用して、標的状態、基準状態、及び中間状態において得られるマイクロ状態のアンサンブルの解析により、基準状態と標的状態との間の自由エネルギー差を計算することと
を含む。
【００１７】
本方法の幾つかの実施形態では、弱結合ポテンシャルは、結合パラメータλの結合伸縮成分λｓｂｓの関数であり、且つ［０，１］内のλｓｂｓの全ての値及びＡａとＡｂとの間の距離ｒの全ての値について、いかなる特異領域も含まない。弱結合ポテンシャルは以下の条件：λｓｂｓが（０，１）内にある場合、ＡａとＡｂとの間の距離が無限に近づくとき、弱結合ポテンシャルは平坦であり、Ａａ及びＡｂが基準状態又は標的状態のいずれかで原子価結合されないとき、弱結合ポテンシャルは、ＡａとＡｂとの間の全ての距離について平坦且つゼロであり、並びにＡａ及びＡｂが標的状態又は基準状態のいずれかで原子価結合されるとき、弱結合ポテンシャルは調和ポテンシャルに戻ることを更に満たす。
【００１８】
本方法の一実施形態では、λｓｂｓに対する弱結合ポテンシャルの一次導関数は、０〜と１のλｓｂｓの全ての値及びＡａとＡｂとの間の距離ｒの全ての値について連続であり、且つ結合される。別の実施形態では、ＡａとＡｂとの間の距離ｒに対する弱結合ポテンシャルの一次導関数及び二次導関数の両方が、λｓｂｓの全ての値について連続であり、
且つ結合され、及びｒが無限に近づくとき、ゼロに近づく。
【００１９】
本方法の幾つかの実施形態では、弱結合ポテンシャルは、（ｒ−ｒ０）２の関数であり、ここで、ｒ０はＡａとＡｂとの間の平衡距離である。一実施形態では、弱結合ポテンシャルは、
【００２０】
【数７】
0006298167.tif 000008
によって表され、式中、ｋは定数であり、並びに関数ｆ、ｇ、及びαはそれぞれ連続関数であり、且つ以下の条件：
ｆ（λｓｂｓ＝０）＝０、
ｆ（λｓｂｓ＝１）＝１、
ｇ（λｓｂｓ＝０）＝１、
ｇ（λｓｂｓ＝１）＝０、
ａ（ｋ，λｓｂｓ＜１）＞０
を満たす。
【００２１】
特定の実施形態では、ｆ（λｓｂｓ）＝λｓｂｓ、ｇ（λｓｂｓ）＝１−λｓｂｓ、及びα（ｋ，λｓｂｓ）は非負の定数である。
幾つかの実施形態では、系状態のそれぞれに対して分子シミュレーションを実行することは、Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において原子価結合しない場合、中間状態のそれぞれでのＡａとＡｂとの間の結合伸縮相互作用エネルギーを計算するために、λｓｂｓＡのスケジュール及び対応する弱結合ポテンシャルを使用することであって、λｓｂｓＡは、基準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、使用することと、Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において原子価結合する場合、ＡａとＡｂとの間の結合伸縮相互作用を計算するために、中間状態のそれぞれ及びλｓｂｓＢに対応する弱結合ポテンシャルについて、λｓｂｓＢのスケジュール及び対応する弱結合ポテンシャルを使用することであって、λｓｂｓＢは、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では０から１まで変化する、使用することとを含む。
【００２２】
本方法の幾つかの実施形態では、系状態のそれぞれに対して分子シミュレーションを実行することは、
（ａ）適切なパラメータを力場の結合角相互作用に使用して、（ｉ）Ａａ及び別の原子Ａｃによって形成される結合と、（ｉｉ）基準状態から標的状態への変換により破断されているか、又は形成されているＡａとＡｂとの結合との間の結合角相互作用を計算することと、
（ｂ）Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において互いに原子価結合しない場合、（ａ）において得られる計算結合角相互作用を結合角結合パラメータλｂａＡで乗算することであって、λｂａＡは、基準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において互いに原子価結合する場合、（ａ）において得られる計算結合角相互作用を結合角結合パラメータλｂａＢで乗算することであって、λｂａＢは、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では０から１まで変化
する、乗算することと、
（ｃ）（ｂ）において得られる結合角相互作用を対応する系状態のシミュレーションステップの総エネルギーに含めることと
を更に含む。
【００２３】
本方法の幾つかの実施形態では、系状態のそれぞれに対して分子シミュレーションを実行することは、
（ａ）力場の二面角相互作用に適切なパラメータを使用して、４つの接続原子群｛Ａｉ，Ａｊ，Ａｋ，Ａｌ｝の二面角相互作用を計算することであって、群はＡａ及びＡｂの両方を含む、計算することと、
（ｂ）Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において互いに原子価結合しない場合、（ａ）において得られる計算二面角相互作用を二面角結合パラメータλｂｄＡで乗算することであって、λｂｄＡは、基準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において互いに原子価結合する場合、（ａ）において得られる計算二面角相互作用を二面角結合パラメータλｂｄＢで乗算することであって、λｂｄＢは、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では０から１まで変化する、乗算することと、
（ｃ）（ｂ）において得られる二面角相互作用を対応する系状態のシミュレーションステップの総エネルギーに含めることと
を更に含む。
【００２４】
上記実施形態では、Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において互いに原子価結合しない場合、結合角相互作用結合パラメータλｂａＡ及び結合二面角相互作用結合パラメータλｂｄＡはそれぞれ、λｓｂｓＡが予め定義される閾値よりも小さい場合、０であるように選択され得、並びにＡａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において互いに原子価結合する場合、結合角相互作用結合パラメータλｂａＢ及び結合二面角相互作用結合パラメータλｂｄＢはそれぞれ、λｓｂｓＢが予め定義される閾値よりも小さい場合、０であるように選択され得る。
【００２５】
本方法の幾つかの実施形態では、状態の全てに対して分子シミュレーションを実行することは、
（ａ）力場の静電相互作用及びファンデルワールスに適切なパラメータを使用して、２つの原子Ａｉ及びＡｊ間の非結合静電相互作用及びファンデルワールス相互作用を計算することであって、対（Ａｉ，Ａｊ）の非結合排除ステータスは、基準状態から標的状態への変換による影響を受ける、計算することと、
（ｂ）Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において原子価結合せず、及びＡｉとＡｊとの間の非結合相互作用が、基準状態では排除されるが、標的状態では排除されない場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電相互作用及びファンデルワールス相互作用を結合パラメータλｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘでそれぞれ乗算することであって、λｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘの両方が、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では０から１まで変化する、乗算すること、
Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において原子価結合し、及びＡｉとＡｊとの間の非結合相互作用が、基準状態では排除されないが、標的状態では排除される場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電相互作用及びファンデルワールス相互作用を結合パラメータλｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘでそれぞれ乗算することであって、λｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘの両方が、基
準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算することと、
（ｃ）（ｂ）において得られる計算非結合静電相互作用及びファンデルワールス相互作用を対応する系状態のシミュレーションステップの総エネルギーに含めることと
を更に含む。
【００２６】
本方法の幾つかの実施形態では、状態の全てに対して分子シミュレーションを実行することは、
Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において原子価結合せず、及びＡｉとＡｊとの間の非結合相互作用が、基準状態では排除されるが、標的状態では排除されない場合、λｖｄｗＡｅｘが中間状態で１よりも小さいとき、λｅｌｅｃＡｅｘがその中間状態で０であるように、基準状態から標的状態までの変換に沿った前記中間状態のそれぞれのスケジュールに従って、λｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘのそれぞれを変更し、並びに
Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において原子価結合し、及びＡｉとＡｊとの間の非結合相互作用が、基準状態では排除されないが、標的状態では排除される場合、λｖｄｗＡｅｘが中間状態で１よりも小さいとき、λｅｌｅｃＡｅｘがその中間状態で０であるように、基準状態から標的状態までの変換に沿った前記中間状態のそれぞれのスケジュールに従って、λｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘのそれぞれを変更すること
を更に含む。
【００２７】
本方法の幾つかの実施形態では、状態の全てに対して分子シミュレーションを実行することは、
（ａ）力場の静電１−４対相互作用及びファンデルワールス１−４対相互作用に適切なパラメータを使用して、別の２つの介在する原子と一緒に、基準状態又は標的状態のいずれかで結合二面角相互作用を形成する２つの原子Ａｉ及びＡｊ間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を計算することと、
（ｂ）Ａａ及びＡｂが、基準状態において互いに原子価結合し、且つ標的状態において原子価結合しない場合、
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれるが、標的状態では含まれない場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＡ１４及びλｖｄｗＡ１４でそれぞれ乗算することであって、λｅｌｅｃＡ１４及びλｖｄｗＡ１４の両方が、基準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算すること、及び
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれないが、標的状態では含まれる場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＢ１４及びλｖｄｗＢ１４でそれぞれ乗算することであって、λｅｌｅｃＢ１４及びλｖｄｗＢ１４の両方が、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算すること、並びに
Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において原子価結合する場合、
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれないが、標的状態では含まれる場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＢ１４及びλｖｄｗＢ１４でそれぞれ乗算することであって、λｅｌｅｃＢ１４及びλｖｄｗＢ１４の両方が、基準状態では０であり、標的状態では１であり、且つ基準状態から標的状態への変換に沿った各中間状態では０から１まで変化する、乗算すること、及び
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれるが、標的状態では含まれない場合、（ａ）において得られるＡｉとＡｊとの間の非結合静電１−４対相互作用及びファンデルワールス１−４対相互作用を結合パラメータλｅｌｅｃＡ１４及びλｖｄｗＡ１４でそれぞれ乗算することであって、λｅｌｅｃＡ１４及びλｖｄｗＡ１４の両方が、基準状態では１であり、標的状態では０であり、且つ基準状態から標的状態への変換に沿った各中間状態では１から０まで変化する、乗算することと、
（ｃ）（ｂ）において得られる計算静電１−４対相互作用及びファンデルワールス１−４対相互作用を対応する系状態のシミュレーションステップの総エネルギーに含めることとを更に含む。
【００２８】
本方法の幾つかの実施形態では、状態の全てに対して分子シミュレーションを実行することは、
（ａ）ＡａとＡｂが、基準状態において互いに原子価結合し、且つ標的状態において原子価結合しない場合、
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まるが、標的状態では含まれない場合、λｖｄｗＡ１４が中間状態で１よりも小さいとき、λｅｌｅｃＡ１４がその中間状態で０であるように、基準状態から標的状態までの変換に沿った中間状態のそれぞれのスケジュールに従ってλｅｌｅｃＡ１４及びλｖｄｗＡ１４のそれぞれを変更し、ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれないが、標的状態では含まれる場合、λｖｄｗＢ１４が中間状態で１よりも小さいとき、λｅｌｅｃＢ１４がその中間状態で０であるように、基準状態から標的状態までの変換に沿った中間状態のそれぞれのスケジュールに従ってλｅｌｅｃＢ１４及びλｖｄｗＢ１４のそれぞれを変更し、並びに
（ｂ）Ａａ及びＡｂが、基準状態において互いに原子価結合せず、且つ標的状態において原子価結合する場合、
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれないが、標的状態では含まれる場合、λｖｄｗＢ１４が中間状態で１よりも小さいとき、λｅｌｅｃＢ１４がその中間状態で０であるように、基準状態から標的状態までの変換に沿った中間状態のそれぞれのスケジュールに従ってλｅｌｅｃＢ１４及びλｖｄｗＢ１４のそれぞれを変更し、
ＡｉとＡｊとの間の非結合１−４対相互作用が、基準状態では含まれるが、標的状態では含まれない場合、λｖｄｗＡ１４が中間状態で１よりも小さいとき、λｅｌｅｃＡ１４がその中間状態で０であるように、基準状態から標的状態までの変換に沿った中間状態のそれぞれのスケジュールに従ってλｅｌｅｃＡ１４及びλｖｄｗＡ１４のそれぞれを変更すること
を更に含む。
【００２９】
本方法の幾つかの実施形態では、ファンデルワールス相互作用の計算は、ソフトコアＬＪ相互作用ポテンシャルを使用することを含み得る。
本方法の幾つかの実施形態では、基準状態及び標的状態のいずれかは、環構造であって、原子Ａａ及びＡｂが互いに結合し、且つ環構造の一部を形成する、環構造を有する分子を含む。
【００３０】
本方法の幾つかの実施形態では、基準状態と標的状態との間の自由エネルギー差を計算することは、結合パラメータλの変更に関連付けられた仕事の特定及び解析により、標的状態、基準状態、及び中間状態において得られるマイクロ状態のアンサンブルの解析を実行することを含む。
【００３１】
本方法の幾つかの実施形態では、基準状態と標的状態との間の自由エネルギー差を計算することは、結合パラメータλが、マイクロ状態の選択されたアンサンブルで瞬時に変更
される際、標的状態、基準状態、及び中間状態において得られるマイクロ状態の適するアンサンブルの熱力学的特性の差の解析により、標的状態、基準状態、及び中間状態において得られるマイクロ状態のアンサンブルの解析を実行することを含む。熱力学的特性の差の解析を実行することは、ＢＡＲ、ＭＢＡＲ、ＷＨＡＭ、ツヴァンツィヒ平均推定量、又は自由エネルギー摂動ファミリ推定量の１つから選択される推定量を適用することを含む。アンサンブルは、例えば、ＮＶＴアンサンブル、ＮＰＴアンサンブル、ＮＶＥアンサンブル、及びμＶＴアンサンブルであることができる。
【００３２】
本方法の幾つかの実施形態では、基準状態と標的状態との間の自由エネルギー差を計算することは、結合パラメータλに関して標的状態、基準状態、及び中間状態で得られるマイクロ状態の適するアンサンブルの熱力学的特性の導関数の熱力学的積分解析を実行することを含む。アンサンブルは、例えば、ＮＶＴアンサンブル、ＮＰＴアンサンブル、ＮＶＥアンサンブル、及びμＶＴアンサンブルであることができる。
【００３３】
本発明は、１つ又は複数のプロセッサと、１つ又は複数のプロセッサに動作可能に結合され、プロセッサによって実行可能な命令を含むメモリとを含む装置も提供し、１つ又は複数のプロセッサは、命令を実行する際、本明細書に記載される方法の様々な実施形態を実行するように動作可能である。本発明は、命令を記憶する非一時的コンピュータ可読媒体を更に提供し、命令は、１つ又は複数のプロセッサによって実行されると、本明細書に記載される方法の様々な実施形態を実行する。
【００３４】
本発明の他の特徴及び利点は、以下の説明及び特許請求の範囲から明らかになる。
本発明は、添付図面を参照することにより、よりよく理解されよう。
【図面の簡単な説明】
【００３５】
【図１ａ】本発明の幾つかの実施形態による、基準系状態と標的系状態との間の例示的な閉環変換を示す図。
【図１ｂ】本発明の幾つかの実施形態による、基準系状態と標的系状態との間の例示的な閉環変換を示す図。
【図１ｃ】本発明の幾つかの実施形態による、基準系状態と標的系状態との間の例示的な閉環変換を示す図。
【図２】本発明の幾つかの実施形態による、結合パラメータλの様々な値での結合の形成又は破断変換に関わる２個の原子間の原子間距離の関数としての弱結合ポテンシャルのプロット。
【図３】本発明の幾つかの実施形態による、結合パラメータλの様々な値での結合の形成又は破断変換に関わる２個の原子間の原子間距離の関数としての弱結合ポテンシャルのプロットと、結合パラメータλに関するポテンシャルの導関数及び粒子間距離に関するポテンシャルの一次導関数及び二次導関数とを示すグラフ。
【図４】自由エネルギー変化が本発明の実施形態の適用に適している、環状構造で生じる例示的な開環変換を示す図。
【図５】自由エネルギー変化が本発明の実施形態の適用に適している、２つの接続された環状構造間の例示的な結合破断を示す図。
【図６】自由エネルギー変化が本発明の実施形態の適用に適している、基準系での非プロリン残基が標的系でプロリン残基に変換する、タンパク質残留物突然変異変換の図。
【図７】自由エネルギー変化が本発明の実施形態の適用に適している、２つの近傍にある残基間の直接接続を生じさせる、基準系での残基（プロリン）が標的系で削除されるタンパク質残基削除変換の図。
【発明を実施するための形態】
【００３６】
本願は、基準系状態と標的系状態との間の自由エネルギー差を計算するコンピュータ実
施方法及びシステムを開示する。特に、環状構造を線形構造に変換する環状構造の結合破断と、線形構造から環状構造に変換する結合形成と（それぞれについては更に後述する）から生じる自由エネルギー差を特定することにおける問題に対処するために、本願に開示される方法及びシステムは、結合伸縮に関数形態を利用し、それにより、自由エネルギー差を計算するために、調和結合関数形態への厳密な接続を錬金術的変換の任意の時点で維持することができる。したがって、本願の方法及びシステムは、有利には、自由エネルギー計算の数値的安定性及び精度を改善することができる。しかし、本明細書に開示されるそのような変更された結合伸縮ポテンシャルを使用する自由エネルギー計算の一般原理が、任意の結合形成及び破断状況に概して適用可能であり、閉環又は開環に限定されないことを認識されたい。
【００３７】
従来の自由エネルギー差計算と同様に、系内の原子は、異なる系状態の系エネルギーを評価するために、異なる群に分類することができる。基準状態及び標的状態の両方が原子ＰＡＢの共通の組を含む。基準状態は原子ＰＡの組を更に含み、及び標的状態は原子ＰＢの組を更に含む。原子ＰＡの組は、基準状態のみに存在し、且つ標的状態には存在せず、及び原子ＰＢの組は、標的状態のみに存在し、且つ基準状態には存在しない。環形成状況では、ＰＡは、２個の末端原子に接続して結合を形成する原子であることができる。変換の過程中、ＰＡ及びＰＢ内の原子は、それ自体の組内の他の原子及びＰＡＢ内の他の原子と相互作用するが、ＰＡ内の原子はＰＢ内のいかなる原子とも相互作用せず、逆も同様である。例えば、結合が２個の原子Ａ１及びＡ２間に形成されようとしている図１ａに示される構造を有する分子の場合、ＰＡＢ＝｛Ａ１，Ａ２，Ａ３，Ａ４，Ａ５，Ａ６，Ａ７，Ａ８，Ａ９｝（原子は、例示のみを目的として原子価４を有するものとして示されている：他の原子又は原子群がこの構造で表現可能なことが理解される）であり、ＰＡ＝｛Ａ１０，Ａ１１，｝であり、ＰＢ＝φである。別の例として、Ａ１２、Ａ１３、Ａ１４を含む分子断片がＡ１とＡ３との間に挿入されて、２つの結合Ａ１−Ａ１２及びＡ３−Ａ１２を形成すべきである、図１ｂに示される骨格構造Ａ１−Ａ３−Ａ２を有する分子の場合、ＰＡＢ＝｛Ａ１，Ａ２，Ａ３，Ａ４，Ａ５，Ａ６，Ａ７，Ａ８，Ａ９｝であり、ＰＡ＝｛Ａ１０，Ａ１１｝であり、ＰＢ＝｛Ａ１２，Ａ１３，Ａ１４｝である。後者の例では、断片を挿入して２つの結合を形成することは、２つのステップを取るものと考えることができる：第１のステップは、Ａ１２を含む断片によりＡ１（又はＡ３）上で鎖Ａ３−Ａ２−Ａ１を線形に成長させることであり、第２のステップは、Ａ１２−Ａ３（又はＡ１２−Ａ１）間で環を閉じることである。閉じられた四員環構造を有する分子が、錬金術的に五員環構造に変換される、図１ｃに示される別の例では、ＰＡＢ＝｛Ａ１，Ａ２，Ａ３，Ａ４，Ａ５，Ａ６，Ａ７，Ａ８，Ａ９，Ａ１０，Ａ１１，Ａ１２｝であり、ＰＡ＝φであり、ＰＢ＝｛Ａ１３，Ａ１４，Ａ１５｝である。
【００３８】
本明細書に記載されるように、環の開閉を含む錬金術的変換の結合ハミルトニアンＨ（λ）は一般に、以下の項：Ｈｂｓ（λｂｓ）、Ｈｂａ（λｂａ）、Ｈｂｄ（λｂｄ）、及びＨｎｂ（λｎｂ）を含むことができ、これらはそれぞれ、結合伸縮項、結合角項、結合二面角項、並びに非結合排除及び１−４対相互作用項に対応する。考察を簡易化するために、以下では、２つの環原子間の原子価結合が形成されている（すなわち、閉環を含む、基準状態から標的状態への変換）場合を詳細に説明する。
【００３９】
初期ステップとして、系内の原子間の結合接続と、ＰＡＢ、ＰＡ、及びＰＢのそれぞれを形成する原子の相対空間構成とを含む、系のトポロジが提供される。
基準状態と標的状態との間の１つ又は複数、例えば複数の中間状態が、異なる値の結合パラメータλによって定義される経路に沿って特定することができ、λの値の増分は、系を基準状態から標的状態に移す。λは、０から１に変化するスカラー変数であることができるが、以下に更に考察する実施形態等の本発明の幾つかの実施形態では、λは、系内の異なるタイプの相互作用に異なる成分を含むベクトルであることができる。分子動力学又
はモンテカルロシミュレーション等であるが、これらに限定されないコンピュータ分子シミュレーションが実行されて、基準状態、標的状態、及び各中間状態でのマイクロ状態のアンサンブルを取得することができる。中間状態のλ値は、基準状態から標的状態への「反応経路」上のそれぞれ隣り合うλ窓間に、分子シミュレーションによってサンプリングされる連続したλ窓内のマイクロ状態間にかなりの重複があるように、既知の技法によって選ぶことができる。
【００４０】
これらの全ての状態に対して分子シミュレーションを実行するに当たり、結合を形成すべき２個の原子Ａａ及びＡｂ（例えば、図１ａではＡ１及びＡ３）間の結合伸縮相互作用エネルギーは、λによって変調される弱結合ポテンシャル（又はその結合伸縮成分）によって定義することができる。λ＝０である（Ａａ及びＡｂが、基準状態で完全に非結合である）場合、弱結合ポテンシャルは、ＡａとＡｂとの間の全ての距離ｒに対して平坦ポテンシャルである。０＜λ＜１である（ＡａとＡｂとの結合が、錬金術的変換において「部分的に形成」中である）場合、ｒ→∞であるとき、弱結合ポテンシャルは平坦ポテンシャルに向かって横ばいになる。すなわち、ＡａとＡｂとの間の距離ｒに対するポテンシャルの偏導関数は、ｒ→∞であるとき、ゼロである。λ＝１である（Ａａ及びＡｂが、標的状態において完全に原子価結合される）場合、弱結合ポテンシャルは調和ポテンシャルに戻る。更に、結合伸縮項のポテンシャルエネルギー関数は、［０，１］内の結合パラメータλの結合伸縮成分λｓｂｓの全ての値及びＡａとＡｂとの間の距離ｒの全ての値に対して、いかなる特異領域も有さない。弱結合ポテンシャル及び弱結合ポテンシャルの幾つかの特性を明らかにする詳細を以下に提供する。
【００４１】
ＯＰＬＳ、ＣＨＡＲＭＭ、及びＡＭＢＥＲ等の普及している分子力学力場では、２個の原子間の結合伸縮相互作用は、以下の形態：
【００４２】
【数８】
0006298167.tif 000009
の調和結合によってモデリングされ、式中、ｋは、使用される力場の結合の強度（又は剛性）を定義する「力の定数」又は「フック定数」であり、ｒは２個の原子間の瞬間距離であり、ｒ０は２個の原子間の平衡距離である。基準状態のハミルトニアンと標的状態のハミルトニアンとの間で結合パラメータを線形にスケーリングする従来の方法では、結合伸縮項は、結合ハミルトニアンにおいて以下の形態：
【００４３】
【数９】
0006298167.tif 000010
を有する。したがって、
【００４４】
【数１０】
0006298167.tif 000011
である。上記式中の積分は、ｒが非常に大きくなるとき、無限に近づく。λが０に近づく限界では、２個の原子間に結合伸縮相互作用はなく、２個の原子間の距離は非常に大きくなることができ、上記積分の計算において特異点及び数値不安定性の問題を生じさせる。実際では、２個の原子間の距離は、シミュレーションボックスのサイズ（シミュレーションが行われ、周期境界条件等の境界条件が適用可能な三次元容積又は単位格子）によって制限され、特異値問題は回避可能である。しかし、上記式中の積分は無制限であるため、自由エネルギーの計算においてなお、数値の不安定性及び不正確性の問題を生じさせるおそれがある。
【００４５】
更に、上記の従来の結合ハミルトニアン関数形態を使用すると、ｒが無限に近づく場合、λ＝０に対してポテンシャルは未定義である。限界値は、λがいかに０に近づくか及びｒがいかに無限に近づくかに依存する。すなわち、
【００４６】
【数１１】
0006298167.tif 000012
である。
【００４７】
２個の原子の結合の消滅及び／又は形成を含む信頼性の高い自由エネルギーを特定可能な安定し効率的なシミュレーションを可能にする経路を得るために、本発明者らは、結合パラメータλが０〜１で変わるとき、２個の原子間の調和相互作用が完全にオン及びオフにされる、２つの物理系を接続する以下の結合ポテンシャル（本明細書では、弱結合ポテンシャルと呼ばれる）
【００４８】
【数１２】
0006298167.tif 000013
を発見し、式中、関数ｆ、ｇ、及びαはそれぞれ、連続関数であり、以下の条件を同時に満たす：ｆ（λ＝０）＝０、ｆ（λ＝１）＝１、ｇ（λ＝０）＝１、ｇ（λ＝１）＝０、α（ｋ，λ＜１）＞０。結合伸縮相互作用について弱結合ポテンシャルに関する本明細書での全ての考察において、式（式９以降）中で使用されるλが、結合パラメータλの結合伸縮成分λｓｂｓを指すことに留意する。
【００４９】
式９によって記述される弱結合ポテンシャルの特定の例では、ｆ（λ）＝λ、ｇ（λ）＝１−λ、及びα（ｋ，λ）＝α＝ｃｏｎｓｔ（定数）である。すなわち、弱結合ポテンシャルは以下の形態を取る。
【００５０】
【数１３】
0006298167.tif 000014
式１０の弱結合相互作用が、λ＝０及びλ＝１である場合、２つの物理的終了状態を正確に復元することを見て取ることができる。
【００５１】
【数１４】
0006298167.tif 000015
式１０の分母へのα（１−λ）（ｒ−ｒ０）２の導入は、λが１よりも小さい場合、調和相互作用を弱結合相互作用（相互作用は、ｒが無限に近づくとき、制限される）に変更し、λ＝１である場合（結合が形成される）、関数は厳密な調和ポテンシャル形態を有する。
【００５２】
上記関数形態は、従来の調和ポテンシャルに関連付けられる特異点及び数値不安定性の問題をなくす。以下、式１０により例示される弱結合相互作用関数形態の幾つかの特性について考察する。
特性１：
弱結合相互作用関数形態は、０〜１のλの全ての値及びｒの全ての値に対して、いかなる特異領域も有さない。上記説明から、すなわち、
【００５３】
【数１５】
0006298167.tif 000016
である。
【００５４】
力の定数ｋ＝２０ｋｃａｌ・ｍｏｌ−１・Å−２と、異なる値のλにおいて弱結合パラメータα＝１とを有するモデル系の場合での弱結合ポテンシャルは、式２に与えられる。ポテンシャルが、λ及びｒの全ての値で連続していることが明確に見て取ることができる。ポテンシャルは、λ＝１での調和ポテンシャルから中間値のλでの弱結合ポテンシャルまでゆっくりと変化し、λ＝０において０になる。
特性２：
【００５５】
【数１６】
0006298167.tif 000017
上記セクションで考察したように、自由エネルギー計算の数値安定性及び正確性を決めるのは、
【００５６】
【数１７】
0006298167.tif 000018
である。上記公式において、λ∈［０，１）である場合、式１３中の一番右側の括弧内で平均される熱力学的特性は、ｒの全ての値で連続しており、制限される。λ＝１である場合、弱結合ポテンシャルは調和ポテンシャルを復元し、したがって、ｒがｒ０に近い位相空間領域のみが、サンプリングされ、アンサンブル平均に入れられる。ｒがｒ０に近い位相空間領域では、式１３中の積分も制限される。したがって、式１３の括弧内の数量は、０〜１のλの全ての値で制限される。
【００５７】
【数１８】
0006298167.tif 000019
は、０〜１のλの全ての値に対していかなる特異領域も有さないため、上記弱結合相互作用関数形態を使用して、正確で信頼性の高い自由エネルギー結果を得ることができる。
【００５８】
例として、力の定数ｋ＝２０ｋｃａｌ・ｍｏｌ−１・Å−２と、弱結合パラメータα＝１とを有するモデル系の場合における、λ＝０．５での結合パラメータλに対する式１０の弱結合ポテンシャルの導関数を図３に描く。
【００５９】
【数１９】
0006298167.tif 000020
が、λ及びｒの全ての値に対していかなる特異領域も有さず、信頼性が高く正確な自由エネルギー計算を可能にすることを見て取ることができる。
特性３：
【００６０】
【数２０】
0006298167.tif 000021
粒子間（又は原子間）距離ｒに対する式１０の弱結合ポテンシャルの一次導関数及び二次導関数の両方が、λ∈［０，１］の全ての値で連続しており、制限され、狭い範囲を有し、ｒ→∞であるとき、０に近づく。λ＝１である場合、弱結合ポテンシャルは調和ポテンシャルに戻り、ｒがｒ０に近い位相空間領域のみが、分子シミュレーションにおいてサンプリングされる。したがって、粒子間距離ｒに対するポテンシャルの一次導関数及び二次導関数も、λ＝１である場合、物理的に許容可能な位相空間領域において連続しており、制限される。したがって、原子への力及び加速度は常に、０〜１のλの全ての値において連続しており、制限され、λの全ての値に対して安定した分子力学シミュレーションの実行を可能にする。例として、力の定数ｋ＝２０ｋｃａｌ・ｍｏｌ−１・Å−２と、弱結合パラメータα＝１とを有するモデル系の場合でのλ＝０．５での粒子間距離ｒに対する弱結合ポテンシャルの一次導関数及び二次導関数を図３に描く。それらの導関数が空間全体にわたって連続しており、安定した分子力学シミュレーションを可能にすることを見て取ることができる。
【００６１】
上記説明から、弱結合ポテンシャルが、従来の方法に関連する特異点及び数値不安定性の問題をなくし、基準系状態、中間系状態、及び標的系状態を含む全ての系状態で、安定した分子力学シミュレーション及びより収束するモンテカルロを実行可能なことが明らかである。本明細書に記載される弱結合ポテンシャルを使用して、原子価結合の破断及び形成に関わる基準系状態と標的系状態との間の自由エネルギー差は、正確で高い信頼性で計算され得る。そのような弱結合ポテンシャルを利用する自由エネルギー計算は、自由エネルギー計算の計算効率及び収束が従来の方法よりも大幅に改善することができる、開環、閉環、又は環再構成を含む錬金術的変換に特に有利である。
【００６２】
上記セクションで述べたように、従来の結合ハミルトニアンでは、ｒが無限に近づくとき、λ＝０である初期状態で限界は存在せず、一方、弱結合ポテンシャルでは、ｒの全ての値に対して限界が存在する。従来の結合ハミルトニアンは、α＝０である弱結合ポテンシャルの特別な事例である。従来の結合ハミルトニアン及び弱結合ポテンシャルは、λ＝１の場合、同じ終了状態に達するが、λ＝０の場合の初期状態は、弱結合パラメータαの値に応じて異なる。以下、初期状態の自由エネルギーが弱結合パラメータαに依存しないことが示される。
【００６３】
以下のポテンシャルエネルギー項（すなわち、λ＝０である式１０）を有するハミルトニアンを考える。
【００６４】
【数２１】
0006298167.tif 000022
α＝０である場合、これは従来の調和ポテンシャルになり、α≠０である場合、これは弱結合ポテンシャルになる。Ｆ（ａ＝０，λ＝０）＝Ｆ（α＝α’＞０，λ＝０）であることを証明する必要がある。なお、
【００６５】
【数２２】
0006298167.tif 000023
であり、式中、
【００６６】
【数２３】
0006298167.tif 000024
である。
【００６７】
アンサンブル平均Ｉ（α）は、
【００６８】
【数２４】
0006298167.tif 000025
であるため、有限数である。
【００６９】
したがって、
【００７０】
【数２５】
0006298167.tif 000026
である。
【００７１】
上記に鑑みて、λ＝０である場合の従来の結合ハミルトニアン及び弱結合ポテンシャルの初期状態は、同じ自由エネルギーを有する。
自由エネルギーの計算において、自由エネルギーの収束は、隣り合う又は連続λ窓間の位相空間領域の重複により影響を受け得る。経験的に、初期状態から最終状態への適する経路は、λの関数としての自由エネルギーの変化が、λの全ての値に対して連続し平滑する場合、すなわち、
【００７２】
【数２６】
0006298167.tif 000027
である場合、達成することができる。
【００７３】
したがって、好ましい実施形態では、弱結合パラメータαの値は、式２１が満たされる場合に得ることができる。
上述したように、開環及び閉環自由エネルギー計算では、結合伸縮相互作用に加えて、結合角項、結合二面角項、及び非結合相互作用を含め、初期状態と最終状態とで異なる他の相互作用エネルギー項があり得る。これらの異なるタイプの相互作用は、変換プロセス中、別様に処理することができる。本明細書に記載されるように、幾つかの実施形態では、結合パラメータλは、異なるタイプの相互作用の成分を含むベクトル形態を取ることができる。例えば、結合ベクトルλは、基準状態と標的状態との間での原子価結合の形成（又は破断）による影響を受ける相互作用に対する以下の成分を含むことができる：結合伸縮エネルギー項の成分λｓｂｓ、結合角エネルギー項の成分λｂａ、結合二面角エネルギー項の成分λｂｄ、非結合排除及び１−４対相互作用のそれぞれのファンデルワールスエネルギー項の成分λｖｄｗ及び静電エネルギー項の成分λｅｌｅｃ。結合ベクトルλは、基準状態と標的状態との間での原子価結合の形成及び／又は破断による影響を受けない他の相互作用の成分λｏｔｈｅｒを更に含むことができる。結合ベクトルλの各成分は［０，１］に属する。基準状態から標的状態への錬金術的変換の過程中、更に後述するように、基準状態で固有の相互作用は、様々なλ成分の第１の組のスケジュールに従ってオフにすることができ、標的状態で固有の相互作用は、様々なλ成分の第２の組のスケジュールに従ってオンにすることができる。
【００７４】
普及している分子力学力場では、結合角相互作用及び結合二面角相互作用は通常、以下のポテンシャルエネルギー形態
【００７５】
【数２７】
0006298167.tif 000028
を有し、式中、θは結合角であり、θ０は平衡結合角であり、ｋθは角力定数であり（θ０及びｋθの両方が、結合角を形成する原子に依存する）、φは二面角であり、ｋφは二面角力定数である（二面角を形成する原子に依存する）。環を開くか、又は閉じる場合、結合の破断又は形成により影響を受ける結合角項及び二面角項は、結合パラメータλの成分λｂａ及びλｂｄのそれぞれにより変調することができる。結合角相互作用項及び二面角相互作用項は、全てのλ値に対して制限されるが、分子が非常に捻れた幾何学的形状である場合、これらの項の絶対値は非常に大きくなり得る。自由エネルギーの計算の精度を改善するために、幾つかの実施形態では、結合伸縮相互作用はまず、相当な程度（例えば、λｓｂｓ＝０．５）までオンにされてから、結合角相互作用及び結合二面角相互作用をオンにすることができる（結合形成中）。このようにして、結合伸縮相互作用は、かなり捻れた幾何学的形状から分子を遠ざけ、自由エネルギー計算において結合角相互作用及び結合二面角相互作用によって生じる不正確性の問題を改善する。
【００７６】
２個の原子間の非結合静電相互作用及びファンデルワールス相互作用は通常、普及している分子力学力場では以下のポテンシャルエネルギー形態
【００７７】
【数２８】
0006298167.tif 000029
によってモデリングされ、式２４及び式２５中、ｒは原子間距離であり、ｑ１及びｑ２は２個の原子の電荷であり、Ｃは定数であり、εはＵｖｄｗ（ｒ）のポテンシャルウェルの深度であり、σは、原子間ポテンシャルＵｖｄｗ（ｒ）がゼロである有限距離である。
【００７８】
ＯＰＬＳ、ＣＨＡＲＭ、ＡＭＢＥＲを含め、多くの力場は、１つ、２つ、又は３つの結合で隔てられた原子間の非結合相互作用を排除するか、又は変更する。特に、２個の原子が３個の結合によって隔てられている（例えば、図１ａに示される２個の原子Ａ１及びＡ８）場合、２個の原子間の、式２４に記述されるような従来の非結合静電相互作用及びファンデルワールス相互作用は計算されない（すなわち、排除される）：その代わり、静電１−４対相互作用及びファンデルワールス１−４対相互作用で置換される。
【００７９】
環を開くか、又は閉じる場合、破断されるか、又は形成される結合にわたる原子の非結合静電相互作用及びファンデルワールス相互作用は、成分λｅｌｅｃ及びλｖｄｗによってそれぞれ変調することができ、これについて更に後述する。本明細書で使用される場合、本発明の幾つかの実施形態では、ファンデルワールス相互作用及び／又は静電相互作用は、全ての構成で対毎の相互作用エネルギーを有限に保つソフトコアレナードジョーンズ（ＬＪ）及び／又はソフトコアクーロンポテンシャルによって記述することができる。例として、ソフトコアＬＪポテンシャルは、以下の関数形態
【００８０】
【数２９】
0006298167.tif 000030
を取り得、式中、それαＬＪはソフトコアパラメータであり、ε及びσは標準レナードジョーンズ相互作用パラメータであり、ｎは、通常、１〜６の値を取る他のパラメータである。ボイトラー（Ｂｅｕｔｌｅｒ）ら「分子シミュレーションに基づく自由エネルギー計算での特異値及び数値不安定性の回避（Ａｖｏｉｄｉｎｇ  ｓｉｎｇｕｌａｒｉｔｉｅｓ
ａｎｄ  ｎｕｍｅｒｉｃａｌ  ｉｎｓｔａｂｉｌｉｔｉｅｓ  ｉｎ  ｆｒｅｅ  ｅｎｅｒｇｙ  ｃａｌｃｕｌａｔｉｏｎｓ  ｂａｓｅｄ  ｏｎ  ｍｏｌｅｃｕｌａｒ  ｓｉｍｕｌａｔｉｏｎｓ）」、ケミカル・フィジックス・レターズ（Ｃｈｅｍ．Ｐｈｙｓ．Ｌｅｔｔ．）、第２２２巻、１９９４年、ｐ．５２９〜５３９を参照のこと。この開示は参照により本明細書に援用される。ソフトコアＬＪ相互作用は、λ＝１である場合、標準レナードジョーンズ相互作用を復元し、λ＝０である場合、０になる。
【００８１】
幾つかの実施形態では、開環及び閉環自由エネルギー計算では、全ての系状態のλスケジュールが使用されて、環構造を開くか、又は閉じることに繋がる結合の形成及び／又は消滅により影響を受ける結合伸縮相互作用、結合角相互作用、結合二面角相互作用、非結合排除相互作用、及び１−４対相互作用を取り扱うことができる。結合の形成又は破断による影響を受ける相互作用の結合パラメータλは、結合破断変換に７つの成分と、基準状態と標的状態への結合形成変換に７つの成分とを含み得る。
【００８２】
結合破断変換に適用可能な結合パラメータλの７つの成分は、以下の項を含む：λｓｂｓＡ、これは、初期状態に存在する（結合は存在し、まだ破断されていない）が、最終状態に存在しない（結合は破断する）結合伸縮相互作用をそれぞれ変調する；λｂａＡ及びλｂｄＡ、これらは、初期状態に存在するが、最終状態に存在しない結合角相互作用及び結合二面角相互作用を変調する；
【００８３】
【数３０】
0006298167.tif 000031
及び
【００８４】
【数３１】
0006298167.tif 000032
、これらは、最終状態では存在するが、初期状態から最終状態への変換中に、初期状態に存在する結合が破断することに起因して、初期状態では排除（非結合排除）される非結合静電相互作用及びファンデルワールス相互作用を変調する；並びに
【００８５】
【数３２】
0006298167.tif 000033
及び
【００８６】
【数３３】
0006298167.tif 000034
、これらは、初期状態では存在するが、初期状態から最終状態への変換中の結合の破断又は形成に起因して、最終状態には含まれない静電１−４対相互作用及びファンデルワールス１−４対相互作用を変調する。
【００８７】
同様に、結合形成変換の７つの成分は、以下の項を含む：λｓｂｓＢ、これは、最終状態に存在する（結合は形成される）が、初期状態に存在しない（結合はまだ形成されていない）結合伸縮相互作用を変調する；λｂａＢ及びλｂｄＢ、これらはそれぞれ、最終状態に存在するが、初期状態に存在しない結合角相互作用及び結合二面角相互作用を変調する；
【００８８】
【数３４】
0006298167.tif 000035
及び
【００８９】
【数３５】
0006298167.tif 000036
、これらは、初期状態では存在するが、最終状態での結合の形成に起因して、最終状態では排除される非結合静電相互作用及びファンデルワールス相互作用を変調する；並びに
【００９０】
【数３６】
0006298167.tif 000037
及び
【００９１】
【数３７】
0006298167.tif 000038
、これらは、最終状態では存在するが、初期状態から最終状態への変換中の結合の形成又
は破断に起因して、初期状態では排除される静電１−４対相互作用及びファンデルワールス１−４対相互作用を変調する。
【００９２】
本明細書で使用される場合、結合パラメータλの成分に関連して使用されるとき、「変調」という用語は、その特定の成分の相互作用エネルギーが、従来の力場のパラメータと、対応する成分結合パラメータとを使用して計算されることを意味する。具体的には、ソフトコアポテンシャルが使用される結合伸縮相互作用及びＬＦ相互作用（又は静電相互作用）では、相互作用エネルギーは、式９及び式２６に従ってそれぞれ計算され、一方、他のタイプの相互作用では、相互作用エネルギーは、従来の力場のパラメータを対応する成分結合パラメータλで乗算したものを使用して計算される。
【００９３】
上述したλスケジュールの一例を以下のスキーム１に示し、スキーム１では、上付き文字（０，１，．．，ｍ，ｍ＋１，．．．，ｎ）は、基準系状態、中間状態、及び標的系状態のインデックスを示し、
【００９４】
【数３８】
0006298167.tif 000039
は、インデックス「ｉ」及び成分「ｃｏｍｐ」を有する状態の各λ成分値である。
【００９５】
【数３９】
0006298167.tif 000040
スキーム１
スキーム１に示されるλｓｂｓＡは、結合破断（例えば、開環）変換にわたり１から０に変化することができる。幾つかの実施形態では、λｓｂｓＡの変化は、変換にわたり線形及び／又は単調であることができる。他の実施形態では、λｓｂｓＡの変化は、変換にわたり非線形及び／又は非単調であることができる。ｍでインデックス付けられる中間系状態ではλｓｂｓＡ＝０．５であることが示されるが、λｓｂｓＡのこのスケジュールは単なる例示であり、非限定的である（例えば、ｍでインデックス付けられる中間系状態に
、０．５よりも小さいか、又は大きい他の値を使用することもできる）。
【００９６】
スキーム１に示されるλｂａＡは、結合破断（例えば、開環）変換にわたり１から０に変化することができる。幾つかの実施形態では、λｂａＡの変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、λｂａＡの変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。
【００９７】
スキーム１に示されるλｂｄＡは、結合破断（例えば、開環）変換にわたり１から０に変化することができる。幾つかの実施形態では、λｂｄＡの変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、λｂｄＡの変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。
【００９８】
更に、特定の実施形態では、分子シミュレーションでのサンプリング効率を改善するために、２個の原子間の結合角相互作用及び結合二面角相互作用は、結合破断変換中、結合伸縮相互作用がオフになる前に、より素早くオフにされて、０になることができる（逆に、結合角相互作用及び結合二面角相互作用は、結合伸縮相互作用が所定の程度までオンになった後でのみ、オンにすることができる）。λｂａＡ／λｂｄＡが、ｍでインデックス付けられた中間系状態で０まで低減可能である（λｂａＡ／λｂｄＡが、λｓｂｓＡの低減よりも素早く低減することにより、初期状態での１からこの中間状態での０に変化可能なことを意味する）ことが示されるが、λｂａＡ／λｂｄＡのスケジュールは単なる例示であり、非限定的である（例えば、λｂａＡ／λｂｄＡは、初期状態からより素早く又はよりゆっくりと０まで低減し得る）。また、λｂａＡ及びλｂｄＡは、それぞれのスケジュールに従って別個に変化することもでき、同期される必要はない。
【００９９】
スキーム１に示される
【０１００】
【数４０】
0006298167.tif 000041
は、結合破断（例えば、開環）変換にわたり０から１に変化することができる。幾つかの実施形態では、
【０１０１】
【数４１】
0006298167.tif 000042
の変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、
【０１０２】
【数４２】
0006298167.tif 000043
の変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。ｍでインデックス付けされる中間系状態において
【０１０３】
【数４３】
0006298167.tif 000044
である（すなわち、
【０１０４】
【数４４】
0006298167.tif 000045
は、初期状態から、ｍでインデックス付けされる中間状態を通して０であることができ、この中間状態から最終状態まで増大することができる）ことが示されるが、このスケジュールは単なる例示であり、非限定的である（例えば、
【０１０５】
【数４５】
0006298167.tif 000046
は、ｍでインデックス付けられる中間状態よりも前又は後の中間状態まで０に保つことができる）。
【０１０６】
スキーム１に示される
【０１０７】
【数４６】
0006298167.tif 000047
は、結合破断（例えば、開環）変換にわたり０から１に変化することができる。幾つかの実施形態では、
【０１０８】
【数４７】
0006298167.tif 000048
の変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、
【０１０９】
【数４８】
0006298167.tif 000049
の変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。ｍでインデックス付けされる中間系状態において
【０１１０】
【数４９】
0006298167.tif 000050
である（
【０１１１】
【数５０】
0006298167.tif 000051
は、値０である初期状態から、値１におけるｍでインデックス付けされる中間状態まで増大し、最終状態まで１であり続けることができる）ことが示されるが、このスケジュールは単なる例示であり、非限定的である（例えば、
【０１１２】
【数５１】
0006298167.tif 000052
は、ｍでインデックス付けられる中間状態よりも前又は後の中間状態において１まで増大することができる）。更に、特定の実施形態では、分子シミュレーションでのサンプリング効率を改善するために、２個の原子間のファンデルワールス相互作用は、変換中、静電相互作用がオンになる前に、完全にオンにすることができる。これはスケジュール１に示されており、
【０１１３】
【数５２】
0006298167.tif 000053
及び
【０１１４】
【数５３】
0006298167.tif 000054
のスケジュールは、
【０１１５】
【数５４】
0006298167.tif 000055
が、中間状態で１よりも小さい場合、
【０１１６】
【数５５】
0006298167.tif 000056
がその中間状態で０であるように調整される。
【０１１７】
スキーム１に示される
【０１１８】
【数５６】
0006298167.tif 000057
は、結合破断（例えば、開環）変換にわたり１から０に変化することができる。幾つかの実施形態では、
【０１１９】
【数５７】
0006298167.tif 000058
の変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、
【０１２０】
【数５８】
0006298167.tif 000059
の変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。ｍでインデックス付けられた中間系状態において
【０１２１】
【数５９】
0006298167.tif 000060
である（すなわち、
【０１２２】
【数６０】
0006298167.tif 000061
が、初期状態における１から、λｓｂｓＡの低減よりも素早い低減によりこの中間状態において０に変化することができる）ことが示されるが、このスケジュールは単なる例示であり、非限定的である（例えば、
【０１２３】
【数６１】
0006298167.tif 000062
は、初期状態からより素早く又はよりゆっくりと０まで低減し得る）。
【０１２４】
スキーム１に示される
【０１２５】
【数６２】
0006298167.tif 000063
は、結合破断（例えば、開環）変換にわたり１から０に変化することができる。幾つかの実施形態では、
【０１２６】
【数６３】
0006298167.tif 000064
の変化は、変換又は変換の一部にわたり線形及び／又は単調であることができる。他の実施形態では、
【０１２７】
【数６４】
0006298167.tif 000065
の変化は、変換又は変換の一部にわたり非線形及び／又は非単調であることができる。ｍでインデックス付けられた中間系状態において
【０１２８】
【数６５】
0006298167.tif 000066
である（
【０１２９】
【数６６】
0006298167.tif 000067
が、初期状態からｍでインデックス付けられる中間状態を通して１に保たれ、次に、その中間状態から最終状態まで０に低減することができる）ことが示されるが、このスケジュールは単なる例示であり、非限定的である（例えば、
【０１３０】
【数６７】
0006298167.tif 000068
は、初期状態から、ｍでインデックス付けられる中間状態の前又は後の中間状態を通して１に保たれ、次に、その中間状態から最終状態まで０に低減することができる）。更に、特定の実施形態では、分子シミュレーションでのサンプリング効率を改善するために、２個の原子間のファンデルワールス相互作用は、変換中、静電相互作用が完全にオフになった後でのみ、オフにされる。これはスケジュール１に示されており、スケジュール１では、
【０１３１】
【数６８】
0006298167.tif 000069
及び
【０１３２】
【数６９】
0006298167.tif 000070
のスケジュールは、
【０１３３】
【数７０】
0006298167.tif 000071
が中間状態で１よりも小さい場合、
【０１３４】
【数７１】
0006298167.tif 000072
が、その中間状態で０であるように調整される。
【０１３５】
結合形成変換のλ成分λｓｂｓＢ、λｂａＢ及びλｂｄＢ、
【０１３６】
【数７２】
0006298167.tif 000073
及び
【０１３７】
【数７３】
0006298167.tif 000074
、
【０１３８】
【数７４】
0006298167.tif 000075
及び
【０１３９】
【数７５】
0006298167.tif 000076
のスケジュールは、結合破断変換のλスケジュールの上記説明に鑑みて当業者によって容易に理解することができ、その理由は、２つの変換プロセスが基本的に、互いを逆にしたものであるためである。しかし、結合形成でのλ成分及び結合破断変換での対応するλ成分（例えば、λｓｂｓＡ及びλｓｂｓＢ）の値が、合計して１になるように見えるが、これは単なる例示であることに留意する。任意の所与の系で、結合は、本願で考慮される同じ錬金術的変換において同時に形成及び破断の両方を行うことはできない。したがって、結合破断変換でのλ成分スケジュールは、結合形成変換でのλ成分スケジュールと調整される必要はない。例として、λｓｂｓＢは、０から１に非線形的に増大するスケジュールを辿ることができ、一方、λｓｂｓＡは、１から０に線形的に低減するスケジュールを辿る。更に、図１ｃに示されるように、本願は、１つの結合（例えば、Ａ１とＡ４との結合）の破断と、別の（１つ又は複数の）結合（例えば、Ａ１とＡ１３との結合及びＡ４とＡ１３との結合）の同時形成を意図する。この場合、結合破断変換のλ成分スケジュール及び結合形成変換のλ成分スケジュールは、別個に選択することができ、互いに依存しない。
【０１４０】
以下、図１ａに示される変換の分子シミュレーションでのエネルギー計算が、上述したλスケジュールの特定の実施形態に関連して示される。
・図１ａ中のＡ１とＡ３との結合形成による影響を受ける結合伸縮相互作用：（Ａ１，Ａ３）、
・図１ａ中のＡ１とＡ３との結合形成による影響を受ける結合角相互作用：
（ｉ）初期状態のみに見られ、最終状態には見られない結合角相互作用：（Ａ２，Ａ１，Ａ１０）、（Ａ４，Ａ１，Ａ１０）、（Ａ５，Ａ１，Ａ１０）、（Ａ２，Ａ３，Ａ１１）、（Ａ８，Ａ３，Ａ１１）、（Ａ９，Ａ３，Ａ１１）、
（ｉｉ）最終状態のみに見られ、初期状態には見られない結合角相互作用：（Ａ２，Ａ１，Ａ３）、（Ａ４，Ａ１，Ａ３）、（Ａ５，Ａ１，Ａ３）、（Ａ１，Ａ３，Ａ２）、（Ａ１，Ａ３，Ａ８）、（Ａ１，Ａ３，Ａ９）、
・図１ａ中のＡ１とＡ３との結合形成による影響を受ける結合二面角相互作用：
（ｉ）初期状態のみに見られ、最終状態には見られない原子が関わる結合二面角相互作用：（Ａ１０，Ａ１，Ａ２，Ａ３）、（Ａ１０，Ａ１，Ａ２，Ａ６）、（Ａ１０，Ａ１，Ａ２，Ａ７）、（Ａ１，Ａ２，Ａ３，Ａ１１）、（Ａ６，Ａ２，Ａ３，Ａ１１）、（Ａ７，Ａ２，Ａ３，Ａ１１）、
（ｉｉ）最終状態のみに見られ、初期状態には見られない原子が関わる結合二面角相互作用：（Ａ３，Ａ１，Ａ２，Ａ６）、（Ａ３，Ａ１，Ａ２，Ａ７）、（Ａ２，Ａ１，Ａ３，Ａ８）、（Ａ２，Ａ１，Ａ３，Ａ９）、（Ａ４，Ａ１，Ａ３，Ａ２）、（Ａ４，Ａ１，Ａ３，Ａ８）、（Ａ４，Ａ１，Ａ３，Ａ９）、（Ａ５，Ａ１，Ａ３，Ａ２）、（Ａ５，Ａ１，Ａ３，Ａ８）、（Ａ５，Ａ１，Ａ３，Ａ９）、（Ａ６，Ａ２，Ａ３，Ａ１）、（Ａ７，Ａ２，Ａ３，Ａ１）、
・図１ａ中の、初期状態では非結合相互作用の計算について排除されるが、最終状態では排除されない原子対：（Ａ１，Ａ１０）、（Ａ１，Ａ１１）、（Ａ２，Ａ１０）、（Ａ２，Ａ１１）、（Ａ３，Ａ１０）、（Ａ３，Ａ１１）、（Ａ４，Ａ１０）、（Ａ５，Ａ１０）、（Ａ６，Ａ１０）、（Ａ６，Ａ１１）、（Ａ７，Ａ１０）、（Ａ７，Ａ１１）、（Ａ８，Ａ１１）、（Ａ９，Ａ１１）、
・図１ａ中の、最終状態では非結合相互作用の計算について排除されるが、初期状態では排除されない原子対：（Ａ４，Ａ８）、（Ａ４，Ａ９）、（Ａ５，Ａ８）、（Ａ５，Ａ９）、
・図１ａ中の、初期状態では含まれるが、最終状態では含まれない１−４原子対：（Ａ１，Ａ８）、（Ａ１，Ａ９）、（Ａ３，Ａ４）、（Ａ３，Ａ５）、
・図１ａ中の、最終状態では含まれるが、初期状態では含まれない１−４原子対：（Ａ４，Ａ８）、（Ａ４，Ａ９）、（Ａ５，Ａ８）、（Ａ５，Ａ９）。
【０１４１】
開環又は閉環変換での結合の形成及び消滅による影響を受けないその他の結合相互作用及び非結合相互作用は、開環及び閉環に関わらない従来の自由エネルギー摂動のおけるように、通常のλスケジュールによって変調することができる（例えば、初期状態から最終状態までのλｏｔｈｅｒの増分であり、ＰＡＢとＰＡとの間の相互作用は１−λｏｔｈｅｒによってスケーリングされ、ＰＡＢとＰＢとの間の相互作用はλｏｔｈｅｒによってスケーリングされる）。更に、初期状態のみに見られるが、最終状態からはなくなっている原子（すなわち、最終状態では「ダミー」原子になる原子）が関わる全ての相互作用も、通常のλスケジュールによって取り扱われる。逆に、逆（結合形成）変換では、最終状態のみに見られ、初期状態にはない原子（例えば、図１ａ中のＡ１１及びＡ１２）が関わる全ての相互作用も、通常のλスケジュールによって取り扱われる。残りの相互作用、すなわち、初期状態のみに見られるか、又は最終状態のみで見られる原子が関わらない開環又は閉環変換での結合の形成及び消滅により影響を受ける原子は、スケジュール１に記載されるものと同様の「特別な」すなわち、図１ａに示される閉環変換でのλスケジュールによって変調することができる。
【０１４２】
・初期状態では含まれるが、最終状態では含まれない結合伸縮相互作用は、本明細書に記載される弱結合相互作用ポテンシャルを使用して、λｓｂｓＡによって変調され（この例では、いずれの相互作用もこのカテゴリに属さない）、
・初期状態では含まれないが、最終状態では含まれる結合伸縮相互作用は、本明細書に記載される弱結合相互作用ポテンシャルを使用して、λｓｂｓＢによって変調され（この例では、原子Ａ１とＡ３との間の結合伸縮相互作用が、このカテゴリに属する）、
・初期状態では含まれるが、最終状態では含まれない結合角相互作用は、λｂａＡによって変調され（この例では、いずれの相互作用もこのカテゴリに属さない）、
・最終状態では含まれるが、初期状態では含まれない結合角相互作用は、λｂａＢによって変調され（この例では、原子（Ａ２，Ａ１，Ａ３）、（Ａ４，Ａ１，Ａ３）、（Ａ５，Ａ１，Ａ３）、（Ａ１，Ａ３，Ａ２）、（Ａ１，Ａ３，Ａ８）、（Ａ１，Ａ３，Ａ９）間の結合角相互作用が、このカテゴリに属する）、
・初期状態では含まれるが、最終状態では含まれない結合二面角相互作用は、λｂｄＡによって変調され（この例では、いずれの相互作用もこのカテゴリに属さない）、
・最終状態では含まれるが、初期状態では含まれない結合二面角相互作用は、λｂｄＢによって変調され（この例では、原子（Ａ３，Ａ１，Ａ２，Ａ６）、（Ａ３，Ａ１，Ａ２
，Ａ７）、（Ａ２，Ａ１，Ａ３，Ａ８）、（Ａ２，Ａ１，Ａ３，Ａ９）、（Ａ４，Ａ１，Ａ３，Ａ２）、（Ａ４，Ａ１，Ａ３，Ａ８）、（Ａ４，Ａ１，Ａ３，Ａ９）、（Ａ５，Ａ１，Ａ３，Ａ２）、（Ａ５，Ａ１，Ａ３，Ａ８）、（Ａ５，Ａ１，Ａ３，Ａ９）、（Ａ６，Ａ２，Ａ３，Ａ１）、（Ａ７，Ａ２，Ａ３，Ａ１）間の結合二面角相互作用が、このカテゴリに属する）、
・初期状態では排除されるが、最終状態では排除されない相互作用は、λｅｌｅｃＡｅｘ及びλｖｄｗＡｅｘによって変調され（この例では、いずれの相互作用もこのカテゴリに属さない）、
・最終状態では排除されるが、初期状態では排除されない相互作用は、λｅｌｅｃＢｅｘ及びλｖｄｗＢｅｘによって変調され（この例では、対（Ａ４，Ａ８）、（Ａ４，Ａ９）、（Ａ５，Ａ８）、（Ａ５，Ａ９）の最終状態において排除される相互作用が、このカテゴリに属する）、
・初期状態では含まれるが、最終状態では含まれない１−４対相互作用は、λｅｌｅｃＡ１４及びλｖｄｗＡ１４によって変調され（この例では、対（Ａ１，Ａ８）、（Ａ１，Ａ９）（Ａ３，Ａ４）（Ａ３，Ａ５）の１−４対相互作用が、このカテゴリに属する）、
・最終状態では含まれるが、初期状態では含まれない１−４対相互作用は、λｅｌｅｃＢ１４及びλｖｄｗＢ１４によって変調される（この例では、対（Ａ４，Ａ８）、（Ａ４，Ａ９）（Ａ５，Ａ８）（Ａ５，Ａ９）の１−４対相互作用が、このカテゴリに属する）。
【０１４３】
初期状態から最終状態までの変換での全ての系状態に適するλスケジュール（スキーム１に示されるλスケジュール等）によって定義されるエネルギー項を用いて、分子シミュレーションが実行されて、λスケジュールに従って基準状態、標的状態、及び中間状態において得られるマイクロ状態のアンサンブルをサンプリングすることができる。λ窓毎に自由エネルギー差が、全ての隣り合うラムダ窓間ΔＦλｉ→λｉ＋１及び／又は基準状態と標的状態との間を含め、ラムダ窓の任意の対間ΔＦｉ→ｊで計算することができる。基準状態と標的状態との間の合計自由エネルギー差は、λスケジュールによって定義される変換路に沿った２つの連続した各状態間の自由エネルギー差を加算することによって得ることができるか、又はサンプリングされた全ての状態からのデータを解析することによって直接得ることができる。
【０１４４】
隣り合うλ窓間又は一般的に、初期状態と最終状態とを含む状態の任意の対間の自由エネルギー差は、様々な方法で計算することができる。例えば、標的状態、基準状態、及び中間状態で得られるマイクロ状態の適するアンサンブルの、例えば内部エネルギー差（自由エネルギー摂動  ＮＶＴアンサンブル）、エンタルピー差（自由エネルギー摂動  ＮＰＴアンサンブル）、又は他の関連する熱力学的特性差（ＮＶＥアンサンブル等の自由エネルギー摂動の他のアンサンブル）の使用により、結合パラメータλとしては、マイクロ状態の選択されたアンサンブルに対して瞬時に変更される。解析は、例えば、バネット受容比（ＢＡＲ：Ｂａｎｎｅｔ  Ａｃｃｅｐｔａｎｃｅ  Ｒａｔｉｏ）、多状態バネット受容比（ＭＢＡＲ：Ｍｕｌｔｉｓｔａｔｅ  Ｂａｎｎｅｔ  Ａｃｃｅｐｔａｎｃｅ  Ｒａｔｉｏ）、加重ヒストグラム解析法（ＷＨＡＭ：Ｗｅｉｇｈｔｅｄ  Ｈｉｓｔｏｇｒａｍ  Ａｎａｌｙｓｉｓ  Ｍｅｔｈｏｄ）、ツヴァンツィヒ平均、又は他の同様の自由エネルギー摂動ファミリ推定量によって更に実行することができる。代替的には、隣り合うλ窓間又は一般には初期状態と最終状態とを含む状態の任意の対間の自由エネルギー差は、得られる微視的状態毎に、結合ベクトルλに対するエネルギーの導関数（ＴＩ  ＮＶＴアンサンブル）、結合ベクトルλに対するエンタルピーの導関数（ＴＩ  ＮＰＴアンサンブル）、又は結合ベクトルλに対する他の関連する熱力学的特性の導関数（ＴＩ  ＮＶＥアンサンブル等のＴＩの他のアンサンブル）解析により計算することができる。他の実施形態では、各λ窓の自由エネルギー差は、例えば、基準状態、標的状態、及び中間状態間で動的に、λ動力学を介して一般性を失わずに遷移することができる動的変数として、結合ベクト
ルλのサンプリングに関連付けられた平均力ポテンシャル（ＰＭＦ：ｐｏｔｅｎｔｉａｌ
ｏｆ  ｍｅａｎ  ｆｏｒｃｅ）の解析により、計算することができ、λ動力学の原理は一般に、ナイト（Ｋｎｉｇｈｔ）ら、λ−動力学自由エネルギーシミュレーション法（λ−ｄｙｎａｍｉｃｓ  ｆｒｅｅ  ｅｎｅｒｇｙ  ｓｉｍｕｌａｔｉｏｎ  ｍｅｔｈｏｄｓ）、ジャーナル・オブ・コンピューテイショナル・ケミストリ（Ｊ．Ｃｏｍｐｕｔ．Ｃｈｅｍ．）、２００９年、第３０巻、ｐ．１６９２〜１７００に記載されており、この開示は参照により本明細書に援用される。λ動力学に基づくサンプリング方法は、λ−モンテカルロ、λ−メタダイナミクス、λ−ＯＳＲＷ、及び他のλ  ＰＭＦサンプリングファミリ法を含むが、これらに限定されない。
【０１４５】
本明細書に記載される自由エネルギー差計算の方法は、幾つかのかなり有用な用途に適用することができ、それらとしては、例えば、
・開環又は閉環を有する同種リガンド間の相対タンパク質−リガンド結合親和性及び／又は相対溶媒自由エネルギー計算（図４参照）、
・大環状化が異なる同種リガンド間の相対タンパク質−リガンド結合親和性及び／又は相対溶媒自由エネルギー計算（図５参照）、
・非プロリン残基からプロリン残基への突然変異又はプロリン残基から非プロリン残基への突然変異のタンパク質熱動力学安定性、タンパク質−リガンド結合親和性、又はタンパク質−タンパク質結合親和性への影響の計算（図６参照）、及び
・残基挿入又は残基削除のタンパク質熱動力学安定性、タンパク質−リガンド結合親和性、又はタンパク質−タンパク質結合親和性への影響の計算（図７参照、図７は、１２９と付番された残基を失うことにより、左側に示される例示的なタンパク質セグメント構造が、右側に示される構造になることを概略的に示す）
が挙げられる。
【０１４６】
開示される主題の自由エネルギー計算の方法の実施形態は、コンピュータプログラムで実施することができ、コンピュータプログラムは、適するハードウェアプラットフォーム、例えば、スタンドアロンコンピュータ、ネットワーク化されたコンピュータ、ネットワークサーバコンピュータ、ハンドヘルドデバイス等のソフトウェア構成要素の形態を取ることができる。開示される方法の様々な態様は、様々なソフトウェアモジュールで実施し得、ソフトウェアがいかに設計されるかに応じて、１つのプロセッサ又は様々なプロセッサによって順次又は並列して実行し得る。プログラムを実行することができる装置は、１つ又は複数のプロセッサ、１つ又は複数のメモリデバイス（ＲＯＭ、ＲＡＭ、フラッシュメモリ、ハードドライブ、光ドライブ等）、入力／出力デバイス、ネットワークインターフェース、及び他の周辺デバイスを含むことができる。プログラムを記憶するコンピュータ可読非一時的媒体も提供される。
【０１４７】
開示される主題の範囲は、本明細書に記載される特定の実施形態によって限定されるべきではない。実際には、本明細書に記載されるものに加えて、開示される主題の様々な変更形態が、上記説明及び添付図から当業者に明らかになるであろう。したがって、そのような変更形態は、添付の特許請求の範囲内にある。
【図１ａ】
0006298167.tif 000078
【図１ｂ】
0006298167.tif 000079
【図１ｃ】
0006298167.tif 000080
【図２】
0006298167.tif 000081
【図３】
0006298167.tif 000082
【図４】
0006298167.tif 000083
【図５】
0006298167.tif 000084
【図６】
0006298167.tif 000085
【図７】
0006298167.tif 000086