Path sampling simulations reveal how the Q61L mutation restricts the dynamics of KRas

The GTPase KRas is a signaling protein in networks for cell differentiation, growth, and division. KRas mutations can prolong activation of these networks, resulting in tumor formation. When active, KRas tightly binds GTP. Several oncogenic mutations affect the conversion between this rigid state and inactive, more flexible states. Detailed understanding of these transitions may provide valuable insights into how mutations affect KRas. Path sampling simulations, which focus on transitions, show KRas visiting several states, which are the same for wild type and the oncogenic mutant Q61L. Large differences occur when converting between these states, indicating the dramatic effect of the Q61L mutation on KRas dynamics. For Q61L a route to the flexible state is inaccessible, thus shifting the equilibrium to more rigid states. Our methodology presents a novel way to predict dynamical effects of KRas mutations, which may aid in identifying potential therapeutic targets. Author summary Cancer cells frequently contain mutations in the protein KRas. However, KRas is a challenging target for anti-cancer drugs, in part because the dynamic behavior of flexible regions in the protein is difficult to characterize experimentally, and occurs on timescales that are too long for straightforward molecular dynamics simulations. We have used path sampling, an advanced simulation technique that overcomes long timescales, to obtain atomistic insight into the dynamics of KRas. Comparing the oncogenic mutant Q61L to the wild type revealed that the mutation closes off one transition channel for deactivating KRas. Our approach opens up the way for predicting the dynamical effects of mutations in KRas, which may aid in identifying potential therapeutic targets.

Ras structure and function. Structure of GTP-bound KRas in the active state 2 (left). The switch regions are highlighted in green for S1 and blue for S2, the α3-helix is highlighted in yellow. The protein is shown as a ribbon with an transparent stick representation for the amino acids in S1 and S2. GTP is shown as solid sticks, with carbon atoms colored in green, oxygen in red, nitrogen in blue and phosphorus in orange. Mg 2+ is shown as a green ball. Schematic representation of the inactive state 1 and the active state 2 of GTP-bound KRas (right). The S1 region is represented in green, the S2 region in blue, and the rest of the protein in grey. State 1 corresponds to the conformational state in which S1 and S2 are more flexible and not bound to GTP. State 2 corresponds to the conformational state in which both S1 and S2 are bound to GTP. State 2 activates downstream effectors like RAF and PI3K by binding them.
transformation, drive oncogenesis and promote tumor maintenance. The Ras family of 5 oncoproteins has been studied extensively for almost three decades, as activation of Ras 6 represents a key feature of malignant transformation for many cancers. In the cancers 7 that contribute most heavily to worldwide mortality, Ras mutations are extremely 8 common [3]. Many Ras-activating mutations have been detected in non-small cell lung 9 cancer (15 to 20% of tumors), colon adenomas (40%) and pancreatic adenocarcinomas 10 (95%) making it the single most common human oncoprotein. Binding of guanosine 11 triphosphate (GTP) activates signal transduction by Ras proteins, while their GTPase 12 function inactivates signal transduction again by hydrolyzing GTP to guanosine 13 diphosphate (GDP). Ras proteins consist of a highly conserved catalytic domain called 14 the G domain and a variable C domain which anchors Ras in the membrane. Several 15 isoforms exist of Ras, which are implicated in different types of cancer [2,3]. A member 16 of this family, KRas-4B, is often found in common and life-threatening cancers, such as 17 lung cancer, colon cancer and pancreatic cancer [3]. 18 In this work we focus on the G domain of the KRas-4B isoform, which contains 166 19 residues and can be considered as the minimal signaling unit. This domain contains the 20 guanine nucleotide binding site and two regions that sense the nature of the bound 21 nucleotide, switch 1, S1, and switch 2, S2. Literature has not reached consensus on 22 which range of residues correspond to each switch region [4,5]. We chose to use a 23 narrow definition that corresponds to the residues that are important for the 24 conformational changes in this study, using residues 30-33 for S1 and residues 60-66 for 25 S2. These regions, highlighted in green (S1) and blue (S2) in figure 1(left), are involved 26 in many interactions between Ras and partners. In the GTP-bound state, Ras interacts 27 with downstream effectors such as the Raf and PI3K kinases [6]. After hydrolysis of As the Q61L mutation may affect the conformations of state 1 and state 2, as well as 48 the transitions between these states, methods that provide high resolution in both space 49 and time are required. Molecular simulations can provide such high detail. Previous 50 molecular simulation studies of Ras proteins focused on the dynamics of the stable 51 states and the effect of oncogenic mutations on those dynamics, as described in the 52 overview by Prakash et al. [16]. It is still an open question how Ras converts from 53 ordered to less ordered conformational states. In particular the effect of mutations on 54 these transitions is unclear. 55 Currently, a brute-force all-atom molecular dynamics (MD) investigation of the 56 mechanism and kinetic aspects of the conformational transitions in Ras is not feasible, 57 as the time scales involved are in the order of hundreds of microseconds, [17,18]. Such 58 long timescales are caused by high free energy barriers separating the stable states. One 59 way of overcoming these barriers is by employing biasing potentials that drive the 60 system towards the barrier region along a predefined reaction coordinate. For HRas, a 61 steered MD study revealed the importance of water molecules in entering the active 62 state 2. [19]. While methods employing additional potentials are well suited for 63 computing free-energy barriers and other thermodynamic properties, they often fail to 64 yield mechanistic insight at ambient conditions, as a poor choice of reaction coordinate 65 may lead to a wrong reaction mechanism, bad sampling and a poor estimation of the 66 rate constants.

