Large-scale, dynamin-like motions of the human guanylate binding protein 1 revealed by multi-resolution simulations

Guanylate binding proteins (GBPs) belong to the dynamin-related superfamily and exhibit various functions in the fight against infections. The functions of the human guanylate binding protein 1 (hGBP1) are tightly coupled to GTP hydrolysis and dimerization. Despite known crystal structures of the hGBP1 monomer and GTPase domain dimer, little is known about the dynamics of hGBP1. To gain a mechanistic understanding of hGBP1, we performed sub-millisecond multi-resolution molecular dynamics simulations of both the hGBP1 monomer and dimer. We found that hGBP1 is a highly flexible protein that undergoes a hinge motion similar to the movements observed for other dynamin-like proteins. Another large-scale motion was observed for the C-terminal helix α13, providing a molecular view for the α13–α13 distances previously reported for the hGBP1 dimer. Most of the loops of the GTPase domain were found to be flexible, revealing why GTP binding is needed for hGBP1 dimerization to occur.


Introduction
Guanosine triphosphate (GTP) binding proteins play essential roles in many cellular 17 processes responsible for the maintenance and regulation of biological functions. Among 18 these proteins are the guanylate binding proteins (GBPs), which belong to the 19 dynamin-related protein family, even though the GTPase domain is the only conserved 20 sequence. They have various functions in the resistance against intracellular pathogens 21 via GTP binding and hydrolysis [1][2][3][4][5]. Generally, an infection is followed by the 22 production of interferons by leukocytes, monocytes and fibroblasts, leading to 23 transcriptional activation of the interferon-stimulated genes. GBPs belong to the 24 vertebrate specific class of interferon-γ induced effector molecules that combat 25 intracellular bacteria, parasites and viruses [6]. The human guanylate binding protein 1 26 (hGBP1) was found to be involved in the defense against viruses, in particular against 27 the vesicular stomatitis virus and the encephalomyocarditis virus, and bacteria [7,8]. 28 hGBP1 was also identified as a marker of various cancer types, such as mammary 29 cancer and for cutaneous lupus erythematosus [9]. 30 hGBP1 is a large, multi-domain GTPase with similar, but quite low nucleotide 31 binding affinities for GTP, GDP (guanosine diphosphate) and GMP (guanosine 32 monoposphate) [10]. It can adopt at least two structural states with different binding 33 affinities to partner proteins. The switch between the two functional states is activated 34 by GTP binding, resulting in an 'active' state (usually the GTP/GDP-bound form) that 35 binds another hGBP1 molecule leading to dimerization or an effector protein for eliciting 36 the desired effect, and a 'silent' state that cannot bind and activate other proteins. The 37 dimerization of hGBP1 occurs through their large GTPase (LG) domains, which 38 stimulates hydrolysis of GTP to GDP and subsequently GMP in two successive cleavage 39 steps [10][11][12][13][14]. The hGBP1 monomer has been shown to be able to also hydrolyze GTP 40 to GDP, but not to GMP [15]. The crystal structure of full-length hGBP1 has been 41 solved in the nucleotide-free, i.e., the apo state (PDB 1DG3) [16] and with the Fig 1. Conformation of the nucleotide-bound hGBP1 crystal structure (PDB 1F5N). The different domains are highlighted in different colors: the LG domain in red with GTP shown in purple, the M domain in green, and the E domain in blue (with different shades used for α12 and α13 for ease of distinction). Loops missing in the crystal structure were modeled. In the enlargement of the LG domain shown on the right, the four GTP-binding site motifs are displayed: the G1-P loop in turquoise, G2-SW1 in blue, G3-SW2 (G3) in magenta, and G4-L2 in orange. The guanine cap (green), loop L1 (yellow), and residues important for dimerization or GTP binding and hydrolysis (shown as sticks with the same color as the corresponding loop) are also highlighted. The aim of the current work is to elucidate the intrinsic dynamics of apo-hGBP1 68 and the hGBP1 dimer. Given the considerable size (67 kDA, 591 residues) and 69 elongated shape of hGBP1, it is to be expected that even without nucleotide binding motion is also present in the hGBP1 dimer, for which a structural model is provided in 81 this study. Other large-scale motions were observed for the C-terminal helix α13, which 82 allows us to explain previously reported experimental data, and for most of the loops of 83 the LG domain loops, which provides a rationale why GTP is required for hGBP1 84 dimerization to occur. Our study provides fundamental insights into the dynamics of 85 both the hGBP1 monomer and dimer with consequences for its function. 86

