Physics-constrained Bayesian inference of state functions in classical density-functional theory.

The swift progression and expansion of machine learning (ML) have not gone unnoticed within the realm of statistical mechanics. In particular, ML techniques have attracted attention by the classical density-functional theory (DFT) community, as they enable automatic discovery of free-energy functionals to determine the equilibrium-density profile of a many-particle system. Within classical DFT, the external potential accounts for the interaction of the many-particle system with an external field, thus, affecting the density distribution. In this context, we introduce a statistical-learning framework to infer the external potential exerted on a classical many-particle system. We combine a Bayesian inference approach with the classical DFT apparatus to reconstruct the external potential, yielding a probabilistic description of the external-potential functional form with inherent uncertainty quantification. Our framework is exemplified with a grand-canonical one-dimensional classical particle ensemble with excluded volume interactions in a confined geometry. The required training dataset is generated using a Monte Carlo (MC) simulation where the external potential is applied to the grand-canonical ensemble. The resulting particle coordinates from the MC simulation are fed into the learning framework to uncover the external potential. This eventually allows us to characterize the equilibrium density profile of the system by using the tools of DFT. Our approach benchmarks the inferred density against the exact one calculated through the DFT formulation with the true external potential. The proposed Bayesian procedure accurately infers the external potential and the density profile. We also highlight the external-potential uncertainty quantification conditioned on the amount of available simulated data. The seemingly simple case study introduced in this work might serve as a prototype for studying a wide variety of applications, including adsorption, wetting, and capillarity, to name a few.

Inhomogeneous solvation theory (IST) and classical density functional theory (DFT) each provide a framework for relating distribution functions of solutions to their thermodynamic properties. As reviewed in this work, both IST and DFT can be formulated in a way that use two "end point" simulations, one of the pure solvent and the other of the solution, to determine the solute chemical potential and other thermodynamic properties of the solution and of subvolumes in regions local to the solute containing hydrating waters. In contrast to IST, where expressions for the excess energy and entropy of solution are the object of analysis, in the DFT end point formulation of the problem, the solute-solvent potential of mean force (PMF) plays a central role. The indirect part of the PMF corresponds to the lowest order (1-body) truncation of the IST expression. Because the PMF is a free energy function, powerful numerical methods can be used to estimate it. We show that the DFT expressions for the solute excess chemical potential can be written in a form which is local, involving integrals only over regions proximate to the solute. The DFT end point route to estimating solvation free energies provides an alternative path to that of IST for analyzing solvation effects on molecular recognition and conformational changes in solution, which can lead to new insights. In order to illustrate the kind of information that is contained in the solute-solvent PMF, we have carried out simulations of β-cyclodextrin in water. This solute is a well studied "host" molecule to which "guest" molecules bind; host-guest systems serve as models for molecular recognition. We illustrate the range of values the direct and indirect parts of the solute-solvent PMF can have as a water molecule is brought to the interface of β-cyclodextrin from the bulk; we discuss the "competition" between these two terms, and the role it plays in molecular recognition.

In general, materials with linear memory constitutive relations are characterized by a relaxation function. This leads to a situation where the free energy for most materials with memory is not unique. There is a convex set of free energy functionals with a minimum and a maximum element. An alternative procedure is proposed which characterizes a material by the kernel of the rate of dissipation functional. Using some recent results, we find that a unique free energy and relaxation function may then be deduced. An example is given for discrete spectrum materials. Also, the new results are used to show that a previously derived general representation of rate of dissipation and free energy functionals is not complete, in the sense that there are valid functionals which cannot be described by this general formula.

