Microsoft PowerPoint pptx

Size: px

Start display at page:

Download "Microsoft PowerPoint pptx"

ゆあももき
5 years ago
Views:

1 平成 24 年 5 月 21 日分子モデリングと分子シミュレーションポテンシャルエネルギー東京大学大学院農学生命科学研究科アグリバイオインフォマティクス教育研究プログラム寺田透

2 講義予定 1. 5 月 21 日 ( 月 ) ポテンシャルエネルギー 2. 5 月 28 日 ( 月 ) 分子動力学法とモンテカルロ法 3. 6 月 04 日 ( 月 ) 分子動力学法の応用 4. 6 月 11 日 ( 月 ) 複合体構造モデリング参考図書 : 岡崎進コンピューターシミュレーションの基礎化学同人

3 本日の講義内容分子軌道法実習課題 1 分子力学法エネルギー最小化分子力学法実習課題 2

4 立体構造とエネルギー物体に力をかけて変形させると物体の持つポテンシャルエネルギーが大きくなる力同様に分子も変形するとその分子が持つポテンシャルエネルギーが変化する

5 分子のポテンシャルエネルギー分子のポテンシャルエネルギーは Schrödinger 方程式を分子軌道法を用いて近似的に解くことで計算できる Hˆ Hˆ E elec elec E elec N i1 E 1 2 elec M A1 M 2 i BA N i1 A1 Z AZ R AB B M Z r 分子のポテンシャルエネルギー A ia N i1 N ji 1 r N: 電子数 M: 原子数 Z A : 原子 A の原子番号 Φ: 波動関数 ij

6 分子軌道法実習 (1) 本実習では量子化学計算ソフトウェア Gaussian 09Wを用いるデスクトップにあるアイコンをダブルクリックしてこのソフトウェアのグラフィックユーザーインターフェイスGaussView 5.0を起動 Current Fragment Control Panel Molecule View Window

7 分子軌道法実習 (2) 1. Control PanelのRing Fragment をクリックし Current Fragmentが benzeneになっていることを確認して Molecule View Windowの中を左クリック 2. Control Panelのメニューから Calculate Gaussian Calculation Setup を選択 3. Job typeを Energy MethodのBasis setを6-31g(d) に設定し Submit 4. インプットファイルを保存するか聞かれるので Save しデスクトップに benzene.gjf として保存する 5. Run Gaussianウィンドウが出るので OK

8 分子軌道法実習 (3) 6. 計算が終わったらGaussian windowを閉じるか聞かれるのではい 7. Gaussian Job Completedウィンドウでは benzene.logファイルを選択し OK 8. Control Panelのメニューから Results Summary を選択 E(RHF) の欄に分子のポテンシャルエネルギー * が表示されている * 1 a.u. = kcal/mol

9 用語の解説 Method:Schrödinger 方程式を解くのに用いる近似法 Hartree-Fock:ab initio 法における基本 Semi-empirical: 半経験的 ( 大規模な系を小さな計算コストで扱えるが結果の信頼性は低い ) DFT: 密度汎関数法 ( 電子相関の効果を比較的小さな計算コストで取り込むことができる ) Basis Set: 各原子におく基底関数系 STO-3G 3-21G 6-31G 6-311G の順により複雑な軌道を表現できる必要に応じて diffuse 関数 (+ ++) 分極関数 ((d) (d,p) など ) を加える

10 分子を変形してみよう Control PanelのModify Bond をクリックした後ベンゼンの隣り合う2つの原子をクリック次のように結合長を変える同様にエネルギーを計算してみよう ( ファイルはデスクトップに benzene2.gjf として保存 )

11 エネルギー最小化分子を変形するとエネルギーが大きくなるエネルギーが小さくなるように変形させることでもとの構造に戻してみよう Job Typeを Optimization に設定し Submit ( ファイルはデスクトップにbenzene3.gjf として保存 )

12 2 原子分子のエネルギー (1) 1. 講義のページから H2.gjf をダウンロードしデスクトップに保存する 2. スタートメニューからすべてのプログラム Gaussian 09W Gaussian 09W を選択し起動 3. メニューの File Open で H2.gjf を開く 4. Existing Job Edit ウィンドウが現れるのでこのメニューから File Exit & Run 5. Output File 名を聞かれるのでデスクトップに H2.out として保存する

