WO2005029385A1

WO2005029385A1 - 分子シミュレーション方法及び装置

Info

Publication number: WO2005029385A1
Application number: PCT/JP2004/013808
Authority: WO
Inventors: Yasushige Yonezawa; Toshikazu Takada; Kazuto Nakata; Toshihiro Sakuma; Haruki Nakamura
Original assignee: Nec Corporation; Osaka University
Priority date: 2003-09-22
Filing date: 2004-09-22
Publication date: 2005-03-31
Also published as: US20070043545A1; JPWO2005029385A1

Abstract

　シミュレーション対象の分子または分子の一部をＱＭ空間とＭＭ空間とに分割し、ＱＭ空間に対して非経験的分子軌道法を適用し、ＭＭ空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレーション方法は、記憶部から、シミュレーション対象の分子または分子の一部を構造データを取り出してＱＭ空間及びＭＭ空間に分割する段階と、ＱＭ空間に関する非経験的分子軌道法における全エネルギー表式の一部を経験的ポテンシャルで置き換える段階と、を有する。

Description

明細書

分子シミュレーション方法及び装置

技術分野

[0001] 本発明は、量子化学の手法によって分子シミュレーションを行うための方法及び装置に関し、特に、量子化学の理論手法のうち、非経験的分子軌道法と分子力場法を融合してひとつの理論体系として扱う QMZMM(Quantum Mechanics/Molecular Mechanics)法による分子シミュレーション方法及び装置に関する。

背景技術

[0002] 量子化学理論の発展や計算機技術の進歩によって、計算により、化合物分子の構造、物性や、分子内の化学結合や分子軌道、電子状態などを精度よくシミュレーションできるようになつてきた。そのような手法の中でも経験的パラメータに原則として依存しない非経験的分子軌道法は、過去数十年にわたる数多くの計算例から、その理論的枠組の正しさが実証されてきた。近年、非経験的分子軌道法の生体分子への適応が、創薬や機能性食品の基礎研究の観点力求められている。例えば、特定のタンパク質の活性部位に高度に選択的に結合するリガンドを発見することは、特定疾患に有効な薬剤の探索の上で極めて有効である。このようなリガンドは、そのタンパク質の活性を阻害するという意味では、インヒビターとして機能する。そこで、リガンドの候補となるべき極めて多数の化合物から実際の候補ィ匕合物を探索するスクリーニングが行われることになるが、その場合、実際に化学反応を起こさせることなぐ計算機シミュレーションによって候補ィ匕合物をスクリーニングすることができれば、全体のスクリー-ングに要する時間の大幅な短縮が望めるほか、大量合成が困難なタンパク質に対するリガンドの探索や、未だに化学合成例のないリガンドを対象にすることができるなど、大きな利点を得ることができる。そこで、非経験的分子軌道法は、計算機によるこのような分子シミュレーションを行うために利用されるようになってきている。

[0003] し力しながら、非経験的分子軌道法では、その全電子計算に必要とされる膨大な計算量から、計算できる分子の大きさは依然として中規模の分子に限定される状況にある。例えば、タンパク質とリガンドの複合体の非経験的分子軌道計算は、上述したように薬候補ィ匕合物の探索に重要な情報をもたらすが、複合体の全系の電子状態計算は、最新のスーパーコンピュータをにもってしても長時間を要するものであって、企業における研究開発の時間スケールでは一般に容認できるものではない。

[0004] ところで一般に、タンパク質では、それが非常にアミノ酸残基数が多いものであっても、リガンドなどと結合する部位は限られた場所であり、したがって、電子状態についての正確なシミュレーションはその部位及びその近傍にっ、てだけ行えばょ、と考えられている。そこで、多数の分子力もなる化学系において、分子もしくは分子の一部を、注目している化学現象の起こる QM(Quantum Mechanics)空間と、それ以外の 2 次的な MM(Molecular Mechanics)空間とに分割し、前者を非経験的分子軌道法など量子力学的取り扱いにより、後者を分子力場法などの経験的ポテンシャルで記述する、 QMZMM法が提唱されている [1]。この方法の利点は：

(1)注目している QM空間での化学現象は、量子力学に基づく方法により記述されるので、計算結果について高い信頼性が期待され；

(2)分子系の大部分を占める MM空間は、計算負荷の少ない経験的ポテンシャルにより記述され、古典力学的取扱いで処理されるので、計算時間が大幅に短縮される；

ということである。 QM/MM法は、多数提案されている分子シミュレーションの理論の中でも、計算結果の信頼性維持と計算時間の短縮という、分子シミュレーションに対する相反する要請を満たす技術として、ここ数年、注目されている。