67
The transition path sampling (TPS) algorithm [20] is another way to address the 68 timescale problem, which avoids these drawbacks. TPS is a Monte Carlo (MC) 69 simulation in the space of trajectories and collects an ensemble of short reactive 70 trajectories connecting a predefined initial and final state, without prior knowledge of 71 the transition state region. The speed-up gained by using TPS and related techniques is 72 tremendous. Assuming a transition rate in the order of 10 s −1 , observing a single 73 transition would require on average 100 ms of MD. In contrast, when using TPS, the 74 barrier region is sampled using MD trajectories of only tens of nanoseconds, thus 75 providing a speed up in the order of several million to a billion. Path sampling methods 76 have been extended from the original two-state formalism to be used with multiple 77 stable states [21]. At a given MC step, a TPS simulation samples a transition between a 78 specific pair of states. However, in the multiple state TPS (MSTPS) approach, it is 79 possible to generate trajectories that sample a different pair of states. The frequency of 80 this switching between transitions depends on the barrier separating different transition 81 channels. Analysis of the switching behavior can provide useful insight into the overall 82 dynamics.

83
For both S1 and S2, three stable states were identified that could be used in MSTPS 84 simulations of the wild type (WT) as well as Q61L. S2  The initial GTP-bound KRas-4B structure was constructed from the crystal structure of 107 GppNHp bound HRas (PDB-code: 4EFL) [17,18]. This was done by first using 108 homology modeling (MODELLER v9.16) [22], using sequential alignment to convert 109 HRas to KRas-4B. Then the GppNHp was manually modified into GTP, by changing 110 the nitrogen into an oxygen and removing the attached hydrogen. Finally, structures of 111 the protein and the GTP were combined into a single file. The initial structure for the 112 mutant (Q61L) was made from this structure by mutating the glutamine (Q) 61 of this 113 final structure into a leucine (L), using MODELLER [22].

