Modeling Protein Conformational Transition Pathways Using

Subscriber access provided by UB + Fachbibliothek Chemie | (FU-Bibliothekssystem)

Article

Modeling protein conformational transition pathways using collective motions and the LASSO method Thomas W. Hayes, and Iain H. Moal J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b01110 • Publication Date (Web): 23 Feb 2017 Downloaded from http://pubs.acs.org on March 1, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Modeling protein conformational transition pathways using collective motions and the LASSO method Thomas W. Hayes and Iain H. Moal∗ European Molecular Biology Laboratory, European Bioinformatics Institute (EMBL-EBI) E-mail: [email protected]

Abstract Many proteins can adopt multiple distinct conformational states which often play different functional roles. Previous studies have shown that the underlying global dynamics through which these states are accessed are, at least in part, encoded by the protein’s topology. In this work we present a method for generating transition pathways between states by perturbing the protein towards a target conformational state along thermally accessible collective motions calculated from the starting conformation. Specifically, the least absolute shrinkage and selection operator (LASSO) is used to identify the most parsimonious route along soft modes calculated using the anisotropic network model. In a survey of 436 conformational changes following protein-protein interaction, we show that such a path exists for most cases and that selected paths are low in energy. We discuss the implications for the atomic modelling of protein recognition, and provide soft energy and parameter bounds which can be employed to efficiently constrain the sampling of such pathways. ∗

To whom correspondence should be addressed

1

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

Global protein dynamics allow a protein to switch functional activity by adopting multiple conformations. While instantaneous thermal fluctuations and longer-scale deformations are both governed by the balance of interatomic forces, low-frequency essential dynamics arise from the rearrangement of the least densely compacted regions. Using this as a basis, elastic network models (ENMs) treat proteins as a system of balls and springs and have proven cheap and convenient heuristics in uncovering flexibility profiles 1–4 and revealing collective degrees of freedom along which alternate states are usually found, such as those co-opted by evolution for catalysis, allosteric regulation and molecular recognition. 5–14 As these collective motions allow the description of motions of many atoms using a small number of low-frequency normal modes, they have seen a wide range of applications. These include modeling biomolecular assemblies by fitting structures to experimental data, such as from electron microscopy 15–17 or x-ray crystal diffraction, scattering or fiber diffraction. 18–20 A requirement for these applications is that a conformation near to the final state is accessible through a linear combination of low-frequency modes calculated from the initial starting structure, but not necessarily the existence of a low energy pathway to the final state. Other applications of ENMs involve sampling along normal coordinates, such as steering molecular dynamics or Monte Carlo simulations, 21–25 with applications in enzyme engineering, 26 calculating ligand binding affinities 27 and exit pathways, 28 and the docking of proteinligand, 29–31 protein-DNA 29,32 and protein-protein 30,33–37 complexes. These methods have the additional requirement of a low-energy transition pathway connecting the start state to the final state and, as only the start state is typically known, can only access configurations derived using modes derived from one configuration. This additional restriction means that the applicability of these methods cannot be assessed by techniques which identify putative transition pathways using two-state models, in which coarse-grained potential energy surfaces of both end states are combined. 38–48 Nevertheless, there have been some studies of transitional paths in which only intra-molecular potential energy surface calculated from 2


Page 2 of 30

Page 3 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


either the start structure, or a structure derived from it, are navigated. Along these lines, Kantarcı-Çarşıbaşı et al. 49 iteratively modeled the transition pathways of hemoglobin and an adenylate kinase by moving towards the end state along a selected low frequency mode of the start state, performing energy minimization, updating the normal modes, and repeating until a configuration near the end state was found. Kürkçüoğlu and Doruker 50 iteratively performed an unbiased search by deforming along the three lowest frequency modes, clustering, minimizing and recalculating modes for each central cluster member, and repeating the process for each resultant structure in each successive iteration. This could generate conformations close to the bound state from which snapshots of a transition path can be obtained by tracing back the generations that lead to it. This work was expanded upon in Kürkçüoğlu et al. 51 . In this paper we present a method of identifying plausible transition pathways by perturbing the start structure towards the end state. We employ the LASSO method, which traces out the path by relaxing a constraint on the extent of perturbation along multiple low-frequency normal modes. As these modes are derived only from the start state, this yields information which can be used in modelling tasks for which only a single-state potential energy landscape is available. As opposed to producing snapshots, this approach creates a continuous path, and it is not constrained to moving along a single mode at a time. There is also no requirement of costly mode recalculation or energy minimisations. In the following sections we give a overview of the method and the concepts upon which it is based, followed by an illustration of the procedure using carboxypeptidase B and the FH2 domain of Bni1p. We then show 14 proteins for which calculated transition pathways pass via intermediate conformations which have been crystallographically isolated. Permutation testing is used to illustrate how the underlying assumptions in the model avoid pathways that proceed via unnatural regions of conformational space. We then survey 436 conformational changes upon molecular recognition events, showing that in most cases a path can be found to a structure within 2Å RMSD of the final state, and that this path is characterized by low ENM energy.

3



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Finally, we discuss the implications for molecular modeling tasks in which sampling along normal coordinates is undertaken, and identify the types of motion for which pathways can and cannot be constructed.

2 2.1

Results and discussion Overview of method

