Knowledge of the Free Energy Landscape topology is the essential key to understanding many biochemical processes. The determination of the conformers of a protein and their basins of attraction takes a central role for studying molecular isomerization reactions. In this work, we present a novel framework to unveil the features of a Free Energy Landscape answering questions such as how many meta-stable conformers there are, what the hierarchical relationship among them is, or what the structure and kinetics of the transition paths are. Exploring the landscape by molecular dynamics simulations, the microscopic data of the trajectory are encoded into a Conformational Markov Network. The structure of this graph reveals the regions of the conformational space corresponding to the basins of attraction. In addition, handling the Conformational Markov Network, relevant kinetic magnitudes as dwell times and rate constants, or hierarchical relationships among basins, completes the global picture of the landscape. We show the power of the analysis studying a toy model of a funnel-like potential and computing efficiently the conformers of a short peptide, dialanine, paving the way to a systematic study of the Free Energy Landscape in large peptides.
Polymers and, more specifically, proteins, show complex behavior at the cellular system level, e. g. in protein-protein interaction networks [1], and also at the individual level, where proteins show a large degree of multistability: a single protein can fold in different conformational states [2]–[4]. As a complex system [5], [6], the dynamics of a protein cannot be understood by studying its parts in isolation, instead, the system must be analyzed as a whole. Tools able to represent and handle the information of the entire picture of a complex system are thus necessary. Complex network theory [7], [8] has proved to be a powerful tool used in seemingly different biologically-related fields such as the study of metabolic reactions, ecological and food webs, genetic regulatory systems and the study of protein dynamics [7]. In this latter context, diverse studies have analyzed the conformational space of polymers and proteins making use of network representations [9]–[12], where nodes account of polymer conformations. Additionally, some studies have tried to determine the common and general properties of these conformational networks [13], [14] looking at magnitudes such as clustering coefficient, cyclomatic number, connectivity, etc. Recently, trying to decompose the network in modules corresponding to the free energy basins, the use of community algorithms over these conformational networks have been proposed [15]. Although this approach has opened a promising path for the analysis of Free Energy Landscapes (FEL), the community based description of the network leads to multiple characterizations of the FEL and thus it is difficult to establish a clear map from the communities found to the basins of the FEL. A similar approach, commonly used to analyze the complex dynamics, is the construction of Markovian models. Markovian state models let us treat the information of one or several trajectories of molecular dynamics (MD) as a set of conformations with certain transition probabilities among them [9], [16], [17]. Therefore, the time-continuous trajectory turns into a transition matrix, offering global observables as relaxation times and modes. In [16]–[18] the use of Markovian models is proposed with the aim of detecting FEL meta-stable states. However, the above approaches to analyze FELs of peptides involves extremely large computational cost: either general community algorithms or large transition matrices. Finally, other strategies to characterize the FEL that have successfully helped to understand the physics of biopolymers, are based on the study of the Potential Energy Surface (PES) [3], [4], [19]–[21]. The classical transition-state theory [22] allows us to project the behavior of the system at certain temperature from the knowledge of the minima and transition states of the PES. This approach entails some feasible approximations, such as harmonic approximation to the PES, limit of high damping, assumption of high barriers, etc. These approximations could be avoided working directly from the MD data. In this article we make a novel study of the FEL capturing its mesoscopic structure and hence characterizing conformational states and the transitions between them. Inspired by the approaches presented in [12], [15] and [16], [17], we translate a dynamical trajectory obtained by MD simulations into a Conformational Markov Network. We show how to efficiently handle the graph to obtain, through its topology, the main features of the landscape: conformers and their basins of attraction, dwell times, rate constants between the conformational states detected and a coarse-grained picture of the FEL. The framework is shown and validated analyzing a synthetic funnel-like potential. After this, the terminally blocked alanine peptide (Ace-Ala-Nme) is studied unveiling the main characteristics of its FEL.
First, we encode a trajectory of a stochastic MD simulation into a network: the CMN. This map will allow us to use the tools introduced henceforth to analyze a specific dynamics of complex systems such as biopolymers. Up to now, we have illustrated the conversion of molecular dynamics data into a graph (the CMN). Now, we show how to efficiently obtain the thermo-statistical data from the mesoscopic description of the CMN.
The alanine dipeptide, or terminally blocked alanine peptide (Ace-Ala-Nme, Figure 3A), is the most simple “biological molecule” that exhibits the common features shown by larger biomolecules. Despite of its simplicity, this system has more than one long-life conformational state with different transition pathways. Since the first attempt by Rossky and Karplus [37] to model this dipeptide solvated, this system has been widely studied in theoretical works [38]–[41]. The alanine dipeptide has been also the appropriate molecule to test tools to explore the FEL [15], [16], [42] and, specifically, to study reaction coordinates [41], [43]. The alanine dipeptide has two slow degrees of freedom, the rotatable bonds () and () (see Figure 3A). The FEL projected onto these dihedral angles let us identify the conformational states that characterize the geometry of biopolymers, namely: alpha helix right-handed (), alpha helix left-handed (), beta strands (,), etc. The number of local minima in the (,) space depends on the effective potential model used to simulate the system. Up to date, electronic structure calculations have identified a total of nine different conformers [44]. Regarding MD simulations different conformational states have been observed: (i) using classical force fields with explicit solvent up to six conformers are detected [16], [38], [39], (ii) at least four stable states by using implicit solvent [15], [38], [40], and (iii) two stable conformers in vacuum conditions [38], [41]. On the other hand, since the angles and seem appropriate to distinguish the metastable states, the kinetics between them is not accurately described with this choice of reaction coordinates, the solvent coordinates and/or other internal degrees of freedom must be taken into account [41], [43]. We have used SSD algorithm to detect the local minima and their corresponding basins for this molecule in the space. For this purpose, a Langevin MD simulation of 250 ns has been performed at a temperature of 400 K (see Text S3 for further details). Additionally the CMN has been built dividing the Ramachandran plot into cells of surface 9°×9° (40×40) and taking dialanine conformations at time intervals of. The resulting CMN have a total of nodes and directed links. The algorithm applied to the CMN network reveals 6 basins. Figure 3B shows the resulting network where nodes belonging to the same basin take the same color. Bringing back this information to the Ramachandran map, these 6 sets of nodes define 6 regions represented in Figure 3C. To better illustrate this division, other representation, where each region has a different color, is shown. By comparing with previous studies on this molecule, we identify the regions in orange, red, yellow and pink with conformers, , and respectively [38], [45]. Besides, region green corresponds to conformer and the blue one to [16], [45]. Remarkably, the basin (one of the less populated state) has been visited 1155 times with a mean stay time of 4. 20 ps. We now look at the coarse-grained picture of the FEL by describing the properties of the 6 basins detected. The different weights of the basins are related to the free energy of the corresponding conformational macro-states. In Table 1 these energy differences are shown, taking the most populated basin as the energy reference. The lowest free energy basins correspond to configurations with φ≤0° (, , and), whereas the two other conformers located in the region φ≥0° have the highest free energy but the largest dwell time. Moreover, we have also analyzed the trapping efficiency of each basin by reporting the mean escape time () as well in Table 1. The FEL can be represented as a dendogram, see Figure 4, where the hierarchical map of the conformers based on Free Energy gives at first glance a global picture of the landscape. Remarkably, the conformer, despite of having one of the highest free energy, looks like the metastable state with longest life. This result is supported by the values of Mean Escape Time shown in Table 1. The alanine dipeptide has been also studied because of its “fast” isomerization and a “slow” transition. Our coarse-grained picture of the FEL also allows us to extract information about these transitions. In Table 2 we show some of the relevant characteristic transition times from a basin a to an adjacent basin b, . [The whole data is shown in the Text S3. ] Transitions between basins with the same sign of are remarkably faster (e. g and). While slow transitions are observed for those hops crossing the line φ = 0° (and), showing them as rare events. Instead, the alanine dipeptide finds more easily paths to go to φ≥0° conformers through and by crossing φ = 180°. To round off the description of the FEL, the dendogram corresponding to the temporal hierarchy is shown in Figure 5. From the figure, it becomes clear that the behavior of the dialanine depends on the time scale used for its observation. Again, the same two different sets of conformers are distinguished from this hierarchy. Additionally, the global minimum conformer is reached in around 100 ps from any basin. Finally, the magnitudes computed here for the alanine dipeptide would allow to construct a first-order kinetic model of 6 coupled differential equations as Eq. (6) (assuming equilibrium intra-basin). This model contains the same information as the kinetic model by Chekmarev et al. for the irreversible transfer of population from [40].
Hierarchical landscapes characterize the dynamical behavior of proteins, which in turn depends on the relation between the topology of the basins, their transitions paths and the kinetics over energy barriers. The CMN analysis of trajectories generated by MD simulations is a powerful tool to explore complex FELs. In this article, we have proposed how to deal with a CMN to unveil the structure of the FEL in a straightforward way and with a remarkable efficiency. The analysis presented here is based on the physical concept of basin of attraction, making possible the study of the conformational structure of peptides and the complete characterization of its kinetics. Note that this has been done without the estimation of the volume of each conformational macro-state in the coordinates space and without the ‘a priori’ knowledge of the saddle points or the transition paths from a local minimum to another. On the other hand, the framework introduced in the article provides us with a quantitative description of the dialanine' s FEL, coming up directly from a MD dynamics at certain temperature. The peptide explores its landscape building the corresponding CMN and the success of extracting the relevant information is up to the ability of dealing with it. Neither the FE basins were defined by the unique criterion of clustering conformations with a geometrical distance [46], nor the rate constants were projected from the potential energy surface [19], [20]. Moreover, the conformers and their properties were computed from the MD with the only limitation of the discretization of time and space. Although we have applied the method to low dimensional landscapes, we expect that high dimensional systems could be also studied, by the combination of this technique with the usual methods to reduce the effective degrees of freedom (like principal component analysis or essential dynamics). In conclusion, the large amount of information obtained by working with the CMN, its potential application to any peptide with a large number of monomers, and the possibility of performing the analysis on top of CMN constructed via several short MD simulations [47], make the approach presented here a promising way to describe the FEL of a protein.
A complete description of complex polymers, such as proteins, includes information about their structure and their dynamics. In particular it is of utmost importance to answer the following questions: What are the structural conformations possible? Is there any relevant hierarchy among these conformers? What are the transition paths between them? These and other questions can be addressed by analyzing in an efficient way the Free Energy Landscape of the system. With this knowledge, several problems about biomolecular reactions (such as enzymatic activity, protein folding, protein deposition diseases, etc.) can be tackled. In this article we show how to efficiently describe the Free Energy Landscape for small and large peptides. By mapping the trajectories of molecular dynamics simulations into a graph (the Conformational Markov Network) and unveiling its structural organization, we obtain a coarse grained description of the protein dynamics across the Free Energy Landscape in terms of the relevant kinetic magnitudes of the system. Therefore, we show the way to bridge the gap between the microscopic dynamics and the macroscopic kinetics by means of a mesoscopic description of the associated Conformational Markov Network. Along this path the compromise between the physical nature of the process and the magnitudes that characterize the network is carefully kept to assure the reliability of the results shown.
computational biology/molecular dynamics
computational biology/macromolecular structure analysis
biophysics/theory and simulation
biophysics/protein folding
Exploring the Free Energy Landscape: From Dynamics to Networks and Back