13 2 原子分子のエネルギー (2) 6. 計算が終了したらメニューから File Exit 7. GaussView 5.0 を起動し Control Panel のメニューの File Open で H2.out を開く ( ファイルの種類を Gaussian Output Files (*.out *.log) にする ) 8. Control Panel のメニューの Results Scan を開く 9. Scan plot ウィンドウの内部を右クリックし Save Data を選択しデスクトップに H2_scan.txt として保存

14 2 原子分子のエネルギー (3) Energy [a.u.] H 2 の共有結合長の変化に伴うエネルギーの変化 y = x x x x Bond length [Å] H2_scan.txtをExcelで開きグラフを書くと 4 次関数でよく近似できることがわかる構造変化に伴うエネルギー変化をあらかじめモデル化しておくことで低い計算コストでエネルギーを求めることができる

15 Energy [a.u.] E E エネルギー関数の近似 H 2 の共有結合長の変化に伴うエネルギーの変化 y = x x x x Bond length [Å] r r Er 0 r r Er kr r 2 0 r E r 実際にはエネルギー最小状態のまわりを熱ゆらぎしている生体分子のシミュレーションで扱う温度は 300 K 程度 kt = 0.6 kcal mol 1 = 10 3 a.u. 実際の結合長の変化は 0.04 Å 程度 2 r 1 Er 2 3 r r O r r r 2 r 0 r が小さければ 2 次式で近似できる

16 参考 : テイラー展開関数 f(x) の x = p のまわりでの展開 x = p からの微小な変位を x とおく f p x f p f p 1 k! x df dx k d f dx df dx x k x p x x 1 2! x 2 2 x 1 d f x 2 3 x x O x x p x p k 2! 2 d f dx dx 2 次までの展開 2 x p x p x 2 誤差項 (x 3 のオーダー )

運動の比較 E F r kr r0 E r r r 2k r 2k, r r r t Acos t r0 0 2 2k 1 r m1m 2 m m 2 r 0 黒実線

17 運動の比較 E F r kr r0 E r r r 2k r 2k, r r r t Acos t r k 1 r m1m 2 m m 2 r 0 黒実線 :ab inito MD 赤破線 : 古典近似解析解ポテンシャルエネルギーを量子化学計算を用いて求める代わりにエネルギー関数を用いて近似してもほぼ同じ軌道が得られる

3 原子分子のエネルギー energy [kcal mol 1 ] 5 4 3 2 1 0 energy [kcal mol 1 ] 0.9 0.91 0.

5 r [Å] 1 95 101 0 90 100 110 120 θ [degree] 113 107 4-5 3-4 2-3 1-2 0-1 r

18 3 原子分子のエネルギー energy [kcal mol 1 ] energy [kcal mol 1 ] r [Å] θ [degree] r [Å] モデル系 : 水分子 r θに関するエネルギー関数の形は rに依存しないエネルギー関数は以下のように近似できる E θ 2 r, k r r k 2 b 0 a 0

分子間相互作用の計算 (1) 1. GaussView 5.0 を起動し Control Panel の Element Fragment をクリックして O を選択 Current Fragment が Oxygen Tetravalent になっていることを確認して Molecule View Window の中を左クリック 2. ややはなれた別の位置を左クリック 3.

19 分子間相互作用の計算 (1) 1. GaussView 5.0 を起動し Control Panel の Element Fragment をクリックして O を選択 Current Fragment が Oxygen Tetravalent になっていることを確認して Molecule View Window の中を左クリック 2. ややはなれた別の位置を左クリック 3. Control Panel の Element Fragment をクリックして C を選択 Current Fragment が Carbon Tetrahedral になっていることを確認して Molecule View Window の H 2 O の水素原子から 1 つ選びクリック CH 3 基に置換 4. 同様にもう 1 つの H 2 O 分子の水素原子を CH 3 基に置換