[0005] このように、従来の全電子計算に対する計算時間の短縮力 QMZMM法の特徴である。 500アミノ酸残基力もなるタンパク質で計算時間を見積もると、全計算時間の 95%以上を QM空間に対する計算に消費されている。 500アミノ酸残基の全体に対して QM法を適用した場合には、その計算には、 5万軌道の程度の原子軌道が必要である。計算時間は、用いる原子軌道の数の 3乗に比例して増減する。ここで QMZ MM法を適用し、活性部位近傍のみを QM空間として扱う場合、計算しなければならない原子軌道の数は 5千軌道程度になり、 1Z1000に計算時間が短縮される。特に、タンパク質 +基質複合体の安定構造や酵素の反応機構の計算には、このような計算を繰り返すことになるので、 QMZMM法以外の手法によって意味のある分子シミユレーシヨンを行おうとすることは、現実的でない。

[0006] 図 1は、タンパク質リガンド複合体における QM空間と MM空間との分割を概念的に示す図である。この図で示すように、リガンド 1がタンパク質 2に結合している活性部位の近傍を QM空間 21とし、それ以外を MM空間 22としている。

[0007] しかしながら、 QMZMM法では、 QM空間に属する原子と MM空間に属する原子との間の相互作用の記述方法に課題が残っている。 QM空間に属する原子と MM空間に属する原子が相互に共有結合を直接形成していない場合には、それらの相互作用をクーロン力やファンデルワールス力などで表現することにより、原理的な困難を生じることなぐ QM空間と MM空間との接続を取り扱うことが可能である。しかしながら、共有結合を形成している場合、その一方の原子が QM空間に属し、他方が M M空間に属することになり、深刻な問題が発生する。例えば，タンパク質の活性部位近傍を QM空間にする場合、共有結合であるタンパク質の主鎖を QM— MM接合面が横切ることになり、隣接する原子間のいわゆる σ (シグマ)結合が人為的に切断されること〖こなる。

[0008] 図 2は、 QM空間と MM空間との分割を説明する図であって、タンパク質におけるペプチド結合の一例を示している。ここでは図示破線で示すように、 N末端（図示左端)力数えて 2番目のアミノ酸残基内にぉ、て主鎖の C N結合の位置で、 QM空間と MM空間とを分けて、る。この QMZMM法に特徴的な空間分割では：

1)共有結合に関与している電子の一方を QM空間の属する原子に帰属させるため、不対電子が人為的にその原子に発生し、その原子がラジカルになってしまう。その結果、系の化学的性質が大きく歪められ、 QM空間の波動関数の正しい記述が得られないことになる；

2) MM空間に属する原子は、電荷をもつ質点として表現されるため、連続的な電子分布をもつ QM空間と質点からなる不連続な MM空間が、ひとつの化学結合を境にして相互作用することになり、系の滑らかで連続的なポテンシャル関数が得られない。このため、 QM空間と MM空間との境界領域近傍の原子の挙動が正しく記述されない危険性がある；

という、理論上の問題を提起する。 [0009] 前者につ!、ては、リンクアトム (Link Atom)法と呼ばれる、不対電子に対して共有結合で結合する水素類似原子を人為的に導入する方法が、伝統的に用いられている [2]。し力しながら、リンクアトム法では、元々の系には存在しない水素類似原子を導入するため、余計なエネルギー項が出現し、その補正が必要となるという問題がある。また、導入された水素類似原子自体が、波動関数の性質に変化をもたらす可能性もある。最近では、 MM空間に帰属された原子を、水素類似原子に置き直す工夫もされて!/、るが、後者の問題の解決策とはなってヽな、。

[0010] さらには、 Bond Orbitalと呼ばれる局在した軌道を共有結合にあてがい、不対電子の発生を防ぐ方法も提案されている [3]。しかしながら、化学結合の軸方向に固定された軌道を用いるため、原子核の移動と共に軌道の回転などの補正が必要となり、理論が複雑になるという難点がある。

[0011] QM空間と MM空間との接合部近傍のポテンシャルが連続的でないという問題は、構造最適化もしくは分子動力学による時間発展の計算を行った時、顕著になる。すなわち、近傍の原子核が、計算上、本来の位置力大きくずれて、例えば、タンパク質の活性部位の構造を正しく維持できなくなる可能性がある。

[0012] この課題を解決する方法として、 ONIOM(our own n- layered integrated molecular orbital + molecular mechanics method)が提唱されている [4]。この方法では、系全体を MM領域としてー且扱うことで、境界領域でのポテンシャルの不連続性の問題は回避しているものの、 QM空間の原子群を全く独立した系として別途計算するため、 MM空間からの影響が QM空間に自己無撞着的に反映されないという問題が指摘されている [3]。

[0013] 以下、本明細書で引用した参考文献を列挙する。