87
To reveal the conformational dynamics of apo-hGBP1 on the microsecond time scale, we 88 performed an all-atom Hamiltonian replica exchange MD simulation (H-REMD) with 30 89 replicas of 400 ns length per replica (see Table 2 [19,20]. To explore the dynamics on 93 the sub-millisecond time scale for apo-hGBP1 and also the hGBP1 dimer, we employed 94 coarse-grained MD simulations using the Martini model [21]. We first present the which form the two-helix bundle α7/8. The residues connecting these two helices are 112 the most flexible, resulting in an RMSF peak of ∼7Å (marked as region 1 in Fig 2).

113
This two-helix bundle is followed by residues L375-T387, which are quite rigid due to 114 their proximity to the LG domain, and the third and longest helix of the M domain, α9, 115 which is characterized by monotonically increasing RMSF values up to one of the two 116 highest RMSF peaks of ∼15Å for residue K429 (region 2 in Fig 2). The next two-helix 117 bundle composed of α10 and α11 has RMSF values between 7 and 15Å, with the lowest 118 values found for the residues connecting these two helices (region 3 in Fig 2). The

119
RMSF of ∼7Å for this turn region is similar to the values for region 2 between α7 and 120 α8, which can be explained by their close geometric proximity. The second RMSF peak 121 of ∼15Å corresponds to the transition between domains M and E (region 4 in Fig 2) conformations that are furthest apart is only 7.7Å, which includes side-chain motions 135 as they were considered during the clustering analysis. Residue R48 is in all clusters 136 solvent-exposed and thus in a position that is not compatible with GTP hydrolysis.

137
This is not too surprising as it is known that only after GTP binding-which was not 138 considered in our simulations-followed by dimerization R48 is positioned toward the 139 γ-phosphate of GTP, stimulating the cleavage of this group by stabilizing the transition 140 state of GTP hydrolysis [22]. Interestingly, K51, which is also crucial for GTPase  However, the two residues S73 and T75, which are important for the hydrolysis 147 reaction [14,22], are in positions different from the ones in the LG domain dimer.

148
Without GTP, they point away from the GTP-binding site. Thus, they must undergo a 149 reorientation for adopting positions supporting GTP hydrolysis upon GTP binding and 150 hGBP1 dimerization. In contrast, residue E99, which is part of the G3-SW2 loop and 151 also of relevance for GTP hydrolysis by forming a composite base together with S73, a 152 bridging water molecule, and GTP itself enabling the transfer of a proton from the 153 water nucleophile to the GTP phosphoryl oxygen [14], remains in a position in 154 agreement with GTP hydrolysis. The stable orientation of the E99 side chain allows its 155 interaction with the Mg 2+ ion (which was not present in our simulations), while 156 residues 101-110 of G3-SW2 are flexible. Loop L1, which is not part of the GTP 157 binding site or dimerization interface, is particularly flexible between residues 156 and 158 165, while the other residues remain in their position. This is understandable as some of 159 these residues belong to an α-helix (see Fig 3), but we nonetheless included them in our 160 analysis because they interact with helix α13 from domain E. We found that the  the open to the closed state, as one might assume from the crystal structures [12,16]. energy surface (FES) of hGBP1 along these two PCs (Fig 4).  Motivated by the dynamics of the E domain observed in our H-REMD simulation and 216 by the hypothesis, derived from size-exclusion chromatography, that hGBP1 can adopt 217 an even more elongated shape than it already has by folding out the E domain (see  Fig S2A). We overlaid the 228 E domain conformations onto the crystal structure of hGBP1, which demonstrates that 229 the E domain is unlikely to fold out as an extended helix if it should detach from the 230 LG domain, as suggested from experiments [24]. The two most stable kinks or turns 231 occured at residues I489 and M505/L506, while the region between Q525 and L542 232 transformed to turns at various positions, but always reversibly as can be seen in the We aimed at understanding the reasons for turn formation that occurred in the 237 isolated α12 helix. This helix is 81 residues long and contains 21 i, i + 3 or i, i + 4 238 combinations of oppositely charged side chains (Fig 5A), which may form salt bridges. 239 While such salt bridges are known to stabilize helices [25,26], we find that in some cases 240 they can also stabilize turn structures. For instance, the interaction E526-K529-E533, 241 which is present in α12, also stabilizes the turn formed at L528 in the second half of the 242 MD trajectoy of the isolated α12 helix. Other salt bridges were newly formed upon turn 243 formation, such as the one between E488 and R493 following the helix → turn transition 244 at I489. In the helix conformation, E488 cannot form a salt bridge with residues 491 or 245 492 as these are valine and glutamate. Also the turn at M505/L506 enabled a novel salt 246 bridge to be formed, namely between E508 and R511, which in the helix is not possible 247 as the side chains of these two residues point in opposite directions. Also between E508 248 and K504 only a weak interaction is possible in α12 due to steric constraints of the helix 249 (the minimum distance between both side chains is above 6Å), while upon turn 250 formation E508 becomes able to interact with R511 instead (minimum distance below 251 2Å). Interestingly, the predicted positions of turn formation coincide with the positions 252 with a reduced α-helical coiled-coil probability of the E domain, which was determined 253 by Syguda et al. [27] using the the COILS program [28]. They performed this analysis 254 to test whether coiled-coil formation may occur between α12 helices. An alternative 255 scenario would be that the isolated E domain forms an intraprotein coiled-coil structure 256 with the first helix extending from L482 to A503 and a second helix from L506 to the 257 end of α12. A third helix with a turn around residue 540 might also be possible, which 258 would be in agreement with both our MD simulations and the COILS predictions. It 259 should be further noted that the helical propensity of the isolated E domain is only 260 slightly reduced upon turn formation. It is 79% averaged over the 500 ns MD simulation, 261 which is close to the value of 82% extracted from the circular dichroism (CD) spectrum 262 of the E domain [27], especially if one considers the errors involved in fitting CD spectra. 263 In the crystal structure (PDB 1DG3) the helical probability is 86% for α12/α13, while 264 it was 84% in the H-REMD simulation of hGBP1. Thus, the presence of the LG and M 265 domains have stabilizing effects on α12, but the overall helicity is not lost. negatively charged residues in the oppositely located α12 helix.

283
In summary, the simulations of full-length hGBP1, the isolated E domain, and α11 284 plus E domain revealed that a folding out of the E domain seems to be an unlikely event. 285 First, it would require the simultaneous breaking of several hydrogen bonds and salt 286 bridges that the E domain forms with both the LG domain and the M domain. Second, 287 even if such detachment of the E domain should happen, the missing tertiary contacts 288 would provoke turn formation at specific residues of α12 on the nanosecond time scale, 289 as demonstrated by our MD simulations. It should also be mentioned that this helix is 290 by a factor of at least four longer than the optimal length of an isolated helix (which is 291 between 9 and 17 amino acids [29]) and thus needs tertiary contacts as present in 292 hGBP1 for it to be stable. Instead of folding out of α12 as an extended helix, a more 293 likely conformational change might be the formation of a coiled-coil motif following 294 fragmentation of α12 into two or three shorter helices. In all cases studied here, α12 has 295 a high tendency of reversible kinking around residues Q525-L542. Even in full-length hGBP1 this kinking is possible, as also the M domain is very flexible in that area, which 297 is where the two-helix bundles α7/8 and α10/11 meet (marked as regions 1 and 3 in 298 Fig 2). We thus conclude that the M and E domains feature a hinge in that area. To explore the conformational flexibility of the apo-hGBP1 monomer on the micro-to 303 millisecond time scale, we performed five continuous MD simulations of lengths between 304 63 µs to 200 µs using the coarse-grained Martini force field (see Table 2). As in the 305 all-atom simulations we observed the kinking motion of the M and E domains in all five 306 coarse-grained simulations. Interestingly, this motion was possible despite the fact that 307 the Martini model preserves the initial secondary structure by applying an elastic 308 network model. Therefore, α12 remained fully helical in these simulations. It further 309 demonstrates that the application of Martini without elastic networks to preserve the 310 tertiary protein structure allows the sampling of conformational changes that also occur 311 in atomistic simulations. with the Martini force field are between 3 to 5 time faster than similar events observed 332 in atomistic simulations [21]. The main transition pathway from A to B goes via states 333 1, 2, 4, and 3 with a probability of 58%. The next most significant transition pathway 334 with a probability of 25% involves only two intermediate states, 2 and 4. The arrows indicate transitions between the states, with the line thickness correlating with the transition fluxes. The FES is colored as gradient from yellow (high free energies) to dark red (low free energies). A representative conformation is shown for each state and colored based on the three domains, while the yellow dots indicate residues Q541 and T581, which are the start and end of the sequence for which the MSM was calculated.  Fig 5B). In order to 339 identify the molecular interactions enabling the motion of α13 beyond this loop, we 340 investigated the interactions between the C-terminal region F565-T590, which includes 341 α13 (F565-M583), and L1 by calculating an average distance map during the simulation 342 interval when this motion occurred. The distance map (Fig S5A) revealed that the 343 strongest contacts are between the sequence 580 QTKMRRRKA 588 of the C-terminus 344 and the sequence 164 EVEDSAD 170 from L1. This indicates that the motion of α13 is 345 facilitated by electrostatic interactions in interplay with the high mobility of the loop L1 346 (see Fig 3), allowing the C-terminal region to move beyond L1 and then further on the 347 surface of the LG domain. We also investigated the cause for the high stability of α13 348 after it has moved beyond L1, yielding Markov states 3, 4 and B, and a very large with the second part of α13 (I576-M583) and the following, flexible residues R584-A588. 356 The strongest contacts are formed between polar residues: E146 and T149-R151 from 357 the LG domain and Q580-T581 from α13. The second interaction hot spot involves 358 sequence N220-F230, i.e., half of helix α4' that extends from S213 to F230, which is in 359 contact with residues S569-Q580 from the first part of α13. Here, the strongest 360 contacts are mainly of hydrophobic nature, involving L224 and C225 from the LG 361 domain and M572 and I576 from α13. In addition, helix α13 has residue L579 exposed 362 to the LG domain, constituting another hydrophobic interaction. Thus, we conclude 363 that the initial motion of α13 is driven by the coaction of the electrostatic interactions 364 between α13 and L1 and the flexibility of L1, and once the helix has moved beyond this 365 loop, it is stabilized by hydrophobic and polar interactions with helices α3 and α4' from 366 the LG domain.