We study a phase-field variational model for the solvaiton of charged molecules with an implicit solvent. The solvation free-energy functional of all phase fields consists of the surface energy, solute excluded volume and solute-solvent van der Waals dispersion energy, and electrostatic free energy. The surface energy is defined by the van der Waals--Cahn--Hilliard functional with squared gradient and a double-well potential. The electrostatic part of free energy is defined through the electrostatic potential governed by the Poisson--Boltzmann equation in which the dielectric coefficient is defined through the underlying phase field. We prove the continuity of the electrostatics---its potential, free energy, and dielectric boundary force---with respect to the perturbation of dielectric boundary. We also prove the $\Gamma$-convergence of the phase-field free-energy functionals to their sharp-interface limit, and the equivalence of the convergence of total free energies to that of all individual parts of free energy. We finally prove the convergence of phase-field forces to their sharp-interface limit. Such forces are defined as the negative first variations of the free-energy functional; and arise from stress tensors. In particular, we obtain the force convergence for the van der Waals--Cahn--Hilliard functionals with minimal assumptions.

The aim of this work is to determine what free energy functionals are expressible as quadratic forms of the state functional \(I^t\) which is discussed in earlier papers. The single integral form is shown to include the functional \(\psi _F\) proposed a few years ago, and also a further category of functionals which are easily described but more complicated to construct. These latter examples exist only for certain types of materials. The double integral case is examined in detail, against the background of a new systematic approach developed recently for double integral quadratic forms in terms of strain history, which was used to uncover new free energy functionals. However, while, in principle, the same method should apply to free energies which can be given by quadratic forms in terms of \(I^t\), it emerges that this requirement is very restrictive; indeed, only the minimum free energy can be expressed in such a manner.

The free energy functions for electron transfer reactions in solution are explored using a previously developed microscopic simulation approach that provides a clear definition of these functions and the variable (the reaction coordinate) used as their argument. The issue of whether the curvatures of the two functions (which correspond to states with nonpolar and polar solutes) are different is given special attention. It is found, in contrast to some previous suggestions, that the curvatures of the two functions are quite similar, even when one would expect differences due to dielectric saturation effects, and that Marcus’ approximation (and, in fact, the linear response theory inherent in this approximation) provides a valid description of the solvent’s role in electron transfer reactions over a wide range of conditions. The present study demonstrates that direct simulations of the reactant and product states do not provide the data needed for determination of the free energy functions in high energy regions. On the other hand, a free energy perturbation approach developed in our early studies does allow one to obtain statistically significant data for these regions of high free energy. Inclusion of the solvent electronic polarizabilities in the potential results in a marked lowering of the solvent reorganization energy, but does not change the ratio of the curvatures of the free energy functions for the nonpolar and polar states. Finally, the effect of solvent tunneling is evaluated with the dispersed polaron method and shown to be rather small, and it is argued that the observed deviations from the classical Marcus theory in the inverted region are due to quantum mechanical tunneling associated with the solute vibrations, rather than to differences in the curvatures of the solvent free energy functions.

Protein molecules in solution have a broad distribution of enthalpy states. A good approximation to the distribution function for enthalpy states can be calculated, using the maximum-entropy method, from the moments of the distribution that, in turn, are obtained from the experimental temperature dependence of the heat capacity. In the present paper, we show that the enthalpy probability distribution can then be formulated in terms of a free energy function that gives the free energy of the protein corresponding to a particular value of the enthalpy. By the location of the minima in this function, the free energy distribution graphically indicates the most probable values of the enthalpy for the protein. We find that the behavior of the free energy functions for proteins falls somewhere between two different cases: a two-state like function with two minima, the relative levels of the two states changing with temperature; and, a single-minimum function where the position of the minimum shifts to higher enthalpy values as the temperature is increased. We show that the temperature dependence of the free energy function can be expressed in terms of a central free energy distribution for a given, fixed temperature (which is most conveniently chosen as the temperature of the maximum in the heat capacity). The nature of this central free energy function for a given protein thus yields all of the thermodynamic behavior of that protein over the temperature range of the denaturation process.

Atomistic Modeling of Excluded-Volume and Soft Interactions in Crowded Environments