非特干文献 1： [1」 J. Gao, Methods and Applications of Combined Quantum Mechanical and Molecular Mechanical Potentials in "Reviews in Computational Chemistry , Vol. 7, K. B. Lipkowitz and D. B. Boyd, Editors, VCH Publishers, Inc. New York, 1996

非特許文献 2 : [2] M. J. Field, P. A.Bash and M. Karplus, J. Comp. Chem., Vol. 11, 700 (1990) 非特許文献 3 : [3] D. M. Philipp, R. A. Friesner, J. Comp. Chem., Vo; 20, 1468 (1999)

非特許文献 4 : [4] T. Vreven, Κ. Morokuma, J. Comp. Chem., Vol. 21, 1419 (2000) 非特許文献 5 : [5]米澤康滋、高田俊和、佐久間俊広、中田一人、中村春木、生物物理、第 43卷、 Supplement 1, B198 (日本生物物理学会第 41会年会講演予稿集) (2003)

非特許文献 6 : [6]電子構造の理論入門:新しい量子化学 (上)、 pl55、東京大学出版社、 1987年

特干文献 7： [7」 P. Puiay, Direct Use of uradient for investigating Molecular Energy Surfaces , in "Modern Theoretical Chemistry 4〃, H. F. Scharfer, Editor, Plenum Press, New York and London, 1977

非特許文献 8 : [8] J. A. Pople, R. Krishnan, H. B. Schlegel and J. S. Binkley, Int. J. Quant. Chem., Vol. 13, 225 (1979)

発明の開示

発明が解決しょうとする課題

[0014] QMZMM法を、生体分子などに対する実用的な分子シミュレーション技術にするためには、前述したように、 1)不対電子による人為的な波動関数の歪曲、 2)不連続なポテンシャルによる分子構造の不安定性、の課題を解決しなければならな、。

[0015] MM法で用いる経験ポテンシャル E(MM)は、（1)式で示すように、結合距離 E 、

bonds 結合角 E 、2面角E 、非共有結合 E で与えられる。

angles torsions nonbonds

[0016] [数 1]

E(MM) = E_bonds + E_angles + E_lorsions + E_nonbonds

ここで、 Ε^ΜΜ(0)は、平衡状態での全 MMエネルギーを表す。このため、全系を経験的ポテンシャルで記述すれば、異常な分子骨格が計算中に出現することはない。他方、 QM空間自体は、量子力学的な手法により厳密に記述されているので、 QM空間の分子構造が異常をきたすことは基本的にありえない。この議論から、 QM空間の内側で、 MM空間に接合する表層部分に、理論的に矛盾なぐ MM空間から経験的ポテンシャルを染み込ませることができれば、分子構造の不安定性の問題を克服できることが分かる。その時、 QM空間の末端の原子に発生する不対電子を、エネルギー表式に基づいて、論理的に隠蔽できることも、解決のための必須の条件である。

[0017] 本発明者らは、 QM領域と MM領域の境界に局所化された分子軌道を配置し、この局所軌道を QM領域に対して Frozen Orbitalとして計算に取り入れることによって、 QM領域と MM領域とのスムースかつ精密な接続を実現できることを既に予想している [5]。

[0018] そこで本発明の目的は、不対電子を隠蔽しつつ、 MM空間に接合する QM空間の表層に、（1)式で与えられる経験的ポテンシャルを染み込ませることのできる、論理的に整合性の取れた分子シミュレーション方法及び装置を提供することにある。

課題を解決するための手段

[0019] 本発明の第 1の様相に従えば、分子シミュレーション方法は、シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、 QM空間に対して非経験的分子軌道法を適用し、 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレーション方法であって、記憶部から、シミュレーション対象の分子または分子の一部を構造データを取り出して QM空間及び MM空間に分割する段階と、 QM空間に関する非経験的分子軌道法における全エネルギー表式の一部を経験的ポテンシャルで置き換える段階と、を有する。

[0020] 本発明の第 2の様相に従えば、分子シミュレーション方法は、シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、 QM空間に対して非経験的分子軌道法を適用し、 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレーション方法であって、記憶部から、シミュレーション対象の分子または分子の一部を構造データを取り出して QM空間及び MM空間に分割し、さらに、 QM空間を、 MM空間に隣接する領域である表層 Q M領域と、表層 QM領域以外の領域である QM領域とに分割する段階と、表層 QM 領域に関して非経験的分子軌道法による全エネルギー表式を得る段階と、 QM領域に関して非経験的分子軌道法による全エネルギー表式を得る段階と、表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換える段階と、

QM空間についての非経験的分子軌道法による全エネルギーを求める段階と、を有する。

[0021] 本発明の第 3の様相に従えば、分子シミュレーション装置は、シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、 QM空間に対して非経験的分子軌道法を適用し、 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレーション装置であって、シミュレーシヨン対象の分子または分子の一部の構造データを格納する記憶部と、記憶部から、シミュレーション対象の分子または分子の一部を構造データを取り出して QM空間及び MM空間に分割し、さらに、 QM空間を、 MM空間に隣接する領域である表層 Q M領域と、表層 QM領域以外の領域である QM領域とに分割する領域分割部と、表層 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換え、 QM空間についての非経験的分子軌道法による全エネルギーを求める第 1の演算部と、を有する。