367
Dimer model from coarse-grained simulations 368 The hGBP1 dimer is formed between GTP/GDP-bound monomers and is thought to be 369 the biologically active state of the protein [13]. Starting from the crystal structure for 370 the LG domain dimer of hGBP1 (PDB 2B92) [12] and the crystal structure of the 371 hGBP1 monomer (PDB 1DG3) [16] we built a complete hGBP1 dimer model as  Table 2). In addition to restraining the β-sheet and the 375 helices of the LG domain to their original positions (see Methods), the residues of the 376 G2-SW1 and GC loops were also restrained. They thus remained in a position of the 377 transition state during GTP hydrolysis as GDP·AlF 3 , which is present in the PDB 378 structure 2B92, is a transition-state mimic.Our Martini simulations of the hGBP1 dimer 379 can thus be considered to represent the nucleotide-bound state even though the 380 nucleotide was not present during the simulations.

381
As for the monomer, we observed the kinking motion of the M and E domains in all 382 five dimer simulations. However, the change in orientation of helix α13 as seen in the 383 coarse-grained simulations of the monomer did not occur. It remained in close proximity 384 to helix α12. Only when the temperature was increased to 320 K, as done for one 385 simulation of 270 µs length (see Table 2), we observed this conformational change, but 386 June 14, 2019 10/24 only in one of the proteins composing the hGBP1 dimer. This indicates that the energy 387 barrier for the this motion must be higher than in the monomer. To find the cause for 388 this, we analyzed whether the interactions between α13 and loop L1 or whether the 389 mobility of that loop would be different in the dimer than in the monomer. Though 390 both cases did not apply. Instead we found that a new salt bridge between α4' and α13 391 involving residues E217 and K567 is formed in the dimer, which became possible by the 392 somewhat different position that α4' adopts in the dimer than in the monomer: it is by 393 ∼6Å further away from the core of the LG domain, bringing this helix closer to α13 394 ( Fig S6A). In each of the dimer simulations at 310 K this salt bridge stayed intact, 395 thereby preventing the 90 • motion of α13. At 320 K, on the other hand, in the protein 396 of the dimer withe the flexible α13 helix this salt bridge broke (Fig S7).