114
The initial structures were put inside a dodecahedral periodic box with a minimum 115 distance between the structures and the side of the box of 1 nm. This resulted in boxes 116 with volumes of 228.154 nm 3 and 230.723 nm 3 for the wild type (WT) and Q61L, 117 respectively. The boxes were filled with TIP3P water [23]. 51 of the waters were  The initial systems were equilibrated in four steps, consisting of energy minimization, an 124 isothermal equilibration, an isothermal-isobaric equilibration and a 1 ns molecular 125 dynamics simulation. The equilibrated structures were used to run four 100 ns 126 molecular dynamics simulations for both WT and Q61L.

128
In the molecular dynamics simulations the atomic interactions were described by the 129 AMBER99SB-ILDN [24] force field, extended with optimized parameters for the 130 triphosphate chain of GTP [25]. Long-range electrostatic interactions were treated via 131 the Particle Mesh Ewald method [26]. The short-range non-bonded interactions (e.g. 132 electrostatics and Van der Waals interactions) were cut off at 1.1 nm.

133
All of the equilibration was performed with GROMACS v.4.6.5. [27]. The leap-frog 134 integrator was used with a time step of 2 fs. Temperature was kept constant at 310 K 135 using the v-rescale thermostat [28] using two temperature coupling groups: the first 136 group consisted of the protein, GTP and Mg 2+ , while the second group consisted of 137 water, Na + , and Cl − . The pressure was kept constant using the Parrinello-Rahman 138 barostat [29] at a pressure of 1 bar. All bond lengths were constrained using the LINCS 139 algorithm [30].
The 100 ns production runs were performed with OpenMM (7.1.0.dev-5e53567) [31]. 141 The constraints were changed to only constrain the bond lengths of bonds to a 142 hydrogen, using SHAKE [32], the integrator was the Velocity Verlet with velocity 143 randomization (VVVR) integrator [33] from OpenMMTools v.0.14 [34] and the barostat 144 was the Monte Carlo barostat [35]. The production simulations were run using the 145 CUDA platform of OpenMM on NVIDIA GeForce GTX TITAN X GPUs. 146 Collective variables and stable states 147 The long molecular dynamics runs were visually analyzed to identify stable states, using 148 VMD [36]. Five types of collective variable functions were used to define the stable 149 states, which are described in table 2 in Appendix 1.

150
For both WT and Q61L the relevant collective variables can be found in table 3 and 151 4, for S1 and S2, respectively. These collective variables are comprised of the collective 152 variable types described in table 2. The stable states for S1 were S1-D33, S1-30-32, and 153 S1-open, and for S2 were S2-GTP, S2-α3 , and S2-open. The definitions of all stable 154 states can be found in table 5.

155
Transition Path Sampling (TPS) 156 In the long molecular dynamics simulations some transitions spontaneously occurred.

157
These transitions were used as the starting transition path for TPS [20,37]. One TPS 158 simulation was performed for S1, starting from the S1-30-32 to S1-open transition. For 159 S2 three TPS simulations where performed, each starting from a different transition.

160
This was done for both WT and Q61L. The initial trajectories were first equilibrated 161 with a TPS simulation until the first decorrelated transition path (a transition path 162 that has no frames in common with the original path) was obtained. This decorrelated 163 path was used as the starting point for the production TPS simulations.

165
The TPS simulations were performed with 166 OpenPathSampling(0.1.0.dev-c192493) [38,39]. Multiple state TPS (MSTPS) [21] was 167 performed with an all-to-all flexible length ensemble, excluding self-transitions. 168 All-to-all means all transitions connecting two states are allowed. A self-transition is a 169 path that starts in a states and returns to that same state after crossing the boundaries 170 set by the state definitions. We used the one-way shooting algorithm [40], with uniform 171 shooting point selection. For the S1 simulations, 1000 shooting trials were performed, 172 while for each of the S2 simulations 2000 shooting trials were performed.

173
Analysis 174 All analysis of the TPS simulations was performed using the tools included in the 175 OpenPathSampling package [38,39], extended with custom Python code.

176
Matplotlib [41] was used for plotting the graphs and triangles. should be similar to the number of paths from B → A, reversed in time, which provides 204 a measure for convergence of the simulation.

