[Q] How to calculate host+guest binding free energy and enthalpy changes from CREGEN output? about crest HOT 3 CLOSED

Hello,
In the calculation of G at temperature T your reference state is the total energy of your system, plus thermostatistical contributions from the respective frequency calculation and solvation free energy terms, so
G(T) = E_gas+δG_solv(T)+G_RRHO(T)
The ensemble free energy is another additive term G_conf(T)=-TS_conf, so in your calculation it should be the opposite sign G_tot(T)=G(T)+G_conf(T). For the enthalpy I don't think it can be calculated like this, because the way it is written above one would end up with G(T) again.

Obtaining good ensemble entropies (and hence free energies) is very difficult. Short runtimes will result in incomplete ensembles and therefore in faulty ensemble free energies. Reproducibility can also be an issue. A part of our research is currently dedicated to these problems and we are hoping to finalize it soon with some precise recommendations for ensemble entropy calculations. This will be the update to version 2.11 of crest.

from crest.

Thanks for this Philipp. Is it legitimate to calculate the guest@host binding enthalpy from the thermochemical output of separate vibrational frequency/hessian runs on the minimum energy ensemble structures obtained from a the above CREST runs? ie H(298K)_guest@host-H(298)_host-H(298)_guest

Your point about short runtimes is well taken. In my systems of interest, using the default CREST -nci setup, luckily the host (cucurbit[8]uril) is quite rigid and by visual inspection of the MTD trajectories I can see the guests leave the CB8 after 2-5ps (depending on the Vbias parameter combination), so I end up with a significant proportion of cavity-bound structures versus surface-bound structures. In replicate runs I have been pleasantly surprised at how well CREST always finds the same minimum energy structures from very different starting poses and although the ensembles can differ a little, the distribution of the dominant low-energy structures is always nearly the same, producing less than 1kcal/mol differences in G_conf(T).

I have also done a few longer runs and it appears that, for these systems, the ensemble entropies converge after 30ps to give free energies that differ by less than 1kcal/mol from the 10ps runs. As even the shorter guest@host runs typically take 24-36 hrs on my old 12-core MacPro, I'm OK with the tradeoff of a small loss of accuracy.

Looking forward to what version 2.11 has to offer.

from crest.

Yes, it is legimate to take G(T) (and related thermostatistical values) from seperate calculations and only use G_conf(T) from the crest output. In fact, this might even be the preferred way to do it since thermochemistry calculations purely at GFN level are not as reliable as, e.g. DFT level. Usually G_conf(T) is much smaller than G(T), so you want to have a good description of the latter. For these calculations and good post-processing of the GFN ensemble at DFT level you could take look at the enso repository. This is a project in which we mainly focus on those kind of things.

from crest.