397
To describe the motion of helix α13 that was observed in the simulation at 320 K, we 398 monitored the distance between the residues Q577 of α13 from the two proteins (Fig 7). 399 As can be seen from the time evolution of the distance and the snapshots associated Förster resonance energy transfer (FRET) studies of the hGBP1 dimer [30], despite the 409 fact that the experimental distances refer to distances between DEER spin labels and 410 FRET dye labels, respectively. The more important aspect is that both experimental 411 techniques predict a change in the Q577-Q577 distance by 50% (FRET) or even more 412 (DEER) upon dimer formation compared to the distance that one obtains from the 413 dimer model that is built based on the crystal structure of the hGBP1 monomer. Motion of helix α13 in the hGBP1 dimer. This motion is monitored by the time evolution of the distance between the two residues Q577 from α13 of the two proteins composing the dimer. Representative snapshots associated with different times are also shown. The two LG domains are shown in red and orange, respectively, for ease of distinction, while the M domains are shown in green and the E domains in blue. The two Q577 residues are indicated by yellow spheres.

415
We studied the conformational dynamics of the nucleotide-free hGBP1 monomer and 416 the hGBP1 dimer using multi-resolution and enhanced MD simulations on the micro-to 417 millisecond time scale. As expected from its highly conserved sequence in the dynamin 418 family, the LG domain is overall very stable in our atomistic H-REMD simulation. The 419 residues involved in the β-sheet and the adjacent α-helices have all fluctuations of less 420 than 2Å, while apart from one case the LG domain loops were found to be flexible 421 (Fig 3 and Table 1). The less flexible loop is the phosphate-binding loop G1-P and those of other enzymes [31] and only partially preorganized prior substrate 430 binding [32,33], which might explain the rather low binding affinity of hGBP1 for GTP 431 (K m = 470 µM) [10]. The LG domain loop with the highest flexibility is the guanine 432 cap (GC), which forms the protein-protein interface in the hGBP1 dimer. Even without 433 nucleotide, the GC can adopt both open and closed conformations and rapidly switch 434 between them. This finding indicates that GTP binding would shift the equilibrium 435 toward the closed GC state, which in turn would facilitate hGBP1 dimer formation via 436 the GC-GC interface requiring a certain stability of the protein recognition motif.