[0022] このように本発明では、 QMZMM法による分子シミュレーションにお!/、て、 QM法が適用される空間のうち MM空間に接続する部分を表層 QM領域とし、この表層 Q M領域には経験的ポテンシャルが染み込むようにして、る。本発明の手法を用いることにより、シミュレーション計算の過程において、タンパク質の活性部位の近傍の構造、すなわち QM空間と MM空間とが接続する領域において、タンパク質の構造が正しく維持される。このように本発明によれば、 QM空間と MM空間とに分割したことに起因する不整合の問題を回避でき、精度の高い分子シミュレーションを行えるようになる。

図面の簡単な説明

[0023] [図 1]タンパク質リガンド複合体における QM— MM空間分割を概念的に示す図である。

[図 2]QM空間と MM空間との分割の一例を説明する図である。 [図 3]本発明の実施の一形態の分子シミュレーション装置の構成を示すブロック図である。

[図 4]QM空間と MM空間との接合部分の分子配列を拡大して示した図である。

[図 5]分子シミュレーションの手順を示すフローチャートである。

発明を実施するための最良の形態

[0024] 本発明に基づく方法は、 QMZMM法によって分子シミュレーションを行う際に、 Q M法が適用される空間の一部にぉ、て、非経験的分子軌道法における全エネルギ一表式の一部を経験的ポテンシャルに置き換えることにより、従来の QMZMM法にお、て発生する QM-MM境界領域の分子構造不安定性を回避するものである。具体的には、 QMZMM法の適用を前提として、シミュレーション対象の分子を QM法が適用される QM空間と、 MM法が適用される MM空間とに分割し、さらに、 QM空間のうち、 MM空間に近接する部分を表層 QM領域とする。そして、表層 QM領域においては、その全エネルギーの一部を経験的ポテンシャルに置き換えている。経験的ポテンシャルにどの程度置き換えるかは、外部からのパラメータで調整可能とすることが好ましい。さらに、表層 QM領域では、その表層 QM領域の一部を構成する分子もしくは分子の一部の波動関数を、局在化分子軌道で表現することが好ましい。このようにして本実施形態では、 QM空間の一部に経験的ポテンシャルを持ち込んでいる。

[0025] まず、本実施形態に基づく方法の理論的側面を詳細に説明する。

[0026] いま、 QM空間と MM空間が接合している表層領域に分布する分子もしくは分子の部分構造に、局在している局在化分子軌道を考える。他方、表層以外の QM空間の原子に属する電子を記述する分子軌道は、その空間全体に拡がる正準軌道

(Canonical Orbital)で与えられるとする。これらの軌道を用いた、 Hartree-Fock (ハートリーフォック）法の全エネルギー E(QM + PQM)は、（2)式で与えられる。

[0027] [数 2]

E(QM + PQM) =∑2H ₊∑∑(υ, -^) +∑ +∑∑(2 - K_y)

(2) 7 7 7

+∑∑(2 - ) ₊∑∑¾^ ₊∑2H ∑∑₍2 - ) +∑ ただし、

[0028] [数 3]

(

= ffiPy (1) ψί (1)― (pj (2) φ_}- (2) άτ_χ άτ₂

Ky = \\_φί (1) (1)― q>j (2) (2) άτ_χ άτ₂ である。ここで、 Qは正準軌道の空間を、 Pは局在化分子軌道の空間を意味している。後述する説明から明らかになるように、 Qが表す空間は、経験的ポテンシャルが持ち込まれずに純粋に QM法によって処理される領域 (QM領域)であり、 Pが表す空間は、 QM法による計算が行われるものの経験的ポテンシャルが持ち込まれる領域（表層 QM領域)である。また、（2)式は、分子軌道法の分野における通常の表記を採用しており、したがって、 φ .は分子軌道であって、

[0029] 画 r

と表される。 R は原子核 Aと Bとの距離、 Z , Zはそれぞれ原子核 A, Bの電荷であ

AB A B

る。 Hartree- Fock法の代わりに、多配置 SCF (自己無撞着場； Self- consistent Field) 法、配置間相互作用法等のポスト Hartree- Fock法への拡張は、局在化分子軌道による部分を、いわゆる" Frozen Core Orbital"と考えれば、同様の取り扱いにより可能である。（2)式の最初の 3項は、 QM領域に属する原子核と電子のみに関係し、次の 3 項は、 QM領域と表層 QM領域との相互作用に関する項であり、（2)式の最後の 3項は、表層 QM領域に属する原子核と電子のみに関係していることが分かる。そこで、表層 QM領域のみに関わる項、すなわち (2)式の最後の 3項に着目し、これらの項に関してパラメータ aを用いて、（3)式のように分割する。ただし 0≤ a≤ 1である。

[0030] [数 5] E(QM ₊ PQM) = 2Hf ₊∑∑(2J_# - ^) +∑¾^ +∑∑(2J,—

p Q Q P 7 7

+ (3) Z_fl

H - a) ∑2 +∑∑(2 - +∑ 十な ∑2H +∑∑(2 —^) +∑

ここでは αと（1— α )を用いて形式的に分割して、るだけであるので、この段階では、全エネルギーの値は、（2)式と全く同等である。いま、（3)式の最後の項すなわち係数 αが乗じられている項について、平衡構造の近傍で、原子核の座標の微小変化 A q に関して級数展開 (テーラー展開)すると、（4)式のように表現される。

