求助柔性多肽分子动力学模拟的收敛问题

dzhen1

      如图所示，针对多肽（九个氨基酸残基）运行了3600ns的分子动力学模拟，观察RMSD发现结构始终未能收敛。根据sob老师关于分子动力学模拟平衡的文章中解释，推测该情况是由于多肽在300K下不断越过势垒，构象频繁变化所导致的。      本人的问题为能否截取1200~2400ns内，RMSD相对平衡区域中的一帧作为该多肽的平衡结构，进行下一步的研究呢。（附图分别为1200ns 1800ns 2400ns 时多肽的结构，可以看出形成的阿尔法螺旋又打开了）

sobereva

不要自己在标题里手写[求助]这种标签，http://bbs.keinsci.com/thread-9348-1-1.html里明确说了。这次给你改了，以后注意

sobereva

看具体之后研究什么

有的小肽本身就没有主导性的稳定结构，但力场参数、模拟设置不当也可能导致构象不稳定

dzhen1

谢谢您的回复，老师。我接下来准备用模拟得到的多肽二级结构生成其表面屏蔽电荷分布模型。力场我使用的是Amber99sb-star-ildn, 进行成品模拟之前依次进行了能量最小化、NVT和NPT平衡，模拟参数我直接用的Gromacs溶菌酶教程的mdp文件（参数如下），想请教您一下这样设置合理吗。
title                   = OPLS Lysozyme NPT equilibration
; Run parameters
integrator              = md        ; leap-frog integrator
nsteps                  = 50000000    ; 2 * 50000000 = 100000 ps (100 ns)
dt                      = 0.002     ; 2 fs
; Output control
nstxout                 = 0         ; suppress bulky .trr file by specifying
nstvout                 = 0         ; 0 for output frequency of nstxout,
nstfout                 = 0         ; nstvout, and nstfout
nstenergy               = 5000      ; save energies every 10.0 ps
nstlog                  = 5000      ; update log file every 10.0 ps
nstxout-compressed      = 5000      ; save compressed coordinates every 10.0 ps
compressed-x-grps       = System    ; save the whole system
; Bond parameters
continuation            = yes       ; Restarting after NPT
constraint_algorithm    = lincs     ; holonomic constraints
constraints             = h-bonds   ; bonds involving H are constrained
lincs_iter              = 1         ; accuracy of LINCS
lincs_order             = 4         ; also related to accuracy
; Neighborsearching
cutoff-scheme           = Verlet    ; Buffered neighbor searching
ns_type                 = grid      ; search neighboring grid cells
nstlist                 = 10        ; 20 fs, largely irrelevant with Verlet scheme
rcoulomb                = 1.0       ; short-range electrostatic cutoff (in nm)
rvdw                    = 1.0       ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype             = PME       ; Particle Mesh Ewald for long-range electrostatics
pme_order               = 4         ; cubic interpolation
fourierspacing          = 0.16      ; grid spacing for FFT
; Temperature coupling is on
tcoupl                  = V-rescale             ; modified Berendsen thermostat
tc-grps                 = Protein Non-Protein   ; two coupling groups - more accurate
tau_t                   = 0.1     0.1           ; time constant, in ps
ref_t                   = 300     300           ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl                  = Parrinello-Rahman     ; Pressure coupling on in NPT
pcoupltype              = isotropic             ; uniform scaling of box vectors
tau_p                   = 2.0                   ; time constant, in ps
ref_p                   = 1.0                   ; reference pressure, in bar
compressibility         = 4.5e-5                ; isothermal compressibility of water, bar^-1
; Periodic boundary conditions
pbc                     = xyz       ; 3-D PBC
; Dispersion correction
DispCorr                = EnerPres  ; account for cut-off vdW scheme
; Velocity generation
gen_vel                 = no        ; Velocity generation is off

sobereva

dzhen1 发表于 2023-1-5 01:31
谢谢您的回复，老师。我接下来准备用模拟得到的多肽二级结构生成其表面屏蔽电荷分布模型。力场我使用的是Am ...

用AMBER19SB好得多（或者14SB，有现成的力场包）
其它没明显问题

dzhen1

好的，谢谢老师。之前用Amber99SB模拟过，结果确实要好的多。