437
Dimer formation, in turn, would further stabilize the closed GC conformation and also 438 the orientations of R48, S73 and T75 in positions supporting GTP hydrolysis, which 439 would explain why the hGBP1 dimer is better able to hydrolyze GTP than the hGBP1 440 monomer is. Our future atomistic simulations of the GTP-bound hGBP1 monomer and 441 dimer will address whether this hypothesis holds. accompanied with local unfolding between residues Q525 and L542, while the contacts 453 to the M domain were found to be vital for overall stability of the long helix α12.

454
Without these tertiary interactions, α12 unfolds on the nanosecond time scale, making 455 it unlikely that the ∼120Å long E domain folds out as intact helix [24], which in  It should be noted that there are also certain differences between the two proteins. 478 First, the closed state is the preferred conformation of BDLP in its apo form [34], while 479 for hGBP1 the crystal structures show that the extended conformation exists for both 480  Another difference between hGBP1 and BDLP is the presence of a paddle domain in 491 BDLP at the trunk tip via which membrane binding is facilitaed, while hGBP1 binds to 492 the membrane following farnesylation at the C-terminus (see Fig 8D and E,   493 respectively). This suggests different membrane binding mechanisms for both proteins. 494 Moreover, also the multimerization seems to be faciliated via different domain 495 interactions in BDLP and hGBP1. The dimer structure of the LG domain obtained by 496 X-ray crystallography [12] implies an elongated geometry with the M and E domains 497 pointing away from the dimer interface (Fig 7), leading to a slightly curved dimer which 498 might induce membrane bending following membrane binding. Such membrane 499 destabilization might be even more likely as a result of the kinking motions of the M 500 and E domains observed in our study. In contrast, multimerization of membrane-bound 501 BDLP is thought to occur via interactions between neighboring neck and trunk helices 502 and a G domain dimer interphase different to that seen for the LG domain dimer of 503 hGBP1 (Fig 8D). Our future simulations will address hGBP1 motions following GTP 504 binding and hydrolysis, its multimerization and lipid binding, which will further 505 highlight possible differences and similarities to other dynamin-like proteins, such as 506 BDLP and MxA.