Our method is based on the notion that conformational changes correspond to motions that proteins are intrinsically predisposed to undergo, in which the predisposed motions are related to the thermally accessible low-frequency normal modes calculated using the anisotropic network model (ANM). 6,7 In the ANM scheme, a normal mode of vibration consists of harmonic linear oscillations of point masses around their equilibrium positions. The oscillations of all points are in phase, and as such reach their greatest extent simultaneously and return to their equilibrium position simultaneously as they vibrate back and forth. The ANM method generates 3N − 6 such normal modes, where N is the number of positions used. Each normal mode has a characteristic frequency, with the lowest frequency modes requiring the least amount of energy in order to be excited. Thus it is the low frequency modes that are the most thermally accessible. Each normal mode can be represented as a vector of length 3N where the first three elements correspond to the displacement of the first point mass along the x, y and z coordinates, with subsequent blocks of three corresponding to Cartesian components of the subsequent points. In this work we use these normal coordinate vectors as a 50-dimensional vector space, comprising of the 50 lowest frequency modes. Each position within this vector space corresponds to a molecular structure, the original structure plus a perturbation specified by a linear combination of normal coordinates; each vector is multiplied by a scalar coefficient, its corresponding coordinate in the vector space, and the resultant vectors are added together and applied as the perturbation to the initial structure. The origin of this vector space corresponds to a structure in which all the atoms are at their 4


Page 4 of 30

Page 5 of 30


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

initial equilibrium coordinates, and all other positions correspond to a perturbation of the initial structure along the direction of one or more low-frequency mode. Thus, the region surrounding the origin maps directly onto the thermally accessible conformations that arise from the application of the ANM on the starting structure. In our proposed method we calculate a path within this vector space from the starting structure such that the corresponding conformation moves towards an alternative molecular arrangement. For any given case, it is highly improbable that a position within the vector space exists in which the coordinates of the perturbed structure exactly match those of the alternative conformation, so our strategy is to calculate a pathway to the position that minimises the RMSD between that position and the alternative state. Thus the method generates a path to the closest possible position to the final state attainable using the linear combination of modes derived from the start state. A single parameter specifies the path from the origin to this structure. Each position along the path is characterised by β, a vector containing the coordinates within the vector space, and which are the coefficients in the linear combination. The path is continuous and yields conformations that are steered towards the target by relaxing an ℓ1 constraint on the magnitude of β and performing a least squares fit of the normal coordinates to the displacements in atomic positions observed in the conformational change. This is done using the least absolute shrinkage and selection operator (LASSO) method. 52 The LASSO uses the gradient of the sum of RMSD and the ℓ1 penalty within the vector space to select the modes down which to travel, identifying the next step in the transition pathway by finding the global minimum using gradient descent, and progressively relaxes the penalty in order to define the whole transition pathway. The gradient descent ensures that out of all the possible routes from the origin to the final position, including those for which the RMSD is monotonically decreasing, the structure is the cloest to the final state that is possible within the constraints of the ℓ1 penalty. Thus it is the shortest and simplest path to the end state, and hence the most parsimonious. For more details, see the Methods section.

6


Page 6 of 30

Page 7 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.2

Illustration of method

To demonstrate the method, Figure 1 shows calculated transition pathways for two cases. The first is the change in the formin homology 2 (FH2) domain of Bni1p 53 upon binding to two actin molecules, mimicking binding to the barbed end of filamentous actin. 54 The C-terminal region of the FH2 contains a well folded and compact region consisting of the knob subdomain, a coiled coil, the post subdomain and a C-terminal helix. The N-terminal region is less structured, with a loop subdomain and an α-helical linker. The calculated transition pathway involves opening of the shallow cavity by elongation of the linker and an increase in coiled-coil/knob/linker and knob/linker/loop angles (Table 1), a flip between cis-staggered configurations, and a rotation of the loop region. This constitutes a large conformational change, with RMSD of 10.3Å between the bound and unbound crystal structures, which drops to 2.7Å between the end of the calculated transition path and the bound structure. Over three quarters of the reduction in RMSD arises from use of only the four lowest frequency modes, showing that the method is preferentially selecting the most thermally accessible motions. Initially, the first two modes are activated, with a 25% reduction in RMSD and the knob/linker/loop angle increasing from 117◦ to 131◦ . This enables the molecule to cross the eclipsed configuration, indicated by the change of sign of the coiled-coil/knob/linker/loop dihedral, which accelerates as the coiled-coil/knob/linker angle increases concomitantly with the moving apart of the knob and linker. The linker/loop distance also increases steadily, until the final stages, when the arc closes back into the new staggered conformation by decreasing the linker/loop distance and coiled-coil/knob/linker angle. This mechanism, which is analogous to the rotation around an sp3 -sp3 chemical bond, would not be possible using methods not based on the underlying global dynamics such as morphing, in which the start state and end state are interpolated. The second example in Figure 1 corresponds to the tick carboxypeptidase inhibitor 55 upon binding to carboxypeptidase B. 56 This motion involves the rearrangement of two globular domains about a flexible linker, as well as rearrangements in the C-terminal domain. The 7



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Table 1: Descriptive distance (d, in Å), angle (θ, in ◦ ) and dihedral (φ, in ◦ ) parameters for the transition path of FH2 domain of Bni1p upon binding to actin, calculated from subdomain geometric centres. φcc-kn-li-lo θcc-kn-li θkn-li-lo dkn-li dli-lo unbound 20.6 132 117.2 22.2 29.7 a 18.6 132.5 130.5 22.5 31 25% 50%a 11 133.6 138.1 23.1 32.2 75%a -1.2 136.5 141.8 24.4 33.7 a -16.6 140 143.1 25.6 34.8 90% -41.1 135.1 151.9 26 33.9 final -50.4 133.6 149.9 26.3 34.3 bound C-terminal knob subdomain: kn, residues 1418-1521. Coiled coil: cc, 1522-1562 and 1640-1719. N-terminal loop: lo, 1350-1400. α-helical linker: li, 1401-1417. a Snapshots are reported as a percentage of the decrease in RMSD. calculated transition path initially selects the first mode, corresponding to torsional motion about the linker, before simultaneously selecting the 3rd and 6th modes, respectively a pincer motion between the two domains, and the concerted rotation of the domains about their hinge points with the linker. This correctly positions the domains relative to one another, after which we see rearrangement of the loop connecting the first β-strand and α-helix of the second C-terminal domain.

