请选择 进入手机版 | 继续访问电脑版

计算化学公社

 找回密码
 现在注册!
查看: 843|回复: 4

[GROMACS] g_mmpbsa能量分解过程中的异常值

[复制链接]

4

帖子

0

威望

23

eV
积分
27

Level 2 能力者

发表于 2018-3-17 04:16:00 | 显示全部楼层 |阅读模式
本帖最后由 510426762 于 2018-3-17 06:47 编辑

模拟了蛋白质和受体在0.15M NaCl水溶液环境中的结合作用,时间长度100ns。rmsd结果稳定了。
g_mmpbsa分解能量,用的一步法分解,参数如下。
1.jpg
  1. ;Polar calculation: "yes" or "no"
  2. polar                = yes

  3. ;=============
  4. ;PSIZE options
  5. ;=============
  6. ;Factor by which to expand molecular dimensions to get coarsegrid dimensions.
  7. cfac                 = 1.5

  8. ;The desired fine mesh spacing (in A)
  9. gridspace         = 0.5

  10. :Amount (in A) to add to molecular dimensions to get fine grid dimensions.
  11. fadd                 = 5

  12. ;Maximum memory (in MB) available per-processor for a calculation.
  13. gmemceil         = 4000

  14. ;=============================================
  15. ;APBS kwywords for polar solvation calculation
  16. ;=============================================
  17. ;Charge of positive ions
  18. pcharge         = 1

  19. ;Radius of positive charged ions
  20. prad                = 0.95

  21. ;Concentration of positive charged ions
  22. pconc           = 0.150

  23. ;Charge of negative ions
  24. ncharge         = -1

  25. ;Radius of negative charged ions
  26. nrad                = 1.81

  27. ;Concentration of negative charged ions
  28. nconc                 = 0.150

  29. ;Solute dielectric constant
  30. pdie                 = 2

  31. ;Solvent dielectric constant
  32. sdie                 = 80

  33. ;Reference or vacuum dielectric constant
  34. vdie                 = 1

  35. ;Solvent probe radius
  36. srad                 = 1.4

  37. ;Method used to map biomolecular charges on grid. chgm = spl0 or spl2 or spl4
  38. chgm            = spl4

  39. ;Model used to construct dielectric and ionic boundary. srfm = smol or spl2 or spl4
  40. srfm            = smol

  41. ;Value for cubic spline window. Only used in case of srfm = spl2 or spl4.
  42. swin                 = 0.30

  43. ;Numebr of grid point per A^2. Not used when (srad = 0.0) or (srfm = spl2 or spl4)
  44. sdens                 = 10

  45. ;Temperature in K
  46. temp                 = 300

  47. ;Type of boundary condition to solve PB equation. bcfl = zero or sdh or mdh or focus or map
  48. bcfl                 = mdh

  49. ;Non-linear (npbe) or linear (lpbe) PB equation to solve
  50. PBsolver         = lpbe


  51. ;========================================================
  52. ;APBS kwywords for Apolar/Non-polar solvation calculation
  53. ;========================================================
  54. ;Non-polar solvation calculation: "yes" or "no"
  55. apolar                = yes

  56. ;Repulsive contribution to Non-polar
  57. ;===SASA model ====

  58. ;Gamma (Surface Tension) kJ/(mol A^2)
  59. gamma           = 0.0226778

  60. ;Probe radius for SASA (A)
  61. sasrad          = 1.4

  62. ;Offset (c) kJ/mol
  63. sasaconst       = 3.84928

  64. ;===SAV model===
  65. ;Pressure kJ/(mol A^3)
  66. press           = 0

  67. ;Probe radius for SAV (A)
  68. savrad          = 0

  69. ;Offset (c) kJ/mol
  70. savconst        = 0

  71. ;Attractive contribution to Non-polar
  72. ;===WCA model ====
  73. ;using WCA method: "yes" or "no"
  74. WCA             = no

  75. ;Probe radius for WCA
  76. wcarad          = 1.20

  77. ;bulk solvent density in A^3
  78. bconc                = 0.033428

  79. ;displacment in A for surface area derivative calculation
  80. dpos                = 0.05

  81. ;Quadrature grid points per A for molecular surface or solvent accessible surface
  82. APsdens                = 20

  83. ;Quadrature grid spacing in A for volume integral calculations
  84. grid            = 0.45 0.45 0.45

  85. ;Parameter to construct solvent related surface or volume
  86. APsrfm          = sacc

  87. ;Cubic spline window in A for spline based surface definitions
  88. APswin          = 0.3

  89. ;Temperature in K
  90. APtemp          = 300