507
In the coarse-grained simulations of both the hGBP1 monomer and dimer, an in our atomistic simulation, even though an enhanced sampling scheme was employed. 513 The motion of the helix α13 has clear implications for the dimerization process. A 514 recent study has suggested that, besides the LG domain interface, the dimerization of 515 hGBP1 also involves an interface between the two helices α13 [30]. The motion of α13 516 observed in our simulations leads to an increased proximity of the two helices in the 517 dimer (Fig 7). While the distance between the two α13 helices that we monitored 518 already agrees with experimental findings, despite only one of the two helices having 519 moved, it is likely that both α13 helices adopt the 90 • rotated position in the dimer.

520
Moreover, it is not implausible that dimers are preferentially formed by monomers 521 which already have both α13 helices rotated, as the α13 motion is more likely to occur 522 in the monomer than in the dimer. This is due to a newly formed salt bridge between 523 α4' and α13 in the dimer, which increases the energy barrier for the motion of α13. The 524 rotated α13 helices would form a protein-protein interface, in addition to the LG 525 domain interface involving the two guanine caps, which would further stabilize the 526 hGBP1 dimer, as had already been suggested by Herrmann and co-workers [27,30].

527
However, based on their DEER and FRET findings they proposed that α12 and α13 528 detach from the LG domain in order to allow for the two α13 helices coming into 529 contact with each other [30]. Our simulation results demonstrate that such a 530 detachment is not necessary to explain their experimental observations. In fact, such 531 detachment would hinder the membrane binding of hGBP1 as our initial simulations of 532 membrane-bound hGBP1 indicate (Fig 8E, unpublished results).

533
In summary, to understand the conformational flexibility of hGBP1 and its All software and web databases used in this work are listed in Table S1. The input 549 needed by the various software and output files created are described in a README file 550 in the Supporting Information. For all atomistic simulations of this study the Amber99SB * -ILDNP force field [37,38] 553 combined with the explicit water model TIP3P [39] were employed. Electrostatic 554 interactions were treated with the particle-mesh Ewald method [40,41] in conjunction 555 with periodic boundary conditions and a real-space cutoff of 12Å. The Lennard-Jones 556 interactions were cut at 12Å. A leapfrog stochastic dynamics integrator was used for The crystal structure of the hGBP1 monomer in ligand-free form with PDB ID 1DG3 561 was used as starting conformation [16]. Missing amino acids from loops in the crystal 562 structure were added with the software ModLoop [44,45]. The final conformation was 563 placed in a dodecahedral box with 12Å between the protein and the box, solvated with 564 108,406 water molecules and 7 Na + ions were added for charge neutrality, resulting in a 565 system with a total number of 335,553 atoms. This particular simulation box was large 566 enough to allow free translation and rotation of the hGBP1 protein without interacting 567 with its periodic images that would otherwise result in simulation artifacts. After 568 energy minimization and equilibration of the system following the same procedure as 569 described in the next paragraph, a Hamiltonian replica exchange MD simulation [46] 570 with 30 replicas was performed. The energy function of hGBP1 including hGBP1-water 571 interactions was modified in each replica by applying biasing factors of 310 K/T with 572 the 30 temperatures T exponentially distributed between 310 and 450 K. This implies 573 one unbiased replica, the so-called target replica at 310 K. The average exchange 574 probability between the replicas was ∼30%. Each replica simulation was 400 ns long, 575 leading to a total of 12 µs for the 30 replicas. The H-REMD simulation was realized 576 with Gromacs 4.5.5 [47] in combination with the PLUMED plugin (version 2.1) [48].

577
The isolated E domain (residues 482-591) and the E domain plus helix α11 of the M 578 domain (residues 456-591) were both simulated at the all-atom level for 500 ns using 579 Gromacs 2016 [49,50]. For the E domain two Cl − ions and for α11 plus the E domain 580 one Na + ion were addded for neutrality of the systems. The energy of both systems was 581 first minimized using a steepest descent algorithm, followed by equilibration of the 582 systems to the desired temperature of 310 K and pressure of 1 atm for mimicking the 583 physiological environment. First, a 0.1 ns N V T equilibration was performed in which 584 the number of atoms (N ), the box volume (V ) and temperature (T ) were kept constant, 585 followed by a 1 ns N pT equilibration to adjust the pressure (p). During equilibration, 586 the protein atoms were restrained with a force constant of 10 kJ mol −1Å−2 allowing the 587 water molecules to relax around the solute. Finally, the 500 ns MD production runs in 588 the N pT ensemble were performed. As no restraints were applied to the protein during 589 the simulations of the E domain (with and without α11), a large cubic box with an edge 590 length of 180Å was created, allowing free rotation and translation of the 120-Å long 591 helix α12 in all directions. The resulting system sizes involved about 573,000 atoms.