[0031] [数 6]

E(QM + PQM) =∑2H +∑∑(2J_ff - ^) +∑¾^ +∑∑(2J,

ί i j A>B A i J

P Q Q P Ύ 7

r— ^ — \ _A£

+

A B P P Z_AZ (4)

+(1—め ∑2^ +∑∑(2^ - ^,) +∑

Λ>Β ^ΑΒ

ふ QE _A Ι E _{A A}

+a

A dq_A B dq_Adq_B もし、 QM空間に属する原子が全て安定構造の位置にあるとすると、言い換えれば平衡点にあるとすると、展開項数を無限に取れば、（3)式と (4)式は、再び等しい全ェネルギーを与えることになる。

[0032] ここで、 QM空間の表層領域に、経験的ポテンシャルによるエネルギー表式を組み込むために、（4)式の最終項を経験的ポテンシャルの項で置き換えて、（5)式のように、全エネルギー表式を変更する。

[0033] [数 7] E(QM + PQM) =∑2H +∑∑(2/, - ^) + ¾^ +∑∑(2J,

i i

MM空間に属する原子との相互作用は，通常の手続きにしたがって、表層 QM領域の原子との間でのみ、結合距離!:、結合角 Θ、 2面角 Ψを用いて記述すれば、全く過不足なぐ滑らかに MM空間のポテンシャルと連続的につながる。（5)式の最後の項は、経験的ポテンシャルの項であるから、パラメータ αは、表層 QM領域において QMによる項と MMによる項とをどのように混ぜる力、言い換えれば、表層 QM領域にお、て経験的ポテンシャルをどの程度"染み込ませる力，を示して、る。ノラメータ aの好適な値は、分子系や、分子シミュレーションの目的によって変化し得るものであるが、例えば、 0. 2程度に設定される。

[0034] このようにして、 QM空間と MM空間を合わせた全エネルギー E(Total)は、（6)式で表される。

[0035] [数 8]

E(Total) = E{QM + PQM) + E(MM) + E(QM + PQM： MM) (6) もちろん、 E(MM)は、 MM空間に対して分子動力学的方法などにより経験的ポテンシャルを適用して得られたエネルギーである。ここで、 QM空間原子と MM空間原子の相互作用を表現する式は、注目している物理的化学的性質により異なるが、ここでは標準的な、 QM空間原子に対する MM空間原子のつくる静電的効果と、 QM空間の電子分布の変化により MM空間原子に誘起される双極子能率で記述される効果とを示してある。

[0036] 以上の説明においては、表層 QM領域が 1つのみ設定されているものとしているが、 1つの分子にぉ、て複数の表層 QM領域が設定されて、る場合にぉ、ても上述の議論を拡張することができる。その場合には、表層 QM領域の相互間でのクーロン力の項及びファンデルワールス力の項も考慮しなければならない。そこで、（2)式を、複数の表層 QM領域が設定されている場合に拡張する。その場合、（2)式の最初の 3項は、 QM領域に属する原子核と電子のみに関する項であるから、変更を要しない。 (2) 式における次の 3項は、 QM領域と表層 QM領域との相互作用に関する項であるから、各表層 QM領域ごとに独立して求めてそれらを加算すればよい。（2)式の最後の 3 項は、上述と同様の手順にしたがって、（7)式のように変形することができる。

[0037] [数 9]

∑2 +∑∑₍₂ - +

=∑ 2H +∑(2J,

Ζ_ΔΖ

∑∑ 2H ₊∑(2J,., - ,.

2 R

i

=∑(l - aJ ∑ 2 十∑(2 - ) Ζ_ΑΖ_Π

2 R

i ζ_λζ_α

∑ 2H +∑(2 - ) +∑∑ 2 R_t

ここでは、表層 QM領域が PQM個存在するものとし、各表層 QM領域ごとに、その表層 QM領域において QMによる項と MMによる項とをどのように混ぜるかを示すパラメータ α が定められている。（7)式中の E ^mは、（5)式中の E^MM(0)に対応するものであつて、 m番目の表層 QM領域についての平衡状態での全 MMエネルギーを表している。 E^mは、（8)式のように表される。

[0038] [数 10]

(8)式は、クーロン相互作用とファンデルワールス相互作用を含む経験的ポテンシャルによる式であり、最初の項は結合距離 rによる寄与を示し、次の項は結合角 Θによる寄与を示し、 3番目の項は二面角 Ψによる寄与を示し、 4番目の項は表層 QM領域の相互間でのファンデルワールス力の項であり、最後の項は表層 QM領域の相互間でのクーロン力の項である。このように、クーロン力及びファンデルワールス力による項を考慮しなければならないことを除けば、表層 QM領域が複数ある場合も、表層 Q M領域が 1つしかな、場合と同様の取り扱、をすることができる。