2.3

Cases for which pathways proceed via isolated intermediates

The method is capable of identifying potential pathways between a start point and a position close to the end point. These pathways are plausible, as they are based on perturbations of the start structure along thermally accessible collective modes, and they are parsimonious, as they correspond to the smallest possible perturbations along these modes required to approach the end point. However, due to the difficulty in observing short-lived intermediate conformations experimentally, it is not possible to validate them by comparison to direct observation of the transition. Nevertheless, occasionally intermediate conformations are isolated crystallographically, being stabilized by crystal contacts, redox state or ligands such as other proteins, small molecules, cofactors, antibodies or ions. While it is not generally

8


Page 8 of 30

Page 9 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


possible to know for certain whether any particular crystal configuration is visited during a transition, it is reasonable to expect alternative configurations to appear along some transition paths. By manually searching a recent survey of alternative conformational states, 57 we identified 14 cases for which an intermediate along the calculated transition pathway is close to an alternative conformation that is distinct from both the start and end state. Each such path can be described by a number of RMSD values which measure the proximity between the relevant structures: the start state (S), end state (E), intermediate state (I), the end of the transition pathway (pfinal ) and the position along the calculated pathway that is closest to I in terms of RMSD (pi ). The S/E RMSD shows the extent of the conformational change, the pfinal shows how close the end point of the transition path is to the final state, and the I/S and I/E RMSDs show how close the experimental intermediate is to the start and end states respectively. The I/pi , pi /S and pi /E RMSDs show how close the position closest to the experimental intermediate is to the intermediate itself, to the start state and to the end state respectively. These values are reported in Table 2, with a visual summary of the terms depicted in the top left panel of Figure 2 alongside the corresponding plots for the 14 examples. Further details of the start state, end state, intermediates and transition pathways are reported in the Supplementary Text.

2.4

Permutation testing

We wished to further investigate the reliability of the calculated pathways by exploring whether the method could produce spurious transitions by attempting to model pathways to contrived end points. If the method can generate pathways to false final conformations then it would be inherently unreliable, as it would be difficult to discern whether a calculated pathway to a genuine end state is produced due the validity of the method, or whether that pathway is also spurious. On the other hand, if the method cannot produce pathways to contrived end points, then it would give us greater confidence in the reliability of the pathways when they are produced. To test this, we employ permutation tests in which we sample 9




Page 10 of 30

Page 11 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 2: Summaries of intermediates in selected pathways, indicating the start point (S), end point (E), final point along transition path (pfinal ), the experimental intermediate (I) and the calculated intermediate closest to it (pi ). See Figure 2 for a visual interpretation of the parameters, and Supplementary Text for more details. Protein Name Fab 17b Cox6c Integrin β3 NikZ DntR " β-alanine synthase Chitinase A " Ski-like protein " " " " NscB1 PKM2 " Mab 4C1 " DAPDH " " ptxD POLB " " "

S 1RZ8_A 1OCC_V 2VDQ_B 4OET_A 2Y7K_A " 2V8D_A 3B8S_A " 3EQ5_C " " " " 3I58_B 3GQY_D " 3RVV_D " 3WBB_A " " 4EBF_A 4DO9_A " " "

E 2NXY_C 2EIK_V 3NIG_B 4OEV_B 2Y84_D " 2V8H_C 3B8S_B " 3EQ5_F " " " " 3I53_B 4JPG_C " 3RVU_D " 3WBF_B " " 4EBF_C 7ICE_A " " "

S/Ea 2.9 1.0 10.3 2.5 1.7 " 4.2 1.6 " 1.9 " " " " 3.3 2.2 " 3.7 " 2.4 " " 2.2 5.1 " " "

pfinal a 0.63 0.2 3.2 0.5 0.5 " 0.4 0.2 " 0.88 " " " " 1.0 0.4 " 0.5 " 0.3 " " 0.6 1.2 " " " a RMSD (Å)

11

I 1RZ8_C 3AG1_I 3ZE2_D 4OEU_B 2Y7K_B 2Y7K_C 2V8D_B 3ARZ_A 3AAR_A 3EQ5_J 3EQ5_H 3EQ5_L 3EQ5_I 3EQ5_K 3I53_A 4B2D_B 4G1N_A 3RVT_D 3RVW_D 3WB9_A 3WB9_B 3WB9_C 4EBF_F 4PHA_A 4NM2_A 2I9G_A 4M47_A

I/Sa 1.5 0.74 5.5 2.0 0.8 1.3 3.3 1.4 1.4 1.0 1.0 1.3 1.5 1.3 1.0 0.9 1.5 0.9 1.1 0.7 0.6 0.6 1.6 1.7 2.8 3.0 2.1


I/Ea 1.5 0.46 5.2 1.0 1.4 2.1 1.0 1.0 1.0 0.7 0.9 0.3 0.3 0.3 2.7 1.7 1.1 3.0 2.7 1.8 1.9 1.9 0.8 4.1 3.2 3.4 3.3

I/pi a 0.61 0.43 2.7 0.80 0.56 0.65 0.48 0.74 0.75 0.64 0.73 0.52 0.58 0.66 0.60 0.59 0.77 0.49 0.58 0.50 0.46 0.50 0.57 1.14 2.17 2.39 1.53

pi /Sa 1.3 0.72 6.4 1.6 0.5 1.2 3.5 1.0 1.0 1.0 0.9 1.5 1.6 1.4 0.8 0.7 1.3 0.8 " 1.9 " " 1.5 1.4 3.1 3.0 2.6