Monte Carlo simulations are used to study the translocation of a polymer into and out of an ellipsoidal cavity through a narrow pore. We measure the polymer free energy F as a function of a translocation coordinate, s, defined to be the number of bonds that have entered the cavity. To study polymer insertion, we consider the case of a driving force acting on monomers inside the pore, as well as monomer attraction to the cavity wall. We examine the changes to F(s) upon variation in the shape anisometry and volume of the cavity, the polymer length, and the strength of the interactions driving the insertion. For athermal systems, the free energy functions are analyzed using a scaling approach, where we treat the confined portion of the polymer to be in the semi-dilute regime. The free energy functions are used with the Fokker-Planck (FP) equation to calculate mean translocation times, as well as translocation time distributions. We find that both polymer ejection and insertion are faster for ellipsoidal cavities than for spherical cavities. The results are in qualitative agreement with those of a Langevin dynamics study in the case of ejection but not for insertion. The discrepancy is likely due to out-of-equilibrium conformational behaviour that is not accounted for in the FP approach.

We present the pressure tensor of an associating Lennard-Jones (LJ) fluid confined in a spherical cavity of hard wall, where a high-order density correlation has been taken into account. To give the two-body association potential for calculating the pressure tensor, an angle-average of site-site attraction over all orientations of two particles is performed. Furthermore, the classical density functional theory is employed to obtain the density profile of the confined fluid, by which the normal and tangential pressure profiles are illustrated under various conditions to show the dependence of the pressure tensor on the association strength, number of associating sites, radius of cavity, and bulk density. As an application, the corresponding surface tension is calculated. It is shown that under a strong association interaction (both association strength and the number of associating sites are large), the pressure profiles are depleted from the wall of the cavity instead of the oscillatory behavior under a weak association interaction. Such a tendency is mainly determined by the competition between association interaction and excluded volume interaction. Therefore, the aggregation state and related properties of an associating LJ fluid within a confinement of nanoscale can be efficiently regulated by the association interaction.

The octanol-water partition coefficient is an important physical-chemical characteristic widely used to describe hydrophobic/hydrophilic properties of chemical compounds. The partition coefficient is related to the transfer free energy of a compound from water to octanol. Here, we introduce a new protocol for prediction of the partition coefficient based on the statistical-mechanical, 3D-RISM-KH molecular theory of solvation. It was shown recently that with the compound-solvent correlation functions obtained from the 3D-RISM-KH molecular theory of solvation, the free energy functional supplemented with the correction linearly related to the partial molar volume obtained from the Kirkwood-Buff/3D-RISM theory, also called the "universal correction" (UC), provides accurate prediction of the hydration free energy of small compounds, compared to explicit solvent molecular dynamics [ Palmer , D. S. ; J. Phys.: Condens. Matter 2010 , 22 , 492101 ]. Here we report that with the UC reparametrized accordingly this theory also provides an excellent agreement with the experimental data for the solvation free energy in nonpolar solvent (1-octanol) and so accurately predicts the octanol-water partition coefficient. The performance of the Kovalenko-Hirata (KH) and Gaussian fluctuation (GF) functionals of the solvation free energy, with and without UC, is tested on a large library of small compounds with diverse functional groups. The best agreement with the experimental data for octanol-water partition coefficients is obtained with the KH-UC solvation free energy functional.