[0039] 次に、表層 QM領域の末端原子にっ、て述べる。末端の原子につ、ては、切断したままでは不対電子の問題が発生するので、共有結合している相手側の、現存する原子を水素類似原子と見なしたリンクアトム法を採用する。この時、基底関数としては、原子価電子の原子軌道のみを用い、全電子計算によるマリケン (Mulliken)電荷を再現するように、 Orbital Exponentを事前に調整しておく必要がある。

[0040] 次に、前述の理論に基づく本実施形態の分子シミュレーション方法及び装置について、具体的に説明する。図 3は、本発明の実施の一形態の分子シミュレーション装置の構成を示すブロック図である。

[0041] この分子シミュレーション装置は、シミュレーション対象となる分子の構造データとして、その分子を構成する各原子の座標データなどを使用する。分子シミュレーション装置は、対象分子の構造データなどを格納する初期データ格納部 11と、初期データ格納部 11に格納されて、る構造データを QM空間と MM空間とに分割し、さらに Q M空間については、 MM空間との接続部である表層 QM領域とそうでない QM領域とに分割する領域分割部 12と、 MM空間について、分子動力学法などの MM法に基づいて分子シミュレーション演算を実行する MM演算部 13と、 QM空間について非経験的分子軌道法による分子シミュレーション演算を実行する QM演算部 14と、上述したパラメータ α (0≤ α≤1)が入力するパラメータ入力部 15と、 MM演算部 13 及び QM演算部 14での演算結果を総合して出力する出力部 16と、を備えている。ここで QM演算部 14は、上述した (5)式に示すように、 QM法による全エネルギー E(Q M + PQM)を算出する。すなわち QM演算部 14は、表層 QM領域については、局在化分子軌道を用い、パラメータ αで指定される割合で、 QM法によるエネルギー成分と経験的ポテンシャルによる成分と重ね合わせて、計算を実行する。なお QM領域に対しては、 QM演算部 14は、通常の QM法と同様に、正準軌道 (拡がった分子軌道）を用いて計算を実行する。なお、 MM演算部 13は、経験的ポテンシャルの計算を主として行い、 QM演算部 14は、主として 2電子積分の計算を行うものであるから、それぞれ、最適なハードウェア構成が異なっている。本実施形態の分子シミュレーション装置は、計算機クラスタとして構成するのに適したものである力その場合であっても、 MM演算部 13及び QM演算部 14は、それぞれ別のハードウア構成の計算機から構成されるようにすることが好まし、。

[0042] 図 4は、 QM空間（QM領域 +表層 QM領域）と MM空間との分割の一例を示している。 QM領域と表層 QM領域との境界は、例えば、図示左から 2番目のアミノ酸残基における C C結合の位置に設定されている。また、表層 QM領域と MM空間との境界は、図示左力も 4番目のアミノ酸残基の C C結合に設定されている。これらの境界の位置は、分子シミュレーションの対象となる分子やシミュレーションの目的によつて適宜に定められるものであるが、計算の精度を高めるためには、一重結合（ σ結合 )の位置に境界を定めることが好ま、。

[0043] 次に、前述した理論的側面に基づき、本実施形態での実際の計算手順について、図 5のフローチャートを利用して説明する。

[0044] まず、ステップ 101において、初期データ格納部 11は、シミュレーション対象の分子の構造データ (座標データ）を初期データ格納部 11から読み出し、その分子を Q Μ空間と MM空間とに分割し、さらに、 QM空間を QM領域と表層 QM領域とに分割し、 MM空間に接合する表層 QM領域に属する分子もしくは分子の一部を仮想分子として切出す。 MM空間についての構造データは MM演算部 13に送られ、ステップ 102において、 MM演算部 13は、 MM空間について、分子動力学法などの、経験的ポテンシャルに基づく手法によって分子シミュレーションを実行する。 QM空間（Q M領域及び表層 QM領域）の構造データは、 QM演算部 14に送られる。なお、 MM 計算を行う際に、表層 QM領域内の、 MM空間に隣接する原子の座標データが必要であるから、その座標データも MM演算部 13に送られる。同様に、 QM演算部 14には、 MM空間内の原子のうち、表層 QM領域に隣接する原子の座標データも送られる。

[0045] 構造データが送られてくると、ステップ 103において、 QM演算部 14は、 QM演算、具体的には Hartree-Fock計算を開始する。このとき、上述したように分子を QM領域、表層 QM領域及び MM空間に分割したことに伴い、分子もしくは分子の一部において共有結合が切られ、その結果、分子もしくは分子の一部の両端または一端に不対電子が発生する。 QM演算部 14は、上述したリンクアトム法により、水素類似原子の導入によりそのような不対電子を隠蔽し、まず、表層 QM領域に関して、

Hartree- Fock計算を行うようにする。

[0046] ステップ 103の結果、正準軌道が得られるから、ステップ 104において、その正準軌道が局在化分子軌道に変換され、ステップ 105において、 MM空間の反対側に位置する水素類似原子の分子軌道の係数を無視して、分子軌道の再規格化が行なわれる。その後、 QM演算部 14は、ステップ 106において、表層 QM領域のこれらの局在化分子軌道の内、 MM空間寄りに局在化した分子軌道を適宜複数選択して Lowdinの直交化 [6]を行い、規格直交化された局在化分子軌道基底を得る。ここで直交化を行うのは、水素類似原子 (リンクアトム)の導入などの影響により、直交性が保証されて!、な！/、からである。

[0047] QM空間内の電子のうち、ステップ 103— 106での処理対象となったもの以外の電子、すなわち QM領域の電子については、 QM演算部 14は、ステップ 107において、非経験的分子軌道法における通常の手順によって、初期分子軌道の係数を求める。 QM演算部 14は、ステップ 106とステップ 107で得た分子軌道を用いて、反対称化された波動関数を作り、上述した (5)式の全エネルギー E(QM + PQM)を求める。このとき、（5)式に含まれるパラメータが必要となるが、このパラメータは、パラメ一タ入力部 15に外部から入力して、パラメータ入力部 15から QM演算部 14に与えられる。その後、 QM演算部 14は、ステップ 109において、変分法に基づく自己無撞着（ SCF)法により、局在化分子軌道以外の分子軌道の係数を最適化し、 QM空間の全エネノレギーを求める。