pi /Ea 1.8 0.34 5.3 1.2 1.3 0.8 0.7 0.8 0.8 1.3 1.3 0.9 0.9 1.0 2.9 1.6 1.0 3.0 " 0.6 " " 0.9 4.2 2.7 2.9 3.3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

200 random permutations of d, the vector of atomic displacements between start and end points, for each of the 96 flexible proteins (RMSD > 2Å) in the protein-protein interaction benchmark 5.0. 58 Each permutation represents a conformational change to an end point with the same magnitude of change as the true conformational change, and the same elements as the true displacement vector, but with those elements randomly shuffled. Consequently, any transition pathway to such a conformation would be necessarily nonsensical. In all 200 permutations of all 96 cases the method fails to find a pathway to a configuration near the contrived end state (Supplementary Fig. 1), greater progress can be made towards the true end state than for any of the spurious ones, and all attempts to find pathways to spurious end states are accompanied by a large and immediate increase in energy. By contrast, most pathways towards true end states are of low energy, increasing only towards the end as indicated by L-shaped energy/RMSD plots (Supplementary Figure 1). The inability of the method to generate pathways to unnatural states illustrates how the underlying model stops the algorithm from forcing a trajectory to the end state when it is prevented from doing so by the topology of the protein and the assumption that the change is driven by the thermally accessible low-frequency motions encoded in the normal modes.

2.5

Survey of pathways taken upon protein-protein interactions

To further investigate the presented method, we applied it to a benchmark of 436 conformational changes occurring upon protein-protein interactions. 58 Figure 3a shows the initial RMSD to the bound state versus the final RMSD expressed as a percentage of the initial RMSD. For cases of intermediate flexibility (initial RMSD values between 2.0Å and 3.5Å), a path to a position less than 2Å from the final state can be identified for 93% of the proteins, and to a position within 1.5Å of the final state for 72%. For the high flexibility cases (> 3.5Å initial RMSD), the method can identify a putative pathway to within 2.0Å for 50% of the proteins. For the proteins with initial RMSD > 2Å, the RMSD to their respective final states along their transition pathways are shown in Figure 3b, with progress along the 12


Page 12 of 30

Page 13 of 30




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

pathway denoted by the normalized sum of absolute coefficients in the linear combination of modes, |β|/N . The initial RMSD values span a range of 2.0 to 10.4Å, with highest density in the 2 to 4Å range. After being perturbed along low-frequency modes towards the final state with an allocated |β|/N of 0.1, RMSD density has lowered to points around 1 to 3Å of the final state. By the time a |β|/N of 0.4 has been allocated, most proteins have reached the end of their transition pathway, with a high density of pathways having converged on a point around 1.2Å RMSD from their bound states. Figure 3c shows how far the transition path can be extended toward the final state as increasing numbers of modes are included, starting from the lowest frequency mode on its own up to a combination of the 50 lowest frequency modes, for the same set of proteins (initial RMSD > 2Å). For many proteins, particularly the most flexible (initial RMSD > 3Å), most of the transition can be modeled using the lowest frequency mode alone. For most cases a position near to the final state can be arrived at using one of, or a combination of, the five lowest frequency modes. Although there are exceptions, few of the paths benefit from the inclusion of the 11th to 50th mode. The drop in RMSD along the pathway, and the corresponding changes in the magnitude of the coefficients in the linear combination (|β|/N ) and the ENM energy are illustrated in Figure 3d and 3e. A gradual increase in |β|/N is typically observed as the RMSD drops. For the energy, however, it is more common to see low energy across most of the transition path, with larger increases only becoming apparent near the terminal point, at which point the higher frequency modes are used to make the final adjustments to the structure (see also Supplementary Figure 1). To give an indication of the sorts of structural perturbations which can be achieved with a given quantity of energy, the inset of Figure 3e shows the extent of conformational changes when the SRP GTPase is thermally activated to 300 or 1000 units of energy by randomly partitioning this energy across random directions in the 50 lowest frequency modes. This shows that, within the ENM model, the energy required to move to a known alternative conformations is not so large that it would result in large

14


Page 14 of 30

Page 15 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and unrealistic perturbations were it to be used to move to random positions of equal energy within the coefficient space. Subsequently the ENM energy, and indeed the number of modes used and |β|/N , can be used to provide soft bounds which, were they to be used to limit conformational sampling in methods that search along normal coordinates such as flexible docking, would ensure that the method could traverse the transition pathway to positions near alternative states while at the same time avoiding highly distorted conformations. For instance, using constraints of |β|/N < 1.0 and E < 1500, and searching the 40 lowest frequency modes, is enough to ensure that sufficient conformational space is accessible in order to sample almost the whole transition path for all of the high flexibility proteins for which a pathway can be modelled. Nevertheless, halving the dimensionality of search space to the lowest 20 modes, and using constraints of |β|/N < 0.7 and E < 500 permits enough conformational freedom to reach the end of the transition path for most cases where such a path is possible, and provides a very large reduction in the difficulty of the search problem. To better understand the factors influencing whether or not a pathway could be identified, individual cases that were particularly successful or unsuccessful were inspected for the types of conformational change involved (Figure 3a upper and lower panels). Unsuccessful cases included those for which the conformational changes are known to be poorly modeled by the anisotropic network model: 59 movement dominated by highly non-linear motion, such as 180◦ rotation of a rigid domain, and non-collective motions such as partial refolding. The method also fails to find transition pathways in some cases where the subunits are weakly coupled, when the conformational change is within a subunit but the low frequency modes are dominated by intersubunit motion. On the other hand, many of the large-scale motions for which pathways were readily calculable involved hinge bending motions, sometimes accompanied by twisting motion. Other tractable transitions include more complex collective motions, such as combinations of twisting, bending and clasping.

15



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

3

Conclusions

