# [gmx-users] different result for entropy with normal mode analysis and schlitter-approximation

oliver.kuhn at uni-duisburg-essen.de oliver.kuhn at uni-duisburg-essen.de
Tue Apr 28 14:29:23 CEST 2009

```Dear Gromacs Users,
I'm trying to calculate entropies from a md trajectory using g_anaeig.
There are two ways to go (question at bottom ;-):

1. NMA and quasi-harmonic approximation: Use a bunch of snapshots (maybe
5-20), minimize each of them to very low maximum forces, calculate the
hessian matrix, diagonalize and use g_anaeig to calculate the entropy from
the resulting eigenvector-matrix assuring that there are no negative
eigenvalues in the eigenvectors 7 to N (first six eigenvectors will not be
part of the calculation). - as follows:

# Energy Minimization
grompp_d -f em_nma.mdp -t md.fitted.trr -time \$t -c md.gro -p protein.top
-o \$t.em.tpr
mdrun_d -v -deffnm \$t.em -table table6-12_4r_doublePrecision.xvg -tablep
table6-12_4r_doublePrecision.xvg

# Hessian Matrix
grompp_d -f nma.mdp -t \$t.em.trr -c md.gro -p protein.top -o \$t.hessian.tpr
mdrun_d -v -deffnm \$t.hessian -table table6-12_4r_doublePrecision.xvg
-tablep table6-12_4r_doublePrecision.xvg

# Diagonalization of the Hessian
g_nmeig_d -f \$t.hessian.mtx -s \$t.hessian.tpr -first 1 -last 10000 -v
\$t.eigenvec.trr

# Entropy calculation - vibrational (without first 6 modes)
g_anaeig_d -v \$t.eigenvec.trr -f \$t.em.trr -s \$t.hessian.tpr -temp 298.15
-nevskip 6 -entropy 2>&1 | tee out.anaeig.Svib.\$t

grep 'The Entropy due to the Quasi Harmonic approximation is'
out.anaeig.Svib.\$t | awk '{print \$10}' >> result/Svib.nma

I use distance-dependent dielectric e=4r, but that doesn't make much
difference.

2. Schlitter approximation based on covariance: Use all snapshots of the
md trajectory, calculate the covariance matrix (g_covar), - diagonalized
matrix will be returned -, and subsequently calculate the entropy with
g_anaeig. - as follows:

# covariance matrix as time average over configurations
g_covar_d -f md\$i.trr -s md.gro -v md\$i.eigenvec.trr -mwa -av
average.\$i.pdb -ascii covar.\$i -xpm covar.\$i -xpma covara.\$i -l covar.\$i
-o md\$i.eigenval.xvg <<- EOF
0
0
EOF

# Analysis of the principal components (and entropy calculation)
g_anaeig -v md\$i.eigenvec.trr -f md\$i.trr -s md.gro -first 1 -last -1
-entropy > out.anaeig.schlitter.\$i

grep 'The Entropy due to the Schlitter formula is' out.anaeig.schlitter.\$i
| awk '{print \$9}' >> result/Svib.schlitter

Somebody before mentioned, he would like to have the undiagonalized
covariance matrix as input for the entropy calculation, I think, that
doesn't make a difference, am I right?

So, practically, I tried to reproduce entropy from Schlitter 1993. A
simulation of a deca-alanine-helix in vacuo in the old gmx force-field
with vdw-cut-off etc. and I could reproduce the value of ca. 700
kJoule/mol K with the Schlitter approximation.
And now the question, why don't I get the same range of values when doing
normal-mode analysis (as described above)?

values of the Schlitter-approximation (for different simulation lengths):
667.365
685.594
681.259
680.269
values given by the quasi-harmonic approximation when calculating from
covariance:
582.731
596.97
590.71
589.07

values from NMA and quasi-harmonic approximation (for 3 snapshots):
21662.9
21674.9
21662.9

So, there is a factor of round 30 between hessian- and covariance-based
entropy!?

I'm totally stuck with this.
If anybody has experience with this phenomenon, any help is appreciated.
Thanx in advance

Oliver Kuhn

```

More information about the gromacs.org_gmx-users mailing list