Calculating the equilibrium composition of systems under given conditions (usually, temperature and pressure) is a core in the majority of the existing computer models of geochemical processes. Recently, Shvarov [1] investigated the thermodynamic consistency of the models of nonideal aqueous solutions that are used in geochemical calculations. He distinguished necessary and sufficient conditions for the thermodynamic consistency of the physicochemical models of real (nonideal) solutions and presented particular examples of errors appearing when these criteria are not met. In addition, Shvarov [1] claimed that the models that are currently used for concentrated electrolyte solutions are thermodynamically inconsistent. Such a pessimistic conclusion is related to a large extent to the fact that Shvarov’s [1] analysis was limited to the Debye‐Huckel model, which is applicable to diluted aqueous solutions, and its empirical modifications aimed at extending to concentrated aqueous solutions. Indeed, such models are thermodynamically inconsistent. However, models have been developed allowing adequate description of the behavior of concentrated aqueous solutions. One of them is the wellknown Pitzer model [2, 3] proposed for the description of the properties of dense electrolyte solutions in a wide range of temperatures, pressures, and compositions. The Pitzer model was developed on the basis of the Gibbs method and is, consequently, thermodynamically sound. If correct numerical algorithms are employed, the Pitzer model allows for the calculation of chemical equilibria in water‐ salt systems within a wide concentration range. The Pitzer model is an example of the rigorous statistical mechanical description of a multicomponent fluid system and the calculation of free energy by the Gibbs method. The derivation of an expression for free energy in the Pitzer model was described in detail in [2, 3], and only its basic principles are shortly discussed here. It is based on the approach developed by McMillan and Mayer [4], who described interactions between dissolved species within the mean field approximation and ignored the detailed description of interactions between individual particles. It should be emphasized that these approximations concern interactions between particles, i.e., a microscopic Hamiltonian. The free energy function is further derived in accordance with the statistical method of Gibbs using the methods of group expansion of Mayer [5, 6] and divergence elimination [6, 7] for the calculation of the configuration integrals of the long-range potentials of electrostatic interactions. Eventually, the excess Gibbs free energy of solution is given in the Pitzer model, after transition from volume concentrations to molalities, as [2]

Empirical free energy as a target function in docking and design: application to HIV-1 protease inhibitors

We describe a method for predicting the three-dimensional (3-D) structure of proteins from their sequence alone. The method is based on the electrostatic screening model for the stability of the protein main-chain conformation. The free energy of a protein as a function of its conformation is obtained from the potentials of mean force analysis of high-resolution x-ray protein structures. The free energy function is simple and contains only 44 fitted coefficients. The minimization of the free energy is performed by the torsion space Monte Carlo procedure using the concept of hierarchic condensation. The Monte Carlo minimization procedure is applied to predict the secondary, super-secondary, and native 3-D structures of 12 proteins with 28-110 amino acids. The 3-D structures of the majority of local secondary and super-secondary structures are predicted accurately. This result suggests that control in forming the native-like local structure is distributed along the entire protein sequence. The native 3-D structure is predicted correctly for 3 of 12 proteins composed mainly from the alpha-helices. The method fails to predict the native 3-D structure of proteins with a predominantly beta secondary structure. We suggest that the hierarchic condensation is not an appropriate procedure for simulating the folding of proteins made up primarily from beta-strands. The method has been proved accurate in predicting the local secondary and super-secondary structures in the blind ab initio 3-D prediction experiment.

Monte Carlo simulations are used to study the conformational behavior of a semiflexible polymer confined to cylindrical and conical channels. The channels are sufficiently narrow that the conditions for the Odijk regime are marginally satisfied. For cylindrical confinement, we examine polymers with a single knot of topology 3_{1}, 4_{1}, or 5_{1}, as well as unknotted polymers that are capable of forming S loops. We measure the variation of the free energy F with the end-to-end polymer extension length X and examine the effect of varying the polymer topology, persistence length P, and cylinder diameter D on the free-energy functions. Similarly, we characterize the behavior of the knot span along the channel. We find that increasing the knot complexity increases the typical size of the knot. In the regime of low X, where the knot/S-loop size is large, the conformational behavior is independent of polymer topology. In addition, the scaling properties of the free energy and knot span are in agreement with predictions from a theoretical model constructed using known properties of interacting polymers in the Odijk regime. We also examine the variation of F with the position of a knot in conical channels for various values of the cone angle α. The free energy decreases as the knot moves in a direction where the cone widens, and it also decreases with increasing α and with increasing knot complexity. The behavior is in agreement with predictions from a theoretical model in which the dominant contribution to the change in F is the change in the size of the hairpins as the knot moves to the wider region of the channel.