205
Kinetics analysis 206 As we assume that the switching samples an equilibrium distribution, the probability P i 207 of sampling a transition i is given by: where n ij is the number of switches from i to j, and t i is the number of MC steps 209 sampling transition i. As the sum of all probabilities is equal to one i P i = 1, 210 equation 1 can be solved for all P i . From the probabilities the switching rate from i to 211 j, k ij , can be calculated by: The values for n and t are taken from the MSTPS simulations. This analysis is adapted 213 from [42].

215
Identification of conformational states 216 The crystal structure of GppNHp bound HRas (PDB: 4EFL) [17,18] was used as a 217 structural template to model the sequence of WT and Q61L KRas-4B with GTP bound. 218 With these two structures we performed four 100 ns MD simulations to explore the 219 conformational space of KRas, for both WT and Q61L. These simulations resulted in 220 the characterization of three stable states for S1, S1-D33, S1-30-32, and S1-open, and 221 three stable states for S2, S2-GTP, S2-α3, and S2-open, shown in figure 2. When S1 is 222 in the S1-D33 state, the side chain of D33 is involved in (water-mediated) hydrogen bond interactions with GTP. For the S1-30-32 state, S1 has shifted along GTP, 224 compared to the S1-D33 state, to form one or more hydrogen bonds between the side 225 chains of residues D30, E31, or Y32 and GTP. The conformations in which S1 has no 226 hydrogen bond interaction with GTP and where it is oriented away from GTP are 227 classified as the S1-open state. The S2-GTP state corresponds to the conformation of 228 S2 where it forms hydrogen bonds with GTP. Two states can occur when S2 is oriented 229 away from GTP. In the S2-α3 state, S2 has multiple interactions between its side chains 230 and the α3-helix. In the S2-open state S2 has no binding interactions with GTP. The 231 parameters for defining these states are listed in Appendix 1. Within the timescale of 232 the path sampling simulations, the conformation of S1 has little effect on the 233 conformation of S2 or vice versa. The same stable states were found for both WT and 234 Q61L, and were stable for at least 100 ns of MD.
235 Fig 2. Stable states of KRas. The stable states found for S1 and S2 are shown with the same coloring as figure 1(left). The S1-D33 state corresponds to the conformation in which D33 in S1 has a hydrogen bond with GTP. The S1-30-32 state corresponds to the conformation in which one or more hydrogen bonds occur between residues 30-32 and GTP. The S1-open state corresponds to the conformation in which S1 has no interactions with GTP and is oriented away from GTP. For S2 the S2-GTP state corresponds to the conformation in which S2 has one or more hydrogen bonds with GTP.
The S2-open state corresponds to a state in which S2 has no interactions with GTP and is oriented away from GTP. The S2-α3 state corresponds to a conformation in which S2 has no interactions with GTP, but instead has 4 interactions with the α3-helix.

Mapping conformational transitions 236
Using MSTPS, we investigated the transitions between the stable states as identified in 237 the MD simulations for S1 and S2 separately. For both S1 and S2 three pairs of 238 transitions can occur: S1-D33 ↔ S1-30-32, S1-D33 ↔ S1-open, and indicate a good acceptance ratio of 34% or higher, and an aggregate simulation time of 245 microseconds.

246
Path sampling simulations for proteins (with stochastic dynamics and diffuse 247 barriers) commonly employ the stochastic, or "one-way" shooting algorithm [37], which 248 improves the acceptance ratio. In this algorithm a trial move replaces only part of the 249 trajectory (forward or backward). Therefore, successive trajectories will have segments 250 with overlapping frames and at least two trials (one forward and one backward) are 251 needed for an accepted trajectory to have no frames in common with the original.