We have presented a method for identifying putative transition pathways between protein conformations using the normal modes of vibration of only one conformation. In doing so we have shown that, for flexible proteins, the accessible conformations encoded in those modes are often sufficient to derive a plausible path to a position near to the desired end state. We have identified 14 cases for which the transition passes close to known intermediate configurations, and used perturbation testing to illustrate how the method will only produce pathways when the conformational change is consistent with the thermally accessible motions inherent in the model, rather than force a spurious pathway between conformers. This result has implications for simulation and optimization tasks involving backbone flexibility such as protein-protein docking. The pre-generation of structural ensembles, either by molecular dynamics, normal modes or other methods, typically fails to produce configurations that are appreciably closer to the final state than the one they were generated from. 60–63 By generating a series of low energy intermediate configurations, the results presented here give credence to methods that instead dynamically search for such as conformation within a search space defined by pre-calculated normal modes. Further, our results provide clear guidelines on the number of modes required, the maximum magnitude of the coefficients required to sample the transition pathway, and limits on the ENM energy beyond which no configurations near to the end state are likely to be found.

4 4.1

Methods Data Sets

To identify proteins with crystallographically isolated intermediate conformations, we started with the alternative conformation data set of Clarke et al. 57 , containing clusters of high resolution (< 2.8Å) structures with 100% sequence identity, covering 1106 proteins. Most of

16


Page 16 of 30

Page 17 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the proteins in the set have only two distinct conformations, often with multiple PDB entries of identical configurations. We searched for pairs of conformers, the most representative of their cluster (cluster member with shortest Euclidean distance to the cluster average), for which a third conformer was available and which was sufficiently distinct from both the representative conformers; either part of a third cluster or on the periphery of a diffuse cluster. In total 35 proteins had an alternative conformer further than 2Å Cα RMSD from both end states, and 184 with RMSD > 1Å. For the conformational changes, the start and end states were arbitrarily attributed to the representative cluster members, and we performed a manual search for cases in which an alternative configuration could be found near the transition path. The second data set, of conformational changes upon binding, was taken from the nonredundant manually curated protein-protein docking benchmark 5.0, 58 which includes 230 structurally superimposed bound/unbound pairs of high-resolution (< 3.25Å) structures. Several proteins were removed, due to a large mismatch in the number of residues in the start and end states (1K4C_l, 1GCQ_r, 2VIS_l, 1H1V_l, 1FC2_l, 1QA9_l, 1ZLI_r), large unfolded protrusions (4GAM_l, 1IRA_r), 12 cases where the unbound antibodies are not available, and the three duplicate cases arising from multiple binding partners, leaving 436 conformational changes upon binding. In these cases, the unbound conformation was considered the start state state and vice versa.

4.2

Normal Mode Calculation

All-atom normal modes were calculated for the start states using the ElNeMo package. 64 Where small molecules are annotated in the benchmark 5.0, these were also included. A 5Å distance cutoff was used such that atom pairs below this threshold were connected by Hookean springs with forces constants of 1 and equilibrium distances equal to that of the starting structures. Hydrogen atoms were discarded, and non-hydrogen atoms treated as points of mass 1. The Hessian matrix for the resultant system was constructed, and diago17



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 30

nalised using the Rotation-Translation-of-Blocks method. 65 All but the 50 lowest frequency non-trivial modes were discarded, with the remaining modes stored as unit vectors, of length N , sorted in ascending order of their corresponding eigenvalues.

4.3

Correspondence Between Residues in the Start and Final States

Corresponding residues in the start and final states were found using Needleman-Wunsch sequence alignment. Unaligned residues were discarded. The atomic coordinates of the remaining C α atoms were converted to vectors u and b for the start/unbound and final/bound structures respectively, with elements ui and bi . These are of the form hx1 , y1 , z1 , x2 . . . zn−1 , xn , yn , zn i, where n is the number of matching residues, and the x, y and z coordinates of the C α atom in the ath matched residue are denoted xa , ya and za respectively. The displacement vector, which represents the conformational change, was calculated as d = b − u.

4.4

Correspondence Between Normal Modes and Conformational Change

The normal modes were calculated using all N atoms, while the conformational change was specified by the coordinates of the n C α atoms appearing in the structures of both the start and final states. To reconcile this disparity the 50 modes were projected into fitting space by discarding the motions of atoms that do not appear in the displacement vector d. The resultant mode vectors are of the form h∆x1 , ∆y1 , ∆z1 , ∆x2 . . . ∆zn−1 , ∆xn , ∆yn , ∆zn i, where ∆ba is the component of the mode for the ath C α atom in the b direction, so that the ith element of this vector corresponds to the same atom and Cartesian coordinate as the ith element in d. The vectors were then collected into the 3n × 50 matrix M , such that its elements, Mij , contain the displacements along the ith coordinate in the jth normal mode. In the following section, the rows of M are denoted as Mi⋆ .

18


Page 19 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


4.5

Energy and Transition Pathway Calculation

The ith coordinate, ci , at any point along the transition path, is given as the starting coordinate plus a linear combination of normal coordinates specified by the parameter vector β, with elements βj :

c i = ui +

50 X

Mij βj

(1)

j=1

and the energy is calculated in procedure-defined units by summing the contribution from each individual mode:

E=

50 X

ωj2 βj2

(2)

j=1

where the frequency ωj2 is the eigenvalue of the jth normal mode. The transition pathway is calculated using the least absolute shrinkage and selection operator (LASSO) method 52 by fitting the normal modes, Mi⋆ ∈ R50 , to the atomic displacements, di ∈ R, i = 1 . . . 3n, while applying an ℓ1 norm constraint. The β parameters are found by solving the objective function 3n 50 X 1 X T 2 min (di − Mi⋆ β) + λ |βj | β∈R50 6n i=1 j=1

(3)

