Dear all, I calculated several 1 microsecond MD simulations for an homodimer in a water box. The structure of the Protein is stable along the trajectories excepted its N-terminal region (that we previously demonstrated by NMR relaxation).
I now try to calculate the configurational entropy of each monomer from these simulations. For this, i first created two trajectories to remove rotation and translation and water molecules) for each monomer (one such trajectory for monomer A and one for monomer B). Then i calculate the configurationnal entropy of each monomer from the covariance matrix (using the coor cova entropy command, andricioaei & Karplus J. Chem Phys. 2001, 115, 6289).
The exact command i used (for monomer A) is:
coor cova firstu 24 nunit 1 begin @startf skip @skips stop @stope - sele segid @sega end sele segid @sega end - entropy unit -1 temp 293 resi
I wrote a script to calculate the configurationnal entropy for increased length of simulation, but apparently the total entropy does not converged (results below).
I wonder if someone tried something like this and what can explain this observation. Alternatively a suggestion to calculate configurationnal entropy from MD trajectories would be appreciated. Best regards. Philippe.
Configuration entropy converges relatively slowly. Are the entropy values on the last line of your table calculated for 1microsecond, or just for the 20 ns from 980 to 1000ns? There is something strange about your data. Each entropy column seems to contain two distinct data sets, one which converges to 9.5 and one to 5.3 (for SA).
I use the QUASiharmonic command within the VIBRAN molecule, roughly: read rtf card name top_all36_prot.rtf (open and close are not necessary) read para card name par_all36_prot.rtf read psf card name my_prot.psf open unit @fu read file name my.trj set T 300 coor dyna firstu @fu nunit @nu sele @solute end
calc nmodes = 3 * ?NSEL bomlev -2 vibran nmode @nmodes open unit 11 write form name @s/quasi-@j.txt outu 11 echu 11
ECHO protein quasi firstu @fu nunit @nu begin @begin stop @stop - sele @prot end - thermo resi temp @T
outu 6
end
Lennart Nilsson Karolinska Institutet Stockholm, Sweden
Thank you for your answer. Yes it is exactly why i am surprised. It seems that, for each monomer, there are two S values depending on the length of the simulation and when plotting these values as a function of the length of the simulation, each of these two values seems to converge to a plateau differing by about a factor 2 (see figure attached). The values of S are calculated for the total range meaning for instance for a stopmax of 1000ns, the set of frames are extracted between 20ns to 1000ns each 20ps. So i was expecting that S converge to a single value. Since my first post, i tried to build a traj for one monomer and to work with it but i found the same 'two values" behaviour. So i will stop with this "coor cova" approach and try to do it with quasi.