[0048] 以上の処理によって、 MM空間の全エネルギーと QM空間の全エネルギーとが求められたことになるから、出力部 16は、ステップ 110において、両方の空間の全エネルギーを (6)式に基づいて組み合わせることにより、シミュレーション対象の分子全体のエネノレギー E(Total)を出力する。

[0049] 分子の安定構造やィ匕学反応トラジェクトリーを求めるのに必要な原子核に働く力は、（6)式を原子核の座標で微分することにより求められる。この時、分子軌道の係数に対する偏微分項が現れるが、正準軌道に対する偏微分項はエネルギー勾配法 [7]により、重なり積分の偏微分に置き換えられる。これに対し局在化分子軌道の係数に対する偏微分項の計算は必要となる力 CPHF(Coupled Perturbed Hartree-Fock)法 [8]により、厳密に計算することができる。次いで、新しい原子核の座標を用いて、表層 QM領域の分子もしくは分子の一部について、局在化分子軌道を新たに求める。この時、前回の局在化分子軌道との重なりを求め、最も重なりの大きい軌道を選択すれば、当初選択された局在化分子軌道に最も近い軌道が自動的に求められる。この局在化分子軌道の計算、 Hartree-Fock法による全エネルギーと分子軌道係数の決定、原子核に働く力による新しい原子核の座標の算出、を繰り返せば、安定構造や化学反応トラジヱクトリーの計算が行える。

[0050] 以上説明した分子シミュレーション装置は、典型的には、計算機クラスタによって実現されるものであり、計算機クラスタには、初期データ格納部 11、領域分割部 12、パラメータ入力部 15及び出力部 16として機能する制御用の計算機と、 MM演算部 13 として機能する計算機または計算機群と、 QM演算部 14として機能する計算機または計算機群が設けられる。これらの計算機は、それぞれその計算機が果たすべき機能を実行するためのプログラムを読み込むことによって、制御用の計算機、 MM演算用の計算機あるいは QM演算用の計算機として機能することになる。そのようなプログラムは、磁気テープ（MT)あるいは CD— ROMなどの記録媒体を介して、あるいはネットワークを介して、計算機に読み込まれる。

[0051] 比較的小規模な分子を対象とした場合には、単一の計算機を用いて、上述した分子シミュレーションを実行することも可能である。その場合には、上述した手順で分子シミュレーションを実行するための計算機プログラムを、スーパーコンピュータやパーソナルコンピュータなどの計算機に読み込ませ、そのプログラムを実行させればよヽ。分子シミュレーションを行うためのプログラムは、磁気テープや CD— ROMなどの記録媒体によって、あるいはネットワークを介して計算機に読み込まれることになる。

[0052] 上述したプログラム、そのようなプログラムを格納した記録媒体、そのようなプロダラムカもなるプログラムプロダクツも、本発明の範疇に含まれる。

産業上の利用可能性

[0053] タンパク質は、生体中で様々な機能を実現している力酵素という観点で捉えると、化学関連産業全般に関わる重要な化学物質である。したがって、本発明に基づいて、タンパク質と基質の相互作用を、量子力学的手法に則って、高い信頼性でシミュレーシヨンできることは、医薬品や機能性食品の製造、化学工業分野全般、さらには、環境保全物質の開発など、広範囲の産業に有効な研究手法を提供することになる。

Claims

請求の範囲

[1] シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、前記 QM空間に対して非経験的分子軌道法を適用し、前記 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレ一シヨン方法であって、