20 分子間相互作用の計算 (2) 5. Control Panel メニューの Calculate Gaussian Calculation Setup を開き Job type を Energy Method の Basis set を 6-31G(d) に設定し Submit ( ファイルはデスクトップに methanol2.gjf として保存 ) 6. 同様に CH 3 OH 1 分子についてもエネルギーを計算 ( ファイルはデスクトップに methanol1.gjf として保存 ) 7. 以下の式を用いて相互作用エネルギーを求める E E AB E A E B E = [2 ( )] = a.u. = 1.15 kcal mol 1

21 分子間相互作用の成分 van der Waals 相互作用斥力 +q q +q 引力 q 静電相互作用

22 静電相互作用エネルギー関数量子化学計算では電子密度 (r) が計算される電子密度 (r) から分子のまわりの静電ポテンシャル(r) が計算できる静電ポテンシャルを再現するように原子に部分電荷をおく静電相互作用を部分電荷どうしの相互作用として計算する E ij q q i i j R R j 電子密度 q +q 静電ポテンシャル +

23 van der Waals 相互作用 5 メタン 4 energy [kcal mol 1 ] HF MP2 メタン分子間に引力が働いている r -1 distance [Å] 無極性分子間に働く引力 ( 分散力 ) は電子相関を考慮した高精度な量子化学計算によって初めて現れる

24 van der Waals 引力の起源無極性分子であっても電子雲のゆらぎによって瞬間的に双極子が生じ近くの分子に双極子を誘起する双極子誘起双極子相互作用は常に引力引力項の起源 <>=0 <E(r)>=0

25 van der Waals エネルギー関数 energy [kcal mol 1 ] distance [Å] 12 6 EvdW r 4 r r 斥力項引力項 r = σ で E vdw = 0 r = 2 1/6 σ で E vdw = ε ( 最小値 ) σ = 3.34 Å ε = 0.36 kcal mol 1

二面角のエネルギーモデル系 : ブタンブタンの二面角の変化に伴うエネルギーの変化 16 Energy [kcal mol 1 ] 14 12 10 8 6 4

26 二面角のエネルギーモデル系 : ブタンブタンの二面角の変化に伴うエネルギーの変化 16 Energy [kcal mol 1 ] Gaussian job file: butane.gjf Dihedral angle [degree]

