Minimum Free Energy Path of Ligand-Induced Transition in Adenylate Kinase
Biological functions of proteins are mediated by dynamical processes occurring on complex energy landscapes 1 ., These processes frequently involve large conformational transitions between two or more metastable states , induced by an external perturbation such as ligand binding 2 ., Time scales of the conformational transition are frequently of order microseconds to seconds ., To characterize such slow events in molecular dynamics ( MD ) trajectories , the free energy profile or the potential of mean force ( PMF ) along a reaction coordinate must be identified ., In particular , the identification of the transition state ensemble ( TSE ) enables the barrier-height to be evaluated , and the correct kinetics would be reproduced if there is only a single dominant barrier ., However , for proteins with many degrees of freedom , finding an adequate reaction coordinate and identifying the TSE is a challenging task placing high demands on computational resources ., The finite-temperature string method 3 , 4 , and the on-the-fly string method 5 find a minimum free energy path ( MFEP ) from a high-dimensional space ., Given a set of collective variables describing a conformational change , the MFEP is defined as the maximum likelihood path along the collective variables ., The MFEP is expected to lie on the center of reactive trajectories and contains only important transitional motions 4 ., Furthermore , since the MFEP approximately orthogonally intersects the isocommittor surfaces ( the surfaces of constant committor probability in the original space ) 4 , the TSE can be identified as the intersection with the isocommittor surface with probability of committing to the product ( or the reactant ) =\u200a1/2 ., The methods and MFEP concepts have been applied to various molecular systems 4–6 including protein conformational changes 7–9 ., With regard to high-dimensional systems like proteins , the quality of the MFEP ( whether it satisfies the above-mentioned properties ) is particularly sensitive to choice of the collective variables ., The collective variables should be selected such that their degrees of freedom are few enough to ensure a smooth free energy surface; at the same time they should be sufficiently many to approximate the committor probability 4 , 9 ., To resolve these contrary requirements , effective dimensional reduction is required ., Large conformational transitions of proteins , frequently dominated by their domain motions , can be well approximated by a small number of large-amplitude principal modes 2 , 10 ., This suggests that the use of the principal components may be the best choice for approximating the committor probability with the fewest number of variables for such large conformational transitions involving domain motions ., A further advantage is the smoothness of the free energy landscape in the space of the large-amplitude principal components ., If the curvature of the MFEP is large , the MFEP may provide a poor approximation to the isocommittor surface since the flux can occur between non-adjacent structures along the path 9 ., The selection of the large-amplitude principal components as the collective variables would maintain the curvature of the MFEP sufficiently small ., Here , we conducted preliminary MD simulations around the two terminal structures of the transition and performed a principal component analysis to obtain the principal components ( see Materials and Methods for details ) ., Following selection of a suitable MFEP , determination of the PMF and characterization of the physical quantities along the MFEP are needed to understand an in-depth mechanism of the transition ., Although the finite-temperature string method yields a rigorous estimate of the gradient of the PMF under a large coupling constant with the collective variables 3–5 ( see Materials and Methods ) , errors in the estimates of the gradients and in the tangential directions of the pathway tend to accumulate during the integration process ., To accurately quantify the PMF and the averages of various physical quantities in a multi-dimensional space , we utilized another statistical method , the multi-state Bennett acceptance ratio ( MBAR ) method 11 , which provides optimal estimates of free energy and other average quantities along the MFEP ., Here , we applied the above proposed methods to the conformational change in Escherichia coli adenylate kinase ( AK ) , the best-studied of enzymes exhibiting a large conformational transition 12–23 ., AK is a ubiquitous monomeric enzyme that regulates cellular energy homeostasis by catalyzing the reversible phosphoryl transfer reaction: ATP+AMP↔2ADP ., According to the analysis of the crystal structures by the domain motion analysis program DynDom 24 , AK is composed of three relatively rigid domains ( Figure 1 ) ; the central domain ( CORE: residues 1–29 , 68–117 , and 161–214 ) , an AMP-binding domain ( AMPbd: 30–67 ) , and a lid-shaped ATP-binding domain ( LID: 118–167 ) ., Inspection of the crystal structures suggests that , upon ligand binding , the enzyme undergoes a transition from the inactive open form to the catalytically competent closed structure 25 ( Figure 1 ) ., This transition is mediated by large-scale closure motions of the LID and AMPbd domains insulating the substrates from the water environment , while occluding some catalytically relevant water molecules ., The ATP phosphates are bound to the enzyme through the P-loop ( residues 7–13 ) , a widely-distributed ATP-binding motif ., The interplay between AKs dynamics and function has been the subject of several experimental studies ., 15N NMR spin relaxation studies have revealed that the LID and AMPbd domains fluctuate on nanosecond timescales while the CORE domain undergoes picosecond fluctuations 12 , 13 ., The motions of these hinge regions are highly correlated with enzyme activity 14 ., In particular , the opening of the LID domain , responsible for product release , is thought to be the rate-limiting step of the catalytic turnover 14 ., Recent single-molecule Förster resonance energy transfer ( FRET ) experiments have revealed that the closed and open conformations of AK exist in dynamic equilibrium even with no ligand present 15 , 16 , and that the ligands presence merely changes the populations of open and closed conformations ., This behavior is reminiscent of the population-shift mechanism 26 rather than the induced-fit model 27 , in which structural transitions occur only after ligand binding ., The population-shift like behaviour of AK has also been supported by simulation studies 17–20 ., Lou and Cukier 17 , Kubitzki and de Groot 18 , and Beckstein et al . 19 employed various enhanced sampling methods to simulate ligand-free AK transitions ., Arora and Brooks 20 applied the nudged elastic band method in the pathway search for both ligand-free and ligand-bound forms ., These studies showed that , while the ligand-free form samples conformations near the closed structure 17–20 , ligand binding is required to stabilize the closed structure 20 ., Despite the success of these studies based on all-atom level models , atomistic details of the transition pathways , including the structures around the TSE , have not been fully captured yet ., In this study , we successfully evaluated the MFEP for both ligand-free and ligand-bound forms of AK using the on-the-fly string method , and calculated the PMF and the averages of various physical quantities using the MBAR method ., Our analysis elucidates an in-depth mechanism of the conformational transition of AK ., The MFEPs for apo and holo-AKs , and their PMFs , were obtained from the string method and the MBAR method , respectively ( see Videos S1 and S2 ) ., The MFEPs were calculated using the same 20 principal components selected for the collective variables ., The holo-AK calculations were undertaken with the bisubstrate analog inhibitor ( Ap5A ) as the bound ligand without imposing any restraint on the ligand ., Figures 2A and 2B show the MBAR estimates of the PMFs along the images of the MFEP ( the converged string at 12 ns in Figures 2A and 2B ) for apo and holo-AK , respectively ., Here , the images on the string are numbered from the open ( ; PDBid: 4ake 28 ) to the closed conformation ( ; PDBid: 1ake 29 ) ., These terminal images were fixed during the simulations to enable sampling of the conformations around the crystal structures ., In the figures , the convergence of the PMF in the string method process is clearly seen in both systems ., Convergence was also confirmed by the error estimates ( Figure S1 ) , and by the root-mean-square displacement ( RMSD ) of the string from its initial path ( Figure S2 ) ., The PMF along the MFEP reveals a broad potential well on the open-side conformations of apo-AK , suggesting that the open form of AK is highly flexible 20 ., This broad well is divided into two regions , the fully open ( ) and partially closed states ( , encircled ) by a small PMF barrier ., In holo-AK ( Figure 2B ) , the MFEP exhibits a single substantial free energy barrier ( ) between the open and closed states , which does not appear in the initial path ., This barrier will be identified as the transition state below ., It is shown in the PMF along the MFEP that the closed form ( tightly binding the ligand ) is much more stable than the open form with loose binding to the ligand ( large fluctuations of the ligand will be shown later ) ., To characterize the MFEP in terms of the domain motions , the MFEP was projected onto a space defined by two distances from the CORE domain , the distance to the LID domain and the distance to the AMPbd domain ( the distance between the mass centers of atoms for the two domains; Figures 2C and 2D ) ., The PMF was also projected onto this space ., The comparison of the two figures shows that ligand binding changes the energy landscape of AK , suggestive that this is not a simple population-shift mechanism ., In apo-AK , the motions of the LID and AMPbd domains are weakly correlated , reflecting the zipper-like interactions on the LID-AMPbd interface 19 ., The MFEP clearly indicates that the fully closed conformation ( ) involves the closure of the LID domain followed by the closure of the AMPbd domain ., The higher flexibility of the LID domain has been reported in previous studies 17 , 19 , 20 ., In holo-AK , the pathway can be described by two successive scenarios , that is , the LID-first-closing followed by the AMPbd-first-closing ., In the open state ( ) , the MFEP is similar to that of apo-AK , revealing that LID closure occurs first ., In the closed state ( ) , however , the AMPbd closure precedes the LID closure ., This series of the domain movements was also identified by the domain motion analysis program DynDom 24 ( Figure S3 ) ., It is known in the string method that the convergence of the pathway is dependent on the initial path ., In order to check whether the MFEP obtained here is dependent on the initial path or not , we conducted another set of the calculations for apo-AK by using a different initial path , which has an AMPbd-first-closing pathway , opposed to the LID-first-closing pathway shown above ., If the LID and AMPbd domains move independently of each other , it is expected that LID-first-closing and AMPbd-first-closing pathways are equally stable ., Despite this initial setup , however , our calculation again showed the convergence toward the LID-first-closing pathway ( see Figure S4 ) ., As described above , this tendency of the pathways would be due to the reflection of the highly flexible nature of the LID domain ., Furthermore , in order to check whether the samples around the MFEP are consistent with the experiments , we compared the PMF as a function of the distance between the Cα atoms of Lys145 and Ile52 with the results of the single-molecule FRET experiment by Kern et al . 16 ( see Figure S5 ) ., The PMF was calculated by using the samples obtained by the umbrella samplings around the MFEP ., In the figure , the stable regions of the PMF for holo-AK are highly skewed toward the closed form , and some population toward the partially closed form was also observed even for apo-AK , which is consistent with the histogram of the FRET efficiency 16 ., To more clearly illustrate the energetics along the MFEP in terms of the domain motions , we separately plot the PMF as a function of the two inter-domain distances defined above ( Figures 3A and 3B ) ., We observe that the PMF of apo-AK has a double-well profile for the LID-CORE distance ( indicated by the blue line in Figure 3A ) , whereas the PMF in terms of the AMPbd-CORE distance is characterized by a single-well ( Figure 3B ) ., The single-molecule FRET experiments monitoring the distances between specific residue pairs involving the LID domain ( LID-CORE ( Ala127-Ala194 ) 15 and LID-AMPbd ( Lys145- Ile52 ) 16 ) revealed the presence of double-well profiles in the ligand-free form ., On the other hand , an electron transfer experiment probing the distance between the AMPbd and CORE domains ( Ala55-Val169 ) 30 showed only that the distance between the two domains decreased upon ligand binding ., Considering the PMF profiles in the context of these experimental results , we suggest that the partially closed state ( ) in apo-AK ( Figure 2A ) can be ascribed to the LID-CORE interactions but not to the AMPbd-CORE interactions ., To elucidate the origin of the stability of the partially closed state , we monitored the root mean square fluctuations ( RMSF ) of the atoms along the MFEP ( see Materials and Methods for details ) ., Figures 3C and 3D show the RMSF along the MFEP for apo and holo-AK , respectively ., In apo-AK ( Figure 3C ) , large fluctuations occur in the partially closed state ( ) around the LID-CORE hinge regions ( residue 110–120 , and 130–140 ) and the P-loop ( residue 10–20 ) ., It has been proposed , in the studies of AK using coarse-grained models , that “cracking” or local unfolding occurs due to localized strain energy , and that the strained regions reside in the LID-CORE hinge and in the P-loop 21 , 23 ., Our simulation using the all-atom model confirmed the existence of “cracking” in the partially closed state , and provided an atomically detailed picture of this phenomenon ., The average structures around the partially closed state revealed that , in the open state , a highly stable Asp118-Lys136 salt bridge is broken by the strain induced by closing motion around ( Figure S6A ) ., This salt bridge has been previously proposed to stabilize the open state while imparting a high enthalpic penalty to the closed state 18 ., Breakage of the salt bridge releases the local strain and the accompanying increases in fluctuation may provide compensatory entropy to stabilize the partially closed state ., A similar partially closed state of the LID domain was also found by the work of Lou and Cukier 31 in which they performed all-atom MD simulation of apo-AK at high temperature ( 500 K ) condition ., In holo-AK , both of the LID-CORE and AMPbd-CORE distances exhibit double-well profiles ( indicated by the red lines in Figures 3A and 3B ) , separating the closed from the open state ., The breakage of the 118–136 salt bridge at around is not accompanied by “cracking” of the hinge region ( Figure 3D ) ., Instead , the hinge region is stabilized by binding of ATP ribose to Arg119 and His134 ( Figure S6B ) , leading to a smooth closure of the LID domain ., This suggests that one role of the salt bridge breakage is rearrangement of the molecular interactions to accommodate ATP-binding 32 ., P-loop fluctuations are also suppressed in holo-AK ( Figure 3D ) ., Consistent with our findings , reduced backbone flexibilities in the presence of Ap5A were reported in the above-mentioned NMR study 13 ., The origin of the double-well profile in holo-AK was investigated via the ligand-protein interactions ., The motion of the ligand along the MFEP was firstly analyzed by focusing on the AMP adenine dynamics , since the release of the AMP moiety from the AMP-binding pocket was observed in the open state ., It is again emphasized that the ligand is completely free from any restraint during the simulations ., PCA was performed for the three-dimensional Cartesian coordinates of the center of mass of AMP adenine , and the coordinates were projected onto the resultant 1st PC in Figure 4A ., The AMP adenine is observed to move as much as 10 Å in the open state ( ; Figure 4B ) , while it is confined to a narrow region of width 1–2 Å ( the binding pocket ) in the closed state ( ) ., Such a reduction of the accessible space of the AMP adenine might generate a drastic decrease in entropy or an increase in the PMF barrier of the open-to-closed transition ., Furthermore , close inspection of the PMF surface reveals the existence of a misbinding event at ( Figure 4B ) , in which the AMP ribose misbinds to Asp84 in the CORE domain , and is prevented from entering the AMP-binding pocket ., This event further increases the barrier-height of the transition ., The MFEP revealed that AMP adenine enters the AMP-binding pocket around , as indicated by a rapid decrease in the accessible area ( Figure 4A ) ., This event is well correlated with the position of the PMF barrier along the MFEP ( Figure 2B ) ., This coincidence between the binding process and the domain closure suggests that the two processes are closely coupled ., Before analyzing the situation in detail , however , it is necessary to assess whether the observed PMF barrier around ( Figure 2B ) corresponds to a TSE , because the PMF barrier is not necessarily a signature of dynamical bottleneck in high-dimensional systems 33 ., TSE validation is usually performed with a committor test 4 , 7 , 9 , 33 ., In principle , the committor test launches unbiased MD simulations from structures chosen randomly from the barrier region , and tests whether the resultant trajectories reach the product state with probability 1/2 ., Here , since limited computational resources precluded execution of a full committor test , 40 unbiased MD simulations of 10 ns were initiated from each of , 33 or 34 , a total of 120 simulations or 1 . 2 , and the distributions of the final structures after 10 ns were monitored 9 ., Figure 5A shows the binned distributions of the final structures assigned by index of the nearest MFEP image ( the blue bars ) ., When the simulations were initiated from the image at ( ) , the distribution biases to the open form-side ( the closed form-side ) relative to the initial structures ., On the other hand , when starting from the image at , the distribution is roughly symmetric around the initial structures ., This result suggests that the TSE is located at ., In other words , it was validated that the TSE was successfully captured in the MFEP , and at the same time , the collective variables were good enough to describe the transition ., A close inspection of the structures around the PMF barrier supported our insufficient committer test and revealed the mechanism of the ligand-induced domain closure ., Figure 5B shows the hydrogen bond ( H-bond ) patterns between the ligand and the protein observed in the average structures at ( before the TSE ) and ( after the TSE ) ., At , Thr31:OG1 ( AMPbd ) forms an H-bond with N7 of AMP adenine , and Gly85:O ( CORE ) forms one with adenine N6 ., These two H-bonds mediate the hinge bending of the AMPbd-CORE domains ., In addition , the H-bond between Gly85:O and adenine N6 helps the enzyme to distinguish between AMP and GMP; GMP lacks an NH2 group in the corresponding position of AMP 34 ., This means that the specificity of AMP-binding operates at an early stage of the ligand binding process ., At , the AMPbd-CORE distance becomes smaller than that at , which allows the formation of 3 additional H-bonds with the ligand: Gln92:OE1 ( CORE ) and adenine N6 , Lys57:O ( AMPbd ) and the ribose O2 , and Arg88:NH1 ( CORE ) and O1 of AMP ., The resulting rapid enthalpy decrease stabilizes the closed conformation ., Gln92:OE1 is also important in establishing AMP specificity; GMP lacks the counterpart atom , adenine N6 ., The strictly conserved Arg88 residue is known to be crucial for positioning AMP so as to suitably receive a phosphate group from ATP 35 ., With regard to the AMPbd closure , our result suggests that Arg88 ( CORE ) , in conjunction with Lys57 ( AMPbd ) , works to block adenine release from the exit channel and to further compact the AMPbd-CORE domains ., A remaining question is how closure of the LID domain follows that of the AMPbd domain ., Unlike the AMP-binding pocket , the ATP-binding sites , including the P-loop , are surrounded by charged residues , which attract interfacial water molecules ., Upon LID closure , most of these water molecules will be dehydrated from the enzyme , but some may remain occluded ., To characterize the behaviors of these water molecules , the 3D distribution function of their oxygen and hydrogen constituents were calculated along the MFEP using the MBAR method ( see Materials and Methods ) ., Figures 6A , 6B , and 6C display the isosurface representations of the 3D distribution functions around the P-loop at , 41 , and 42 , respectively ., The surfaces show the areas in which the atoms are distributed four times as probably as in the bulk phase ., At , the ATP phosphates are not yet bound to the P-loop because an occluded water molecule ( encircled ) is wedged between the phosphate and the P-loop , inhibiting binding of ATP and and bending of the side-chain of “invariant lysine” ( Lys13 ) , a residue that plays a critical role in orienting the phosphates to the proper catalytic position 36 ., This occluded water molecule may correspond to that found in the crystal structure of apo-AK ( PDBid: 4ake ) ( Figure 6D , encircled ) ., Figures 6B and 6C clearly demonstrate that , upon removal of this water molecule , the ATP phosphates begin binding to the P-loop ., These observations were confirmed by plots of the PMF surface mapped onto a space defined by the LID-CORE distance versus the index of image ( Figure 6F ) , which shows that the PMF decreases discontinuously upon dehydration followed by LID domain closure ., Interestingly , compared with the crystal structure ( PDBid: 1ake ) ( Figure 6E ) , the position of the ATP moiety is shifted to the AMP side by one monophosphate unit ., This may be a consequence of early binding of the AMP moiety ., At a later stage ( around ) , this mismatch was corrected to form the same binding mode as observed in the crystal structure ., This reformation of the binding mode may be induced by the tight binding of ATP adenine to the LID-CORE domains , and will not occur in the real enzymatic system containing ATP and AMP instead of the bisubstrate analog inhibitor , Ap5A ., In this study , we have applied the on-the-fly string method 5 and the MBAR method 11 to the conformational change of an enzyme , adenylate kinase , and successfully obtained the MFEP ( Figures 2A and 2B ) ., The MFEP yielded a coarse-grained description of the conformational transitions in the domain motion space ( Figures 2C and 2D ) ., At the same time , the atomistic-level characterization of the physical events along the MFEP provided a structural basis for the ligand-binding and the domain motions ( Figures 3–6 ) ., This kind of multiscale approach used here is expected to be useful generally for complex biomolecules since the full space sampling can be avoided in an efficient manner ., We have shown that in the TSE of holo-AK , the conformational transition is coupled to highly specific binding of the AMP moiety ., Our results have been validated by unbiased MD simulations ., The mechanism of the AMPbd domain closure is consistent with that proposed by the induced-fit model ( Figure S7A ) , and follows a process similar to that of protein kinase A , previously investigated by a coarse-grained model 32:, ( i ) the insertion of the ligand into the binding cleft initially compacts the system;, ( ii ) additional contacts between the ligand and non-hinge region further compact the system ., The closure of the LID domain is more complicated ( Figure S7B ) ., It was shown that apo-AK can exist in a partially closed state , stabilized by the “cracking” of the LID-CORE hinge and the P-loop , even with no ligand present ., The cracking of the hinge region enables rearrangement of molecular interactions for ATP-binding , which induces a smooth bending of the hinge ., Along with the LID closure , ATP is conveyed into the P-loop , with removal of an occluded water molecule ., The closure of the LID domain follows the “population-shift followed by induced-fit” scenario discussed in Ref ., 37 , in which a transient local minimum is shifted toward the closed conformation upon ligand binding ., This two-step process of the LID domain closure is similar to the two-step mechanism reported in recent simulation studies of the Lysine- , Arginine- , Ornithine-binding ( LAO ) protein 38 and the maltose binding protein 39 ., In holo-AK , AMPbd domain closure occurs early ( at ) , while the LID domain closes at later stages ( ) ., An interesting question is whether an alternative pathway is possible in the presence of the real ligands ( ATP and AMP ) instead of Ap5A ., Ap5A artificially restrains the distance between the ATP and AMP moieties ., During the process with real ligands , the dynamics of the LID and AMPbd domains is expected to be less correlated ., Nevertheless , for full closing of the LID domain , we conjecture that the AMPbd domain should be closed first , enabling the interactions on the LID-AMPbd interface to drive the dehydration around the P-loop ., This suggests that full recognition of ATP by the LID-CORE domains occurs at a later stage of the conformational transition ., This conjecture may be related to the lower specificity of E . coli AK for ATP compared with AMP 34 ., Nonspecific AMP-binding to the LID domain has previously been suggested to explain the observed AMP-mediated inhibition of E . coli AK at high AMP concentrations 40 ., A missing ingredient in the present study is the quantitative decomposition of the free energy in each event , such as the ligand binding and the interactions on the LID-AMPbd interface ., For enhanced understanding of the conformational change , our methods could be complemented by the alchemical approach 41 ., Varying the chemical compositions of the system during the conformational change would enable us to elucidate the effects of ligand binding , cracking , and dehydration in a more direct manner ., We prepared three systems from the following initial structures:, ( i ) “apo-open system” , X-ray crystal structure of the open-form without ligand ( PDBid: 4ake 28 ) ,, ( ii ) “holo-closed system” , crystal structure of closed-form with Ap5A ( PDBid: 1ake 29 ) ,, ( iii ) “apo-closed system” , structure created by removing Ap5A from the holo-closed system ., The protonation states of the titratable groups at pH 7 were assigned by PROPKA 42 , implemented in the PDB2PQR program package 43 , 44 ., The apo-open and apo-closed systems yielded identical assignments , which were used also for the holo-closed system ., These systems were solvated in a periodic boundary box of water molecules using the LEaP module of the AMBER Tools ( version 1 . 4 ) 45 ., A padding distance of 12 Å from the protein surface was used for the apo-open system ., For the apo-closed and holo-closed systems , a longer padding distance of 20 Å was used to avoid interactions with periodic images during the closed-to-open transition ., Two Na+ ions were added to neutralize the closed-apo and open-apo systems , while seven Na+ ions were required to neutralize the closed-holo system ., The systems were equilibrated under the NVT condition at 300 K by the following procedure: First , the positions of solvent molecules and hydrogen atoms of the protein ( and Ap5A ) were relaxed by 1 , 000 step minimization with restraint of non-hydrogen atoms ., Under the same restraints , the system was gradually heated up to 300 K over 200 ps , followed by 200 ps MD simulation under the NVT condition at 300 K while gradually decreasing the restraint forces to zero , but keeping the restraints on atoms needed in the string method ., The system was further equilibrated by 200 ps MD simulation under the NPT condition ( 1 atm and 300 K ) , adjusting the density of the water environment to an appropriate level ., The ensemble was finally switched back to NVT , and subjected to additional 200 ps simulation at 300 K , maintaining the restraints ., The equilibration process was conducted using the Sander module of Amber 10 45 , with the AMBER FF03 force field 46 for the protein , and TIP3P for water molecules 47 ., The parameters for Ap5A were generated by the Antechamber module of AMBER Tools ( version 1 . 4 ) 45 using the AM1-BCC charge model and the general AMBER Force Field ( GAFF ) 48 ., Covalent bonds involving hydrogen atoms were constrained by the SHAKE algorithm 49 with 2 fs integration time step ., Long-range electrostatic interactions were evaluated by the particle mesh Ewald method 50 with a real-space cutoff of 8 Å ., The Langevin thermostat ( collision frequency 1 ps−1 ) was used for the temperature control ., The production runs , including the targeted MD , the on-the-fly string method , the umbrella sampling , and the committor test , were performed with our class library code for multicopy and multiscale MD simulations ( which will soon be available ) T . Terada et al . , unpublished , using the same parameter set described above ( unless otherwise noted ) ., Protein structures and the isosurfaces of solvent density were drawn with PyMOL ( Version 1 . 3 , Schrödinger , LLC ) ., The calculations were performed using the RIKEN Integrated Cluster of Clusters ( RICC ) facility ., It has been shown that normal modes or principal modes provide a suitable basis set for representing domain motions of proteins 2 , 10 ., In particular , it has been argued that the conformational change in AK can be captured by a set of principal modes of apo-AK 31 , 51 ., In this study , we have defined the collective variables for the on-the-fly string method using the principal components of apo-AK ., The PCA was carried out in the following manner: After the equilibration process , 3 ns MD simulations were executed at 300 K without restraint for both apo-open and apo-closed systems ., The obtained MD snapshots from both systems were combined in a single PCA 52 , removing the external contributions by iteratively superimposing them onto the average coordinates 53 , 54 ., The PCA was then conducted for the Cartesian coordinates of the atoms ., It was found that the first principal mode representing the largest-amplitude merely represents the difference between the open and closed conformations ., The fluctuations in the two structures were expressed in the principal modes of smaller amplitudes ., The cumulative contributions of these modes ( ignoring the first ) are shown in Figure S8 ., As expected , the principal modes represent the collective motions of the LID and AMPbd domains ( Figure S9 ) ., The first 20 principal components ( 82% cumulative contribution , ignoring that of the first ) were adopted as the collective variable of the string method ., These components were sufficient to describe the motions of three domains in AK for which at least degrees of freedom are required in the rigid-body approximation ., The additional eight degrees of freedom were included as a buffer for possible errors in the estimation of the principal modes ., The sum of the canonical correlation coefficients between the two sets of the 20 principal components , one calculated using the samples of the first half ( 0–1 . 5 ns ) snapshots and the other using the last half ( 1 . 5–3 ns ) snapshots , was 11 . 8 ( ∼12 ) , suggesting that the subspace of the domain motions was converg
Large-scale conformational changes in proteins involve barrier-crossing transitions on the complex free energy surfaces of high-dimensional space ., Such rare events cannot be efficiently captured by conventional molecular dynamics simulations ., Here we show that , by combining the on-the-fly string method and the multi-state Bennett acceptance ratio ( MBAR ) method , the free energy profile of a conformational transition pathway in Escherichia coli adenylate kinase can be characterized in a high-dimensional space ., The minimum free energy paths of the conformational transitions in adenylate kinase were explored by the on-the-fly string method in 20-dimensional space spanned by the 20 largest-amplitude principal modes , and the free energy and various kinds of average physical quantities along the pathways were successfully evaluated by the MBAR method ., The influence of ligand binding on the pathways was characterized in terms of rigid-body motions of the lid-shaped ATP-binding domain ( LID ) and the AMP-binding ( AMPbd ) domains ., It was found that the LID domain was able to partially close without the ligand , while the closure of the AMPbd domain required the ligand binding ., The transition state ensemble of the ligand bound form was identified as those structures characterized by highly specific binding of the ligand to the AMPbd domain , and was validated by unrestrained MD simulations ., It was also found that complete closure of the LID domain required the dehydration of solvents around the P-loop ., These findings suggest that the interplay of the two different types of domain motion is an essential feature in the conformational transition of the enzyme .
Conformational transitions of proteins have been postulated to play a central role in various protein functions such as catalysis , allosteric regulation , and signal transduction ., Among these , the relation between enzymatic catalysis and dynamics has been particularly well-studied ., The target molecule in this study , adenylate kinase from Escherichia coli , exists in an open state which allows binding of its substrates ( ATP and AMP ) , and a closed state in which catalytic reaction occurs ., In this molecular simulation study , we have elucidated the atomic details of the conformational transition between the open and the closed states ., A combined use of the path search method and the free energy calculation method enabled the transition pathways to be traced in atomic detail on micro- to millisecond time scales ., Our simulations revealed that two ligand molecules , AMP and ATP , play a distinctive role in the transition scenario ., The specific binding of AMP into the hinge region occurs first and creates a bottleneck in the transition ., ATP-binding , which requires the dehydration of an occluded water molecule , is completed at a later stage of the transition .
computational chemistry, molecular dynamics, biophysic al simulations, chemistry, biology, computational biology