記憶部から、前記シミュレーション対象の分子または分子の一部を構造データを取り出して前記 QM空間及び前記 MM空間に分割する段階と、

前記 QM空間に関する非経験的分子軌道法における全エネルギー表式の一部を経験的ポテンシャルで置き換える段階と、

を有する、分子シミュレーション方法。

[2] シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、前記 QM空間に対して非経験的分子軌道法を適用し、前記 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレ一シヨン方法であって、

記憶部から、前記シミュレーション対象の分子または分子の一部の構造データを取り出して前記 QM空間及び前記 MM空間に分割し、さらに、前記 QM空間を、前記 MM空間に隣接する領域である表層 QM領域と、前記表層 QM領域以外の領域である QM領域とに分割する段階と、

前記表層 QM領域に関して非経験的分子軌道法による全エネルギー表式を得る段階と、

前記 QM領域に関して非経験的分子軌道法による全エネルギー表式を得る段階と前記表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換える段階と、

前記 QM空間についての非経験的分子軌道法による全エネルギーを求める段階とを有する、分子シミュレーション方法。

パラメータを入力する段階をさらに有し、前記表層 QM領域の全エネルギ一部を経験的ポテンシャルに基づく項に置き換える割合が前記パラメータによって指定される、請求項 2に記載の方法。

[4] 複数の前記表層 QM領域が設定され、前記各表層 QM領域ごとに前記パラメータが入力され、表層 QM領域相互間のクーロン相互作用及びファンデルワールス相互作用を含む経験的ポテンシャルが用いられる、請求項 3に記載の方法。

[5] 前記表層 QM領域における波動関数を局在化分子軌道で表現する段階をさらに有する、請求項 2乃至 4のいずれか 1項に記載の方法。

[6] 分子構造最適化及び時間発展計算を行う段階をさらに有し、前記分子構造最適化及び時間発展計算に際し前記局在化分子軌道を更新する、請求項 5に記載の方法

[7] シミュレーション対象の分子または分子の一部を QM空間と MM空間とに分割し、前記 QM空間に対して非経験的分子軌道法を適用し、前記 MM空間に対しては経験的ポテンシャルに基づく方法を適用して分子シミュレーションを行う分子シミュレ一シヨン装置であって、

前記シミュレーション対象の分子または分子の一部の構造データを格納する記憶部と、

前記記憶部から、前記シミュレーション対象の分子または分子の一部の構造データを取り出して前記 QM空間及び前記 MM空間に分割し、さらに、前記 QM空間を、前記 MM空間に隣接する領域である表層 QM領域と、前記表層 QM領域以外の領域である QM領域とに分割する領域分割部と、

前記表層 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、前記 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、前記表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換え、前記 QM空間についての非経験的分子軌道法による全エネルギーを求める第 1の演算部と、

を有する、分子シミュレーション装置。

[8] パラメータを入力するパラメータ入力部をさらに有し、前記第 1の演算部が前記表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換える割合が前記パラメータによって指定される、請求項 7に記載の装置。

[9] 複数の前記表層 QM領域が設定され、前記各表層 QM領域ごとに前記パラメータが入力され、表層 QM領域相互間のクーロン相互作用及びファンデルワールス相互作用を含む経験的ポテンシャルが用いられる、請求項 8に記載の装置。

[10] 前記 MM空間に関して経験的ポテンシャルによる方法に基づいて全エネルギーを求める第 2の演算部と、

前記第 1の演算部で求めた全エネルギーと前記第 2の演算部で求めた全エネルギ一とから前記シミュレーション対象の分子または分子の一部の全エネルギーを算出する出力部と、

をさらに有する、請求項 7乃至 9のいずれか 1項に記載の装置。

[11] コンピュータを、

シミュレーション対象の分子または分子の一部の構造データを格納する記憶部、前記記憶部から、前記シミュレーション対象の分子または分子の一部を構造データを取り出して前記 QM空間及び前記 MM空間に分割し、さらに、前記 QM空間を、前記 MM空間に隣接する領域である表層 QM領域と、前記表層 QM領域以外の領域である QM領域とに分割する領域分割部、

前記表層 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、前記 QM領域に関して非経験的分子軌道法による全エネルギー表式を求め、前記表層 QM領域の全エネルギー表式の一部を経験的ポテンシャルに基づく項に置き換え、非経験的分子軌道法によって前記 QM空間についての全エネルギーを求める第 1の演算部、

として機能させるプログラム。

[12] 前記コンピュータを、さらに、

前記 MM空間に関して経験的ポテンシャルによる方法に基づ、て全エネルギーを求める第 2の演算部、

前記第 1の演算部で求めた全エネルギーと前記第 2の演算部で求めた全エネルギ一とから前記シミュレーション対象の分子または分子の一部の全エネルギーを算出する出力部、として機能させる、請求項 11に記載のプログラム。

コンピュータが読み取り可能な記録媒体であって、請求項 11または 12に記載のプログラムを格納した記録媒体。