252
These no-overlap trajectories are referred to as "decorrelated", and are required for    To quantify the relative frequency or "population" of each transition and the switching 295 rate between transitions, we applied a kinetics analysis approach as developed for between S2 and residues H95, Y96, Q99, and R102 of the α3-helix (see Appendix S1for 329 definitions of these distances). Note that path densities do not show stable states, as 330 the trajectories are stopped when reaching one. Comparing the two path density plots 331 indicate that the WT simulations sample a larger region on the y-axis, as the WT path 332 density extends to above 1.25 nm, while the Q61L path density is more confined to the 333 region below 1.25 nm in S2 − α3 distance.  The channel far away from the α3-helix represents a mechanism involving water 354 molecules solvating S2, resulting in S2 extending into the solvent, away from both GTP 355 and the α3-helix. The channel close to the α3-helix represents S2 moving along a 356 hydrophobic pocket on the α3-helix. In this reaction mechanism, S2 can either enter the 357 S2-α3 state by forming four contacts between S2 and the α3-helix, or by sliding along 358 Fig 7. Schematic overview of the effect of the Q61L mutation on the dynamics of Kras. (left) WT (right) Q61L. The S1 region is represented in green, the S2 region in blue for WT and red for Q61L, the α3-helix in yellow, and the rest of the protein in grey. Downstream effectors are also shown in grey. Assuming only the S2 GTP-bound state triggers the downstream effectors, the Q61L mutation alters the conformational space such that one channel to reach the open state becomes very unlikely. This would lead to either a shift in the equilibrium distribution between the open and GTP-bound state or to transitions occurring more frequently. Both of these effects would lead to an increased probability to encounter downstream effectors while in the GTP-bound state, which would trigger the downstream signaling networks.
the helix until entering S2-open. The Q61L mutation changes a hydrophilic residue to a 359 hydrophobic one, thus lowering the affinity of S2 for water. Therefore, this latter, 360 solvated, channel, which is easily accessible for the WT protein, becomes much less 361 likely for Q61L. Moreover, the mutated S2 has stronger interactions with the α3-helix, 362 as shown by the higher path density in the channel close to α3-helix ( Figure 6  binding other proteins [7], suggesting that these structured sub-states are more similar 378 to the active state 2 than the inactive state 1.

379
The more open and flexible sub-states of the S2-open state are less likely to be 380 recognized by downstream effectors. A β-strand in the PI3 kinase interacts with KRas 381 via both S1 and S2 [6], which can only occur when both S1 and S2 are in a closed 382 conformation. The more open conformations are harder to reach in Q61L, and indeed, 383 we only observe these flexible sub-states in the WT simulations, thus providing an 384 explanation for the increased probability of Q61L to bind a downstream effector. This is 385 consistent with the fact that the Q61L mutation is more prevalent in tumors involving 386 HRas than KRas. Residue 95 is a glutamine in KRas, whereas it is a histidine in HRas, 387 rendering the α3 helix slightly less hydrophilic and thus potentially enhancing the 388 interaction with the mutated S2. This prediction may be tested by repeating the NMR 389 experiment as performed by Geyer et al. [11], comparing the effect of the Q61L 390 mutation on the switching frequency. Alternatively, the lifetime of the S2-α3 state could 391 be measured using 15 N NMR spectroscopy, by labeling nitrogen atom NE2 in the Q95 392 side chain, located in the α3-helix. Finally, the α3-helix identified as important for the 393 transitions between the ordered, active state 2 and the flexible, inactive state 1 might 394 provide a new target for the development of compounds that could ameliorate the effect 395 induced by the Q61L mutation.

397
In this work we investigated the conformational space and dynamic behavior of KRas in 398 complex with GTP using path sampling. The loops in KRas that interact with GTP 399 each visit three different conformational states. Surprisingly, these conformational states 400 do not change upon introducing the Q61L mutation, located in region S2. However, the 401 mutation has a significant effect on the transitions between the conformational states of 402 region S2. While the WT protein frequently changes from one transition to another, the 403 mutant hardly changes at all. Closer examination of the various transitions revealed 404 that S2 in the WT protein is more likely to be solvated than in the Q61L mutant. The 405 Q61L mutation prevents direct solvation of S2, which is an accessible route for the WT 406 Appendix S1: Stable state definitions The stable states are defined by ranges in 418 collective variables. This appendix provides a guide to these stable state definitions. 419 Table 2 gives the types of collective variables, while Tables 3 and 4 list the collective 420 variables used to define the stable states for S1 and S2, respectively. Table 5 gives the 421 ranges in collective variable space for the stable states found for S1 and S2. 422 Table 2. List of the different collective variable types.
CV type Description Minimum distance The smallest distance between two groups of atoms. Using MDTraj [43], distances of every atom pair were calculated. The lowest is the minimum distance. Circular mean center of mass (cCOM) * The circular mean center of mass (cCOM) is a center of geometry calculation that allows for periodicity. The system is first mapped onto a cube, followed by the calculation of the center of mass, using the procedure of ref. [44]. Then the cCOM is mapped back onto the original axes.

Number of hydrogen bonds
The number of hydrogen bonds is calculated by counting how many of the possible donor-acceptor pairs form a hydrogen bond. A hydrogen bond in this code is defined by having a H donor -acceptor distance smaller than 0.25 nm and having an X donor -H donor -acceptor angle larger than 2 3 π rad. Number of water mediated hydrogen bonds The number of water mediated hydrogen bonds is calculated by first selecting all water oxygens that are within a distance of 0.35 nm from both input groups. If a water forms a hydrogen bond with both groups, it is counted as a water mediated hydrogen bond. The hydrogen bond calculation is done as described above.

Number of bonds
The number of bonds is the number of pairs, from a given list of pairs, for which the minimum distance is smaller than 0.35 nm. The minimum distance is defined in the minimum distance cv type. * Note that this circular mean center of mass does not give the actual center of mass. In the CV type column the name of the collective variable type as used in table 3 and 4 are shown, with their description in the Description column. Table 3. List of the relevant collective variables for the stable state definitions of S1.
CV Description d GTP asp30 a The minimum distance between the C γ of aspartic acid 30 and the heavy atoms of GTP, including Mg 2+ . d GTP glu31 a The minimum distance between the C δ of glutamic acid 31 and the heavy atoms of GTP, including Mg 2+ . d GTP tyr32 a The minimum distance between the side-chain oxygen of tyrosine 32 and the heavy atoms of GTP, including Mg 2+ . d GTP asp33 a The minimum distance between the C γ of aspartic acid 33 and the heavy atoms of GTP, including Mg 2+ . n hbonds GTP asp30 c The number of hydrogen bonds between the side-chain oxygens of aspartic acid 30 and the hydroxyl groups on the ribose of GTP. n hbonds GTP tyr32 c The number of hydrogen bonds between the hydroxyl oxygen of tyrosine 32 and the hydroxyl groups on the ribose of GTP. n hbonds tyr32 GTP c The number of hydrogen bonds between the oxygens of GTP and the hydroxyl group of tyrosine 32. n hbonds ile55 tyr40 c The number of hydrogen bonds between the backbone carbonyl of isoleucine 55 and the backbone amide of tyrosine 40. n hbonds GTP S1 c The number of hydrogen bonds between the backbone carbonyls of valine 29 and aspartic acid 30 and the hydroxyls on the ribose of GTP. n h med bonds GTP asp33 d The number of water mediated hydrogen bonds between all oxygens of aspartic acid 33 and all oxygens of GTP. n h med bonds MG asp33 d The number of water mediated hydrogen bonds between all oxygens of aspartic acid 33 and Mg 2+ . n h med bonds GTP nnb d The number of water mediated hydrogen bonds between the side-chain oxygens of aspartic acid 30, glutamic acid 31 and tyrosine 32 and all oxygens of GTP. n h med bonds MG nnb d The number of water mediated hydrogen bonds between the side-chain oxygens of aspartic acid 30, glutamic acid 31 and tyrosine 32 and Mg 2+ . a This collective variable uses the minimum distance as described in These definitions apply to both WT and Q61L. The CV column shows the collective variable names, with their description in the Description column. The minimum distance between the heavy atoms of glycine 12 and the heavy atoms of glycine 60. d gly12 gln61 a,wt The minimum distance between the heavy atoms of glycine 12 and the side-chain heavy atoms of glutamine 61. d gly12 leu61 a,Q61L The minimum distance between the heavy atoms of glycine 12 and the side-chain heavy atoms of leucine 61. d GTP glu62 a The minimum distance between the C δ of glutamic acid 62 and the heavy atoms of GTP, including Mg 2+ . d GTP glu63 a The minimum distance between the C δ of glutamic acid 63 and the heavy atoms of GTP, including Mg 2+ . d cCOM GTP S2 a,b The minimum distance between the circular mean center of mass of all atoms of residues 61 to 66 and the circular mean center of mass of all atoms of GTP, including Mg 2+ . n S2 α3 e The number of combinations between the sets of {histidine 95, tyrosine 96, glutamine 99, arginine 102} and {{61} * , glutamic acid 62, glutamic acid 63, tyrosine 64} for which the minimal distance between the side-chain heavy atoms of the residue from the first set and all heavy atoms of the residue from the second set is smaller than 0.35 nm. a This collective variable uses the minimum distance as described in These definitions apply to both WT and Q61L. The CV column shows the collective variable names, with their description in the Description column. Table 5. List of the stable state definitions for KRas. * {61} is glutamine 61 for the wild type and leucine 61 for the Q61L mutant. a One or more of the conditions must be true. b All conditions must be true.
The State column are the names of the stable states. Every stable state is build by combining the Constraints and Logic columns. For example in set notation the S2-GTP state corresponds to

435
This simulation enters a different transition channel than the other simulations, which 436 would explain these altered numbers.

437
The Least Changed Path (LCP) is the sequence of frames that, together, represent 438 all accepted trajectories of a path sampling simulation. When running backwards in the 439 simulation (from the last MC step to the first) a sequence of frames that are in between 440 the shooting point of the latest trajectory and the shooting point of the trajectory in 441 which that shooting point of the latest trajectory is replaced, on the last accepted 442 trajectory before this next trajectory is added to the LCP. This is continued until the 443 first MC step is reached. These LCPs represent the barrier region that is sampled 444 during the TPS simulation. [45].  This could be attributed to a lack of convergence of the simulations. The number of 499 switches that occur between the transitions is similar for both the WT and the Q61L 500 simulation.

501
Representative trajectories of all sampled transitions of both the WT and the Q61L 502 were visually compared. No distinct differences in transition mechanisms between WT 503 and the Q61L were observed. In conclusion, these results suggest that the mutation in 504 S2 has little effect on the dynamical behavior of S1. transition for WT. The frames in these movies are rendered with the switch regions 507 highlighted in green for S1 and blue for S2. The α3-helix is highlighted in yellow. The 508 protein is shown as a ribbon with an transparent stick representation for the amino acids 509 in S1 and S2. GTP is shown as solid sticks, with carbon atoms colored in green, oxygen 510 in red, nitrogen in blue and phosphorus in orange. Mg 2+ is shown as a green ball.

511
S6 Video Movie of a typical trajectory of the S2-GTP and the S2-open 512 transition for Q61L. The frames in these movies are rendered with the switch regions 513 highlighted in green for S1 and blue for S2. The α3-helix is highlighted in yellow. The 514 protein is shown as a ribbon with an transparent stick representation for the amino acids 515 in S1 and S2. GTP is shown as solid sticks, with carbon atoms colored in green, oxygen 516 in red, nitrogen in blue and phosphorus in orange. Mg 2+ is shown as a green ball.