复制代码
计算出来的能量,发现energy_MM里有异常值。主要集中在s4 s5 s6,Protein-MOL VdW Energy接近于0,VdW = s4 - s0 - s2
  1. @ s0 legend "Protein VdW Energy"
  2. @ s1 legend "Protein Elec. Energy"
  3. @ s2 legend "MOL VdW Energy"
  4. @ s3 legend "MOL Elec. Energy"
  5. @ s4 legend "Protein-MOL VdW Energy"
  6. @ s5 legend "Protein-MOL Elec. Energy"
  7. @ s6 legend "Protein-MOL Total Energy"
  8.           0.000          0.000          0.000          0.000          0.000         -0.001         -0.342         -0.343
  9.         100.000          0.000          0.000          0.000          0.000         -0.001         -0.437         -0.438
  10.         200.000          0.000          0.000          0.000          0.000       -278.876        -37.203       -316.079
  11.         300.000          0.000          0.000          0.000          0.000       -286.627        -43.203       -329.830
  12.         400.000          0.000          0.000          0.000          0.000         -0.001         -0.025         -0.025
  13.         500.000          0.000          0.000          0.000          0.000       -307.671          3.073       -304.598
  14.         600.000          0.000          0.000          0.000          0.000       -285.399        -26.875       -312.274
  15.         700.000          0.000          0.000          0.000          0.000         -0.001          0.014          0.013
  16.         800.000          0.000          0.000          0.000          0.000         -0.001         -0.132         -0.133
  17.         900.000          0.000          0.000          0.000          0.000         -0.001          0.037          0.036
  18.        1000.000          0.000          0.000          0.000          0.000         -0.001          0.319          0.318
  19.        1100.000          0.000          0.000          0.000          0.000         -0.001          0.156          0.155
  20.        1200.000          0.000          0.000          0.000          0.000         -0.001          0.322          0.321
  21.        1300.000          0.000          0.000          0.000          0.000       -260.126         -5.969       -266.095
  22.        1400.000          0.000          0.000          0.000          0.000         -0.001         -0.495         -0.496
  23.        1500.000          0.000          0.000          0.000          0.000         -0.001         -0.023         -0.023
  24.        1600.000          0.000          0.000          0.000          0.000         -0.001         -0.141         -0.142
  25.        1700.000          0.000          0.000          0.000          0.000         -0.001         -0.035         -0.035
  26.        1800.000          0.000          0.000          0.000          0.000         -0.001         -0.068         -0.069
  27.        1900.000          0.000          0.000          0.000          0.000         -0.001          0.088          0.088
  28.        2000.000          0.000          0.000          0.000          0.000       -264.009         -3.774       -267.782
  29.        2100.000          0.000          0.000          0.000          0.000         -0.001         -0.227         -0.227
  30.        2200.000          0.000          0.000          0.000          0.000         -0.001          0.137          0.136
  31.        2300.000          0.000          0.000          0.000          0.000         -0.001         -0.138         -0.139
  32.        2400.000          0.000          0.000          0.000          0.000         -0.001          0.105          0.105
  33.        2500.000          0.000          0.000          0.000          0.000         -0.001          0.039          0.038
复制代码
4.jpg

这个算作异常值么?还是由于这个时刻蛋白质和配体的确没有VdW力作用了呢?


大概情况如上

4

帖子

0

威望

23

eV
积分
27

Level 2 能力者

 楼主| 发表于 2018-3-17 05:10:05 | 显示全部楼层
好像可能解决了。
看到这篇帖子提示了。http://bbs.keinsci.com/thread-6822-1-1.html
g_mmpbsa里也说明了
WARNING: Trajectory should be PBC corrected and molecule should not be PBC broken. To make molecule whole in trajectory, please follow these links: PBC and trjconv.

63

帖子

0

威望

409

eV
积分
472

Level 3 能力者

发表于 2018-3-26 14:55:13 | 显示全部楼层
本帖最后由 HarrisLin 于 2018-3-26 14:58 编辑

是PBC的問題沒錯 要用noPBC分析

4

帖子

0

威望

23

eV
积分
27

Level 2 能力者

 楼主| 发表于 2018-5-30 11:56:19 | 显示全部楼层
510426762 发表于 2018-3-17 05:10
好像可能解决了。
看到这篇帖子提示了。http://bbs.keinsci.com/thread-6822-1-1.html
g_mmpbsa里也说明 ...

是的,之后仔细阅读,发现漏掉了这句。
用这种方式做了noPBC
echo "0"      | gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_whole.xtc -pbc whole -n index.ndx
echo "1 0"   | gmx trjconv -s md_0_1.tpr -f md_0_1_whole.xtc -o md_0_1_mol_center.xtc -center -pbc mol -n index.ndx

60

帖子

0

威望

338

eV
积分
398

Level 3 能力者

发表于 2018-9-14 18:28:23 | 显示全部楼层
学习了,并收藏之。
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

手机版|北京科音自然科学研究中心|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949-1号 )

GMT+8, 2018-11-15 08:56 , Processed in 0.135207 second(s), 26 queries .

快速回复 返回顶部 返回列表