592
The velocity rescaling thermostat was employed to regulate the temperature in the 593 N V T simulations, while the Nosé-Hoover thermostat [51,52] and the isotropic 594 we replaced the helix involving residues T133-F175 in the dimer with the same helix 612 from the monomer (Fig S7B). For completing the full dimer model, we combined the 613 LG domain dimer up to V288 with the monomer conformation starting from N289. The 614 resulting dimer model was first subjected to energy minimization at atomistic resolution 615 and in explicit water, and then converted to the Martini model. The coarse-grained 616 hGBP1 dimer was inserted in a rectangular box with edge lengths of 284, 111 and 108Å 617 filled with Martini water and 14 Na + ions, amounting to a final system with 32,836 618 particles. Similarly to the Martini monomer simulation, the stable parts of the LG 619 domain were restrained for both proteins in order to avoid the dimer to rotate and 620 translate, and also to preserve the conformation of the LG domains. In addition, also 621 the G2-SW1 and GC loops were restrained to keep the LG domains in its dimer-specific 622 conformation as present in PDB 2B92 [12]. As listed in Table 2, for the dimer system

626
To create pictures of the 3D protein structures, the Visual Molecular Dynamics (VMD) 627 software [55] and PyMol [56] were used. If not stated otherwise, for the analysis of the 628 H-REMD simulation the data collected by the target replica was used. To quantify the 629 stability and flexibility of hGBP1 during the atomistic MD simulations, the root mean 630 square fluctuations (RMSF) of the C α around their average positions was calculated.

631
The RMSF calculated for different time intervals of the H-REMD simulation was used 632 to demonstrate that this simulation had converged within 400 ns per replica (Fig S8). given in the Results.

640
The main structural changes of the M and E domains sampled by the target replica 641 of the H-REMD simulation were identified based on a principal component analysis 642 (PCA) [59]. Given the existence of many different structural fluctuations in hGBP1, we 643 found that applying the PCA to the entire protein is not the best way for separating the 644 different large-scale motions of the protein from each other. Therefore, we only 645 considered the Cartesian coordinates of the M domain and α12 during that analyis as 646 the RMSF analysis had revealed that these regions exhibit the largest flexibility. Helix 647 α13 of the E domain was not included as it was very stable throughout the H-REMD 648 simulation. We projected the conformations from the H-REMD target replica onto the 649 first two principal components (PCs), calculated two-dimensional histograms and then 650 the 2D free energy surface along the two PCs.

651
For the analysis of the coarse-grained simulations of the hGBP1 monomer, we 652 applied the Markov state model (MSM) approach using the PyEmma software [60]. 653 First, the five trajectories were subjected to the time-lagged independent component 654 analysis (TICA) [61], a method well suited for dimensionality reduction and recently 655 applied with success in the field of MD simulations [62][63][64][65]. The variance of the first 656 two time-lagged independent components (TICs) amounted to 27% of the total variance, 657 and they described best the conformational change involving helix α13. The MSM was 658 then built by clustering the trajectories projected onto the first two TICs using the 659 uniform time clustering algorithm with 300 microstates. We estimated the implied time 660 scales from the MSM for 10 different lag times, based on which we selected a lag time of 661 450 ns for caclculating the MSM of our system. For the identification of metastable 662 Markov states we applied the fuzzy spectral clustering method PCCA+ [66,67] and 663 used transition path theory [68][69][70], which is implemented in the PyEmma software, to 664 calculate the reactive fluxes yielding the mean first passage times between the states. between E217 and K567. In the monomer this salt bridge is never formed (black) as 717 these two residues are too far away from each other. In the dimer, helix α4' adopts a 718 slightly different position than in the monomer (see Fig S7), allowing a salt bridge being 719 formed between E217 and K567. This salt bridge impedes the motion of α13 in the 720 dimer. Only if the temperature is raised to 320 K, α13 starts moving in one of the 721 monomers composing the dimer (blue) while it remains intact in the other one (red). At 722 50 µs the salt bridge in the monomer with the flexible α13 is completely broken, which 723 corresponds to the time when this helix has adopted the 90 • rotated position (see Fig 724  7). It should be noted that the distances shown here are between the centres of the