二面角のエネルギー関数周期関数の重ね合わせと1 番目と4 番目の原子間の van der Waals 相互作用の和で表される E 0 k 1 cos n i i 16 14 i i 4 r14 12 r14 6 Energy [kcal mol 1

27 二面角のエネルギー関数周期関数の重ね合わせと1 番目と4 番目の原子間の van der Waals 相互作用の和で表される E 0 k 1 cos n i i i i 4 r14 12 r14 6 Energy [kcal mol 1 ] van der Waals 0 周期関数 Dihedral angle [degree] 1 番目と 4 番目の原子間のみ van der Waals 相互作用を考慮

28 課題 1 ブタンの二面角の変化に伴うエネルギーの変化を表すポテンシャルエネルギー関数のパラメータ n δ ε σ を決定せよ講義のページにあるbutane_scan.xlsxをダウンロードして利用すること 1 番目と4 番目の原子間距離はGauss View 5.0 のScan Plotウィンドウのメニューの Plot Plot Molecular Property を開き Bond 1 4 を指定して求めること

29 分子力学法量子化学計算は計算コストが大きい構造変化に伴うエネルギー変化をあらかじめ関数でモデル化しておくことで低い計算コストでエネルギーを求める E 1 結合長 k r r k k n r 1 cos b ij b b b ij 4 ij rij 12 a 2 結合角 a ij rij a 6 4 非共有結合相互作用 a q q i r ij j d d 3 二面角 d d d

30 力場パラメータの決定力場パラメータとは? ポテンシャル関数で用いられるパラメータ ( 平衡結合長ばね定数部分電荷など ) 非経験的パラメータ量子化学計算の結果からパラメータを求める経験的パラメータ構造や熱力学量などの実験値を再現するようにパラメータを決める

31 問題点と解決法 (1) 生体高分子は多数の原子からなる全体について量子化学計算を行うのは困難化学的 ( 混成軌道の種類置換基の種類 ) に類似した原子は同じ原子種とみなし同じ力場パラメータを割り当てる同じ原子種を含む小さなモデル化合物についてパラメータを決定する

32 原子種の例 (CHARMm) HA CT 同じ原子種 HA CT OT HO 同じ原子種 HA HA HA CT CT C6R 同じ原子種 HO OT

33 問題点と解決法 (2) 凝縮相 ( 液相など ) では分子が接近しているため第 3 の分子の位置が 2 つの分子の相互作用に影響を与える気相で決めたポテンシャルエネルギー関数をそのまま適用できない気相凝縮相

E 有効ポテンシャルエネルギー r, r, r Er, r Er, r Er, r Er, r, r 1 2 3 E 1 eff eff eff

分子系の相互作用エネルギー E eff (r 1,r 2 ): 有効 2 体間相互作用エネルギー厳密な多体間相互作用の効果を 2

34 E 有効ポテンシャルエネルギー r, r, r Er, r Er, r Er, r Er, r, r E 1 eff eff eff r, r E r, r E r, r 1 E(r 1,r 2,r 3 ):3 分子系の相互作用エネルギー E(r 1,r 2 ):2 分子系の相互作用エネルギー E eff (r 1,r 2 ): 有効 2 体間相互作用エネルギー厳密な多体間相互作用の効果を 2 体間相互作用のエネルギー関数に取り込む凝縮相このエネルギー関数のパラメータを実験値を再現するように決める

35 水分子のモデル (1) r(oh) HOH r* q H SPC TIP3P r(oh) [Å], HOH [degree] r* [Å 12 kcal mol -1 ], 2r* 6 2 [Å 6 kcal mol -1 ] q O = 2q H r* van der Waals 相互作用は酸素原子間のみ計算する Jorgensen et al. J. Chem. Phys. 79, 926 (1983)

36 水分子のモデル (2) 密度蒸発熱定圧比熱膨張率圧縮率 SPC TIP3P 実験値密度 [g cm 3 ] 蒸発熱 [kcal mol 1 ] 定圧比熱 [cal mol 1 K] 膨張率 [10 5 K 1 ] 圧縮率 [10 6 atm 1 ] いずれも 25 1 atm における値 Jorgensen et al. J. Chem. Phys. 79, 926 (1983)

37 生体高分子の力場パラメータポテンシャルエネルギー関数のパラメータ ( 力場パラメータ ) は分子シミュレーションのソフトウェアと共に配布されている AMBER CHARMM GROMOS, GROMACS

38 エネルギー最小化 (1) 立体構造 ( 座標 ) を変化させてエネルギー関数の値が最小になるようにすること立体構造最適化とも呼ばれる分子動力学シミュレーションを行う際には原子同士のぶつかりを排除したり構造のゆがみを正すために最初に必ず行う

39 エネルギー最小化 (2) 1 次のアルゴリズム Steepest descent( 最急降下 ) 法単純だが収束までに多段階を要することがある Conjugate gradient( 共役勾配 ) 法エネルギー関数が N 次元の 2 次形式で近似できる場合 N 回の操作で極小に到達する 2 次のアルゴリズム Newton-Raphson 法収束は早いが Hessian( 2 E) の計算に膨大な時間がかかる

40 Discovery Studio 3.0 Client の起動 Discovery Studio 3.0 Clientのアイコンをダブルクリックして起動メニューの View Explores Tools を選択し Toolsタブを表示させておくメニューの View Toolbars Sketching を選択し Sketching tool barを表示させておく

41 低分子化合物の生成 (1) 1. メニューの File New Molecule Window を選択し新しいウィンドウを開く 2. SketchingツールバーからRing を左クリック Molecule Windowの真ん中付近を左クリック 3. ViewツールバーのRotate Sketchingモードを解除を左クリック 4. メニューの View Hierarchy を選択し Hierarchy Windowを表示 5. 炭素原子を1つ左クリックし選択黄色でマークされる 6. Hierarchy Windowのツリーを展開し対応する原子がマークされていることを確認

42 低分子化合物の生成 (2) 7. Molecule Window の中で右クリックしてメニューを出し Attributes of C5 を選択 ( 原子名は選択した原子によって異なる ) Hybridization が Sp3 になっていることを確認 8. Molecule Window の中で何もないところを左クリックして選択を解除 9. Ctrl キーと A を同時に押して全原子を選択 10. メニューの Chemistry Bond Aromatic を選択 Attributes の Hybridization が Sp2 になっていることを確認せよ

43 低分子化合物の生成 (3) 11. Simulation ボタンを左クリック 12. Tools タブの Change Forcefield を左クリックして展開する 13. Forcefield を CHARMm Partial Charges を Momany-Rone とし Apply Forcefield を左クリック 14. 水素原子が付加され Forcefield Status が Molecule 1 typed with CHARMm となっていることを確認 15. Attributes で部分電荷や原子種がアサインされていることを確認せよ

44 参考 :Molecule Window の操作 (1) Molecules Windowの中を左クリックしてアクティブにしてから以下の操作を行う回転をクリックしてから 3D Window の中で左ドラッグ並進をクリックしてから 3D Window の中で左ドラッグズームをクリックしてから 3D Window の中で左ドラッグ

45 参考 :Molecule Window の操作 (2) 選択 Molecule Window 上で原子をクリックその原子が選択され黄色でマークされる Molecule Window 上で原子をダブルクリックその原子を含む残基が選択されマークされる Hierarchy Window でもチェイン残基原子グループ (backbone など ) 単位で選択できる何もないところをクリックすると選択を解除できる Hierarchy Window では Ctrl キーを Molecule Window では Shift キーを押しながらクリックすると複数選択ができる属性 (attribute) 選択した後右クリックで表示できる

46 参考 :Molecule Window の操作 (3) Home 最初の向き位置に戻す Fit to Screen ( 選択した ) 構造を Window にフィットするように並進拡大縮小 Center Structure ( 選択した ) 構造の中心が Window の中心に来るように並進

47 参考 : 力場パラメータ Discovery Studio には CHARMM が統合されているシミュレーションの対象に応じて力場パラメータを使い分けるのが良い CHARMm:General-purpose Momany and Rone allatom forcefield that also provides automatic parameter estimation charmm22:academic all-atom forcefield used for simulating protein systems charmm27:academic all-atom forcefield used for simulating DNA and protein systems

48 エネルギー最小化 (1) 1. メニューの Window Close All で現在出ているWindowを閉じる (Saveするか聞かれるが No で良い) 2. メニューの File New Molecule Window を選択し新しいウィンドウを開く 3. Sketchingツールを使ってdibenzo-p-dioxin O を作れ ( 原子の変更はメニューの Chemistry Element ) O

49 エネルギー最小化 (2) 4. Forcefield には CHARMm を用いる 5. Tools タブの Run Simulation のツリーを展開し Minimization をクリック設定ウィンドウが現れる 6. 右下のように設定して Run ( 計算が失敗する場合は Clear Forcefield して再度 Apply Forcefield ) 7. メニューの View Data Table でエネルギー値を確認せよ

50 SMILES による化合物の表現低分子化合物の構造の表現方法の 1 つの SMILES がある新しい Molecule Window を開きメニューの File Insert From SMILES を選択 Smiles string に入力すると対応する分子が生成される (aromatic は小文字 ) benzene: c1ccccc1 dibenzo-p-dioxin: c13ccccc1oc2ccccc2o3 alanine: [N+][C@@H](C)C(=O)[O-]

51 課題 2 以下の分子を作成しエネルギー最小化計算をせよ Oseltamivir: インフルエンザ治療薬 ( 商品名 :Tamiflu) Force field は CHARMm とすることエネルギー最小化後の構造の図を PNG 形式で保存すること

52 課題の提出課題 1については作成したExcelファイル課題 2についてはエネルギー最小化後の構造の図のファイルを添付してメールで寺田宛に送ること課題 2についてはエネルギー最小化前と後のエネルギー値を本文に記載すること件名は分子モデリング課題とし本文に氏名と学生証番号を必ず明記すること

Microsoft PowerPoint pptx

Microsoft PowerPoint pptx 平成 23 年 5 月 25 日分子モデリングと分子シミュレーションポテンシャルエネルギー東京大学大学院農学生命科学研究科アグリバイオインフォマティクス教育研究プログラム寺田透講義予定 1. 5 月 25 日 ( 水 ) ポテンシャルエネルギー 2. 6 月 01 日 ( 水 ) 分子動力学法とモンテカルロ法 3. 6 月 08 日 ( 水 ) 分子動力学法の応用 4. 6 月 15 日 (