The position along the transition path is given by the coordinate λ. The starting structure corresponds to λ = λmax , the smallest value of λ for which |β| = 0, and decreases to λ = 0, corresponding to an unconstrained fit, yielding the β vector which minimizes the RMSD of coordinates c against the final state b and is equivalent to the fitting used by Moal and Bates 33 and Oliwa and Shen 66 . Starting from large λ, the pathway is calculated using coordinate descent, 67 soft thresholding the partial residual using the update 3n X 1 (j) β˜j ← S Mij (di − d˜i ), λ 3n i=1

19


(4)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 30

where

(j) d˜i =

X

Miℓ β˜ℓ

(5)

ℓ6=j

with the soft thresholding operators given by     z−γ     S(z, γ) = z + γ       0

if z > 0 and γ < |z| if z < 0 and γ < |z|

(6)

if γ ≥ |z|

The update is applied until convergence, at which point λ is decreased and the process is repeated using the previous β vector as a warm start. Calculations were performed in R using the glmnet package, 67 using up to 1000 samples of λ logarithmically spaced in the [λmax , λmax × 10−4 ] range, and terminating when the fractional change in RMSD deviance between λ samples reaches less than 10−5 , for cases in which subsequent λ samples are unnecessary to reach the final position. In order to compare the extent of deformation between proteins of different size, we normalize β by N .

Acknowledgement The authors thank Paul Bates, Emanuele Monza, Juan Fernández-Recio and Nick Goldman for providing useful comments, and Nick Goldman for hosting T.W.H during his Summer internship. This work was funded by the European Molecular Biology Laboratory and the Biotechnology and Biological Sciences Research Council (I.H.M.: Future Leader Fellowship BB/N011600/1, and T.W.H.: Research Experience Placements scheme BB/M011194/1).

Supporting Information Available The following files are available free of charge.

20


Page 21 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


• SupplementaryFigure1.pdf: RMSD and energy of calculated pathways towards true (red dotted) and contrived (green) final states. • SupplementaryText.pdf: Details of protein transition pathways that pass via experimentally determined intermediate configurations. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Tirion, M. M. Large Amplitude Elastic Motions in Proteins from a Single-Parameter, Atomic Analysis. Phys. Rev. Lett. 1996, 77, 1905–1908. (2) Bahar, I.; Atılgan, A. R.; Erman, B. Direct evaluation of thermal fluctuations in proteins using a single-parameter harmonic potential. Fold. Des. 1997, 2, 173–181. (3) Haliloğlu, T.; Bahar, I.; Erman, B. Gaussian Dynamics of Folded Proteins. Phys. Rev. Lett. 1997, 79, 3090–3093. (4) Yang, L.-W.; Eyal, E.; Chennubhotla, C.; Jee, J.; Gronenborn, A. M.; Bahar, I. Insights into equilibrium dynamics of proteins from comparison of NMR and X-ray data with computational predictions. Structure 2007, 15, 741–749. (5) Hinsen, K. Analysis of domain motions by approximate normal mode calculations. Proteins 1998, 33, 417–429. (6) Doruker, P.; Atilgan, A. R.; Bahar, I. Dynamics of proteins predicted by molecular dynamics simulations and analytical approaches: application to alpha-amylase inhibitor. Proteins 2000, 40, 512–524. (7) Atılgan, A. R.; Durell, S. R.; Jernigan, R. L.; Demirel, M. C.; Keskin, O.; Bahar, I. Anisotropy of fluctuation dynamics of proteins with an elastic network model. Biophys. J. 2001, 80, 505–515. 21



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(8) Tama, F.; Sanejouand, Y. H. Conformational change of proteins arising from normal mode calculations. Protein Eng. 2001, 14, 1–6. (9) Krebs, W. G.; Alexandrov, V.; Wilson, C. A.; Echols, N.; Yu, H.; Gerstein, M. Normal mode analysis of macromolecular motions in a database framework: developing mode concentration as a useful classifying statistic. Proteins 2002, 48, 682–695. (10) Alexandrov, V.; Lehnert, U.; Echols, N.; Milburn, D.; Engelman, D.; Gerstein, M. Normal modes for predicting protein motions: a comprehensive database assessment and associated Web tool. Protein Sci. 2005, 14, 633–643. (11) Rueda, M.; Chacón, P.; Orozco, M. Thorough validation of protein normal mode analysis: a comparative study with essential dynamics. Structure 2007, 15, 565–575. (12) Dobbins, S. E.; Lesk, V. I.; Sternberg, M. J. E. Insights into protein flexibility: The relationship between normal modes and conformational change upon protein-protein docking. Proc. Natl. Acad. Sci. U.S.A. 2008, 105, 10390–10395. (13) Meireles, L.; Gür, M.; Bakan, A.; Bahar, I. Pre-existing soft modes of motion uniquely defined by native contact topology facilitate ligand binding to proteins. Protein Sci. 2011, 20, 1645–1658. (14) Gür, M.; Zomot, E.; Bahar, I. Global motions exhibited by proteins in micro- to milliseconds simulations concur with anisotropic network model predictions. J. Chem. Phys. 2013, 139, 121912. (15) Tama, F.; Miyashita, O.; Brooks, C. L. Normal mode based flexible fitting of highresolution structure into low-resolution experimental data from cryo-EM. J. Struct. Biol. 2004, 147, 315–326. (16) Tama, F.; Miyashita, O.; Brooks, C. L. Flexible multi-scale fitting of atomic structures

22


Page 22 of 30

Page 23 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


into low-resolution electron density maps with elastic network normal mode analysis. J. Mol. Biol. 2004, 337, 985–999. (17) Hinsen, K.; Reuter, N.; Navaza, J.; Stokes, D. L.; Lacapère, J.-J. Normal mode-based fitting of atomic structure into electron density maps: application to sarcoplasmic reticulum Ca-ATPase. Biophys. J. 2005, 88, 818–827. (18) Tirion, M. M.; ben Avraham, D.; Lorenz, M.; Holmes, K. C. Normal modes as refinement parameters for the F-actin model. Biophys. J. 1995, 68, 5–12. (19) Wu, Y.; Ma, J. Refinement of F-actin model against fiber diffraction data by long-range normal modes. Biophys. J. 2004, 86, 116–124. (20) Delarue, M.; Dumas, P. On the use of low-frequency normal modes to enforce collective movements in refining macromolecular structural models. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 6957–6962. (21) Zhang, Z.; Shi, Y.; Liu, H. Molecular Dynamics Simulations of Peptides and Proteins with Amplified Collective Motions. Biophys. J. 2003, 84, 3583–3593. (22) Miloshevsky, G. V.; Jordan, P. C. The open state gating mechanism of gramicidin a requires relative opposed monomer rotation and simultaneous lateral displacement. Structure 2006, 14, 1241–1249. (23) Isin, B.; Schulten, K.; Tajkhorshid, E.; Bahar, I. Mechanism of signal propagation upon retinal isomerization: insights from molecular dynamics simulations of rhodopsin restrained by normal modes. Biophys. J. 2008, 95, 789–803. (24) Yan, Q.; Murphy-Ullrich, J. E.; Song, Y. Structural insight into the role of thrombospondin-1 binding to calreticulin in calreticulin-induced focal adhesion disassembly. Biochemistry 2010, 49, 3685–3694.

23



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(25) Gür, M.; Zomot, E.; Cheng, M. H.; Bahar, I. Energy landscape of LeuT from molecular simulations. J. Chem. Phys. 2015, 143, 243134. (26) Monza, E.; Lucas, M. F.; Camarero, S.; Alejaldre, L. C.; Martínez, A. T.; Guallar, V. Insights into Laccase Engineering from Molecular Simulations: Toward a Binding-Focused Strategy. J. Phys. Chem. Lett. 2015, 6, 1447–1453. (27) Takahashi, R.; Gil, V. A.; Guallar, V. Monte Carlo Free Ligand Diffusion with Markov State Model Analysis and Absolute Binding Free Energy Calculations. J. Chem. Theory Comput. 2014, 10, 282–288. (28) Kotev, M.; Lecina, D.; Tarragó, T.; Giralt, E.; Guallar, V. Unveiling Prolyl Oligopeptidase Ligand Migration by Comprehensive Computational Techniques. Biophys. J. 2015, 108, 116–125. (29) Lindahl, E.; Delarue, M. Refinement of docked protein-ligand and protein-DNA structures using low frequency normal mode amplitude optimization. Nucleic Acids Res. 2005, 33, 4496–4506. (30) May, A.; Zacharias, M. Protein-ligand docking accounting for receptor side chain and global flexibility in normal modes: evaluation on kinase inhibitor cross docking. J. Med. Chem. 2008, 51, 3499–3506. (31) Borrelli, K. W.; Cossins, B.; Guallar, V. Exploring hierarchical refinement techniques for induced fit docking with protein and ligand flexibility. J. Comput. Chem. 2010, 31, 1224–1235. (32) Zacharias, M.; Sklenar, H. Harmonic modes as variables to approximately account for receptor flexibility in ligand-receptor docking simulations: Application to DNA minor groove ligand complex. J. Comput. Chem. 1999, 20, 287–300.

24


Page 24 of 30

Page 25 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(33) Moal, I. H.; Bates, P. A. SwarmDock and the use of normal modes in protein-protein docking. Int. J. Mol. Sci. 2010, 11, 3623–3648. (34) Li, X.; Moal, I. H.; Bates, P. A. Detection and refinement of encounter complexes for protein-protein docking: taking account of macromolecular crowding. Proteins 2010, 78, 3189–3196. (35) Mashiach, E.; Nussinov, R.; Wolfson, H. J. FiberDock: Flexible induced-fit backbone refinement in molecular docking. Proteins 2010, 78, 1503–1519. (36) Venkatraman, V.; Ritchie, D. W. Flexible protein docking refinement using posedependent normal mode analysis. Proteins 2012, 80, 2262–2274. (37) Shen, Y. Improved flexible refinement of protein docking in CAPRI rounds 22-27. Proteins 2013, 81, 2129–2136. (38) Kim, M. K.; Jernigan, R. L.; Chirikjian, G. S. Efficient generation of feasible pathways for protein conformational transitions. Biophys. J. 2002, 83, 1620–1630. (39) Kim, P. M.; Lu, L. J.; Xia, Y.; Gerstein, M. B. Relating three-dimensional structures to protein networks provides evolutionary insights. Science 2006, 314, 1938–1941. (40) Miyashita, O.; Onuchic, J. N.; Wolynes, P. G. Nonlinear elasticity, proteinquakes, and the energy landscapes of functional transitions in proteins. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 12570–12575. (41) Best, R. B.; Chen, Y.-G.; Hummer, G. Slow protein conformational dynamics from multiple experimental structures: the helix/sheet transition of arc repressor. Structure 2005, 13, 1755–1763. (42) Okazaki, K.-i.; Koga, N.; Takada, S.; Onuchic, J. N.; Wolynes, P. G. Multiple-basin energy landscapes for large-amplitude conformational motions of proteins: Structure-

25



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

based molecular dynamics simulations. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 11844– 11849. (43) Chu, J.-W.; Voth, G. A. Coarse-grained free energy functions for studying protein conformational changes: a double-well network model. Biophys. J. 2007, 93, 3860– 3871. (44) Zheng, W.; Brooks, B. R.; Hummer, G. Protein conformational transitions explored by mixed elastic network models. Proteins 2007, 69, 43–57. (45) Yang, S.; Roux, B. Src kinase conformational activation: thermodynamics, pathways, and mechanisms. PLoS Comput. Biol. 2008, 4, e1000047. (46) Yang, Z.; Májek, P.; Bahar, I. Allosteric transitions of supramolecular systems explored by network models: application to chaperonin GroEL. PLoS Comput. Biol. 2009, 5, e1000360. (47) Zhu, F.; Hummer, G. Gating transition of pentameric ligand-gated ion channels. Biophys. J. 2009, 97, 2456–2463. (48) Das, A.; Gür, M.; Cheng, M. H.; Jo, S.; Bahar, I.; Roux, B. Exploring the conformational transitions of biomolecular systems using a simple two-state anisotropic network model. PLoS Comput. Biol. 2014, 10, e1003521. (49) Kantarcı-Çarşıbaşı, N.; Haliloğlu, T.; Doruker, P. Conformational transition pathways explored by Monte Carlo simulation integrated with collective modes. Biophys. J. 2008, 95, 5862–5873. (50) Kürkçüoğlu, Z.; Doruker, P. Ligand Docking to Intermediate and Close-To-Bound Conformers Generated by an Elastic Network Model Based Algorithm for Highly Flexible Proteins. PLoS ONE 2016, 11, e0158063.

26


Page 26 of 30

Page 27 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(51) Kürkçüoğlu, Z.; Bahar, I.; Doruker, P. ClustENM: ENM-Based Sampling of Essential Conformational Space at Full Atomic Resolution. J. Chem. Theory Comput. 2016, 12, 4549–4562. (52) Tibshirani, R. Regression Shrinkage and Selection Via the Lasso. J. R. Stat. Soc. B 1994, 58, 267–288. (53) Xu, Y.; Moseley, J. B.; Sagot, I.; Poy, F.; Pellman, D.; Goode, B. L.; Eck, M. J. Crystal structures of a Formin Homology-2 domain reveal a tethered dimer architecture. Cell 2004, 116, 711–723. (54) Otomo, T.; Tomchick, D. R.; Otomo, C.; Panchal, S. C.; Machius, M.; Rosen, M. K. Structural basis of actin filament nucleation and processive capping by a formin homology 2 domain. Nature 2005, 433, 488–494. (55) Arolas, J. L.; Popowicz, G. M.; Lorenzo, J.; Sommerhoff, C. P.; Huber, R.; Aviles, F. X.; Holak, T. A. The three-dimensional structures of tick carboxypeptidase inhibitor in complex with A/B carboxypeptidases reveal a novel double-headed binding mode. J. Mol. Biol. 2005, 350, 489–498. (56) Pantoja-Uceda, D.; Arolas, J. L.; García, P.; López-Hernández, E.; Padró, D.; Aviles, F. X.; Blanco, F. J. The NMR structure and dynamics of the two-domain tick carboxypeptidase inhibitor reveal flexibility in its free form and stiffness upon binding to human carboxypeptidase B. Biochemistry 2008, 47, 7066–7078. (57) Clarke, D.; Sethi, A.; Li, S.; Kumar, S.; Chang, R. W. F.; Chen, J.; Gerstein, M. Identifying Allosteric Hotspots with Dynamics: Application to Inter- and Intra-species Conservation. Structure 2016, 24, 826–837. (58) Vreven, T.; Moal, I. H.; Vangone, A.; Pierce, B. G.; Kastritis, P. L.; Torchala, M.; Chaleil, R.; Jiménez-García, B.; Bates, P. A.; Fernández-Recio, J.; Bonvin, A. M.;

27



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Weng, Z. Updates to the Integrated Protein-Protein Interaction Benchmarks: Docking Benchmark Version 5 and Affinity Benchmark Version 2. J. Mol. Biol. 2015, 427, 3031–3041. (59) Yang, L.; Song, G.; Jernigan, R. L. How Well Can We Understand Large-Scale Protein Motions Using Normal Modes of Elastic Network Models? Biophys. J. 2007, 93, 920– 929. (60) Grünberg, R.; Leckner, J.; Nilges, M. Complementarity of structure ensembles in protein-protein binding. Structure 2004, 12, 2125–2136. (61) Smith, G. R.; Sternberg, M. J. E.; Bates, P. A. The relationship between the flexibility of proteins and their conformational states on forming protein-protein complexes with an application to protein-protein docking. J. Mol. Biol. 2005, 347, 1077–1101. (62) Kuroda, D.; Gray, J. J. Pushing the Backbone in Protein-Protein Docking. Structure 2016, 24, 1821–1829. (63) Pallara, C.; Rueda, M.; Abagyan, R.; Fernández-Recio, J. Conformational Heterogeneity of Unbound Proteins Enhances Recognition in Protein-Protein Encounters. J. Chem. Theory Comput. 2016, 12, 3236–3249. (64) Suhre, K.; Sanejouand, Y.-H. ElNemo: a normal mode web server for protein movement analysis and the generation of templates for molecular replacement. Nucleic Acids Res. 2004, 32, W610–614. (65) Tama, F.; Gadea, F. X.; Marques, O.; Sanejouand, Y. H. Building-block approach for determining low-frequency normal modes of macromolecules. Proteins 2000, 41, 1–7. (66) Oliwa, T.; Shen, Y. cNMA: a framework of encounter complex-based normal mode analysis to model conformational changes in protein interactions. Bioinformatics 2015, 31, i151–160. 28


Page 28 of 30

Page 29 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(67) Friedman, J.; Hastie, T.; Tibshirani, R. Regularization Paths for Generalized Linear Models via Coordinate Descent. J. Stat. Softw. 2010, 33, 1–22.

29




Page 30 of 30

Modeling Protein Conformational Transition Pathways Using

Recommend Documents