High pressure inhibits signaling protein binding to the flagellar motor and bacterial chemotaxis through enhanced hydration

High pressure below 100 MPa interferes inter-molecular interactions without causing pressure denaturation of proteins. In Escherichia coli, the binding of the chemotaxis signaling protein CheY to the flagellar motor protein FliM induces reversal of the motor rotation. Using molecular dynamics (MD) simulations and parallel cascade selection MD (PaCS-MD), we show that high pressure increases the water density in the first hydration shell of CheY and considerably induces water penetration into the CheY-FliM interface. PaCS-MD enabled us to observe pressure-induced dissociation of the CheY-FliM complex at atomic resolution. Pressure dependence of binding free energy indicates that the increase of pressure from 0.1 to 100 MPa significantly weakens the binding. Using high-pressure microscopy, we observed that high hydrostatic pressure fixes the motor rotation to the counter-clockwise direction. In conclusion, the application of pressure enhances hydration of the proteins and weakens the binding of CheY to FliM, preventing reversal of the flagellar motor.


Results and Discussion
CheYp active form is stable even at high pressure. First, we examined the pressure effects on phosphorylation-dependent stability of monomeric CheY active form by MD simulation. Phosphorylation of CheY at the Asp57 sidechain induces a conformational change of CheY from the inactive (Fig. 1, orange) to the active ( Fig. 1, green) form that is mainly characterized by a transition of the sidechain χ 1 angle of Tyr106 from ~ 60° to ~ −150°3 9,40 , resulting in an increase in the binding affinity of CheY for the N-terminal segment of FliM (FliM N ) 39,41 . To model the combination of the inactive (Tyr106 χ 1 ~ 60°) or active (~−150°) form, phosphorylated (CheYp) or non-phosphorylated at Asp57 (CheY), four different models of CheY (aCheYp, iCheYp, aCheY, and iCheY) were constructed and simulated for 1 μs by MD at 0.1 and 100 MPa.
In the MD simulations of CheYp started from the active form (aCheYp: χ 1 = −150°), χ 1 remained at approximately −150° at both 0.1 and 100 MPa (upper panel in Fig. 2a). Even when the MD of CheYp was started from the inactive form (iCheYp: χ 1 = 60°), χ 1 angle made a transition from 60° to −150° and remained at approximately −150° at both pressure conditions (second upper panel in Fig. 2a). The calculated free energy of the inactive form compared to the active form (ΔG M ) was higher than that of the active form by a few kcal/mol when CheY is phosphorylated (Table S1). These results clearly indicated that, for CheYp, the active form is more favorable than the inactive form in the pressure range ≤100 MPa. The observed stabilization of phospho-CheY in the active form at ambient pressure is consistent with previous experimental observations [42][43][44] . It is worth mentioning that the time scale of the active-inactive transition is slowed at 100 MPa. In the case of iCheYp, the transition at 100 MPa occurred after 120 ns, while this transition was observed within 20 ns at 0.1 MPa. The same trend was also seen in the other cases. This slowing presumably reflects the pressure-induced slowdown of molecular motion and a consequent stiffening [45][46][47] .
In the case of non-phosphorylated CheY, the active form is more destabilized. Non-phosphorylated CheY (lower two panels of Fig. 2a) tended to make conformational transitions between the active and inactive forms, with the exception of iCheY at 0.1 MPa. Starting from the active form without phosphorylation (aCheY) at 0.1 MPa, CheY exhibited conformational transitions between the two forms and eventually shifted predominantly to the inactive form, while at 100 MPa, CheY mostly remained in the active form; this pattern presumably resulted the slowed conformational change at high pressure. MD of non-phosphorylated CheY started from the inactive (iCheY) underwent no such transition and retained the inactive form at 0.1 MPa; few transitions were Initial MD structures of the phosphorylated CheY in the active (aCheYp: green) and inactive (iCheYp: orange) forms. Y106 and the phosphorylated D57 residues are shown in stick models and those of aCheYp are colored on an atomic-color basis. A close-up view of Y106 is shown on the right. The FliM N structure complexed with aCheYp is shown in magenta. In this paper, the molecular structure was visualized using VMD 82 .
observed at 100 MPa. The free energy difference ΔG M was significantly dependent on the initial form unlike that of the phosphorylated CheY (Table S1). Overall, the MD simulation of monomeric CheY clearly indicated that activated phospho-CheY is not affected by high pressure below 100 MPa.
Pressure causes significant conformational change of CheY at 100 MPa. In the case of aCheYp at 100 MPa, we observed an interesting conformational transition of the N-and C-terminal segments during MD (Fig. 2b); this transition began at approximately 0.5 μs and was completed at approximately 0.74 μs. The N-terminal conformational change involved the first five residues of CheY, a region that is far from the FliM N -binding site; the C-terminal conformational change occurred as a significant bending of the helix caused by large mainchain dihedral changes in Lys119 and Leu120, as revealed by Dihedral Transition Analysis (DTA) 48 . Since CheY residues 122, 123, and 126 are part of the interface residues of the complex, the binding affinity may be weakened by loss of some interactions if this conformational transition occurs upon binding; however, the bending occurred in the direction opposite to the FliM N -binding site (compare Figs. 1a and 2b), and this conformation is not expected to sterically interfere with FliM N binding. In iCheY, a similar N-terminal change in the first five residues was observed at 100 MPa. In addition, we also identified a "sidechain flip" 48 of Asp74, situated in at the end of the α-helix and far from the FliM N -binding interface. www.nature.com/scientificreports www.nature.com/scientificreports/ Pressure induces earlier detachment of FliM n . Next, we investigated pressure effects on the aCheYp-FliM N complex by 1-μs MD at 0.1, 50, and 100 MPa. During the 1-μs MD, the complex structure was stable and its dissociation was not observed. Experimentally, dissociation of some protein complexes has been shown to occur at ~ 60 MPa 14,49 , but simulating the dissociation process using standard MD is very difficult because the time scale of the complex dissociation can be much longer than the MD time scale 26 . Some of the residues between residues 83-125 of CheY make contacts with FliM N . The complex structure is stabilized by a salt bridge between Lys119 of CheY and Asp12 of FliM N (Fig. S1) that was maintained during 90% of the simulation time, and also by a few mainchain-sidechain hydrogen bonds. Since the aforementioned conformational transition at 100 MPa starts from Lys119, that change may weaken the binding if the transition occurs. However, such a conformational change was not observed in the simulation of the complex at high pressure nor in the following dissociation simulations.
The dissociation process of the aCheYp and FliM N complex was observed by PaCS-MD at 0.1, 50, and 100 MPa. In the present work, we focused on examining the pressure dependence of the binding affinity; five PaCS-MD trials were sufficient for this purpose, as will be shown below. In all the trials, dissociation up to 50 Å in the inter-center of mass distance, d, was observed within 400 cycles (Fig. 3a). This computational cost corresponds to a total simulation time of 40 ns (0.1 ns MD per replica × 400 cycles) and total computational cost of 400 ns (40 ns × 10 replicas). Relatively large variations in the number of PaCS-MD cycles required to cause the dissociation were consistent with the results of the earlier paper 31 , which showed that this number varied in different trials when 10 replicas were used. The generated dissociation pathways are shown in Fig. 3b,c, and in Movies S1-S3. Although the direction of dissociation varied, the helical segment of the C-terminal segment of FliM N always detached from aCheYp first, with subsequent unbinding of the N-terminal region (Fig. 3c-e). The former step, consisting of detachment of the FliM N C-terminal helix, occurred at d = 25 ± 2, 25 ± 2, and 27 ± 2 Å at 0.1, 50, and 100 MPa (mean ± standard deviation (SD)), respectively, indicating a lack of clear pressure www.nature.com/scientificreports www.nature.com/scientificreports/ dependence. This process correlated with breaking of the aforementioned salt bridge (Lys119: CheY-Asp12: FliM N ), which occurred slightly earlier at d = 20 -21 Å. The subsequent unbinding of the N-terminal region occurred at d = 45 ± 3, 41 ± 3, and 35 ± 4 Å at 0.1, 50, and 100 MPa, respectively, suggesting that higher pressure induces earlier complete detachment of FliM N . We regard the phase before the helix detachment as the bound state (Phase 1), the phase before the complete dissociation as the partially bound state (Phase 2), and the last phase as the unbound state (Phase 3). This pressure-dependent dissociation is investigated further below.
In addition, we examined pressure effects on monomeric FliM N by analyzing MD trajectories generated by two different procedures. In the first procedure, MD simulations were started from the FliM N structure taken from the crystal structure of the complex (PDB ID: 1F4V 39 ) and solvated as a monomer. MDs lasting 1 μs were conducted at 0.1, 50, and 100 MPa and analyzed. In the second procedure, the FliM N structures completely dissociated from aCheYp (d > 40 Å) by PaCS-MD were analyzed. FliM N in the former tended to maintain the helical region formed in the CheY-bound state, whereas in the latter, FliM N was slightly more unstructured in the helical region, a change that was induced during the detachment process.
Thermodynamic properties of the proteins are conserved up to 100 MPa. Various properties of FliM N and aCheYp in monomeric and complex states were calculated at different pressures (Table S2). The root-mean-square deviation (RMSD) of α-carbons (C α ) of the MD-representative structure at high pressures compared to that at 0.1 MPa was small (~1 Å) for monomeric CheY and the aCheYp-FliM N complex, except for monomeric aCheYp and iCheY at 100 MPa, in which the aforementioned conformational changes occurred. The RMSD of monomeric FliM N is larger, mainly due to the structural fluctuation of the unstructured segment at the N-terminus, which is also very flexible in the complex form. In the helical region of FliM N (residues 8-15), the RMSD was much smaller. These results indicated that hydrostatic pressures up to ~50 MPa do not significantly change the structures of the monomers and complex within the range of the simulations, but suggested that 100-MPa pressure can change the structure of monomeric CheY. The solvent accessible surface area (SASA) and the excluded volume (V ex ), the cavity volume (V cav ), and the isothermal compressibility (κ T ) also showed no significant change between the two pressure conditions (Table S2). It is worth mentioning that the TIP3P model 50 used in the current simulation has been reported to reproduce experimental compressibility up to 200 MPa 51 . Since CheY is a small protein, the V cav of CheY was very small even at 0.1 MPa and shrank slightly at 100 MPa; however, the amount of shrinkage was equivalent only to the size of one water molecule, which is smaller than the observed cavity and volume fluctuations. The pressure-independent behavior of κ T is consistent with the fact that protein structure does not significantly change in this pressure range, indicating that high pressure mainly compresses water. This result is considered to be reasonable, given that the κ T of bulk water is ~ 2-fold larger than that of typical proteins 10,52 . This observation also is consistent with the Young's moduli of hydrated proteins, which are in the range of 1 ~3 GPa for actin, tubulin, and flagellin 53 . The Young's modulus of CheY should be equal to or greater than those of these proteins, given that CheY is a smaller and more-compact globular protein.
Expected elastic deformation at 100 MPa should be 3.3% if Young's modulus is 3 GPa, which suggests very small deformation comparable to thermal fluctuations (~1Å).

Pressure increases protein hydration and induces water penetration into the interface. Next,
we analyzed the hydrogen bond network, including water molecules. The applied pressure did not alter the number of hydrogen bonds connecting protein-protein (H PP ) that sums up both intra-and inter-protein hydrogen bonds (Table S2). The result is consistent with the aforementioned result that high pressure ≤ 100 MPa does not change protein structure in most cases. Similarly, the number of protein-water-protein hydrogen bond network (H PWP ) was almost pressure-invariant. On the other hand, the number of protein-water hydrogen bonds (H PW ) in monomeric CheY and the complex increased as the pressure increased, which agrees with the expectation that high pressure results primarily in the compression of water. Similar changes in the hydrogen bonding network at the protein-water interface previously have been reported to be induced by pressure 23 .
Next, we investigated water molecules in the first and second hydration shells around the proteins, namely the first and second shell surface waters, which are defined as the water molecules within 3 and 3-6 Å (respectively) from any protein atoms. We counted the number of surface waters in the first and second shells, N 1SW and N 2SW , respectively. Interestingly, significant increases in hydrated water molecules at high pressure were observed (Fig. 4a), consistent with the aforementioned increase in H PW (Table S2). Compared to 0.1 MPa, ~ 30 more water molecules were found in the first hydration shell of monomeric CheY and the aCheYp-FliM N complex at 100 MPa, which is equivalent to ~ 8% increase in the ratio (r 1SW ). This increase was slightly higher than the average increase in the density of bulk water at the increased pressure (~5%). Bulk water densities at 0.1, 50, and 100 MPa were calculated to be 3.29 ± 0.70, 3.37 ± 0.69, and 3.45 ± 0.67 × 10 −2 Å −3 , respectively (see Methods for detail). These results indicated that the pressure increase induces an increase of water density around the protein surfaces compared to that of bulk water. This observation is consistent with previous theoretical and experimental studies of proteins in solution at similar pressures 1,3,11,21 . This tendency of hydration increase at higher pressure also was seen in the second hydration shell (~6% increase in r 2SW ) but the increase was nominally less than that seen in r 1SW . The increase of hydration water was exceptionally higher for aCheYp in the monomernic state at 100 MPa, in which (as mentioned earlier) a conformational transition of the C-terminal helix was seen. The same tendency also was observed in the number of interface waters, which was analyzed as described below. This observation implied that the conformational transition promotes greater hydration of aCheYp at higher pressure.
We further focused on hydration around the complex interface and conducted similar analysis of the number of interface water in the first (N 1IW ) and second (N 2IW ) hydration shells around the interface residues (Table S2). These residues are defined as the residues involved in CheYp-FliM N contacts during more than 80% of the time of the second half of the 1-μs MD trajectory at 0.1 MPa. As expected, the increase in N 1IW was significant in all cases. For example, the water density of the interface in the first shell in the complex increased 6% when comparing (2020) 10:2351 | https://doi.org/10.1038/s41598-020-59172-3 www.nature.com/scientificreports www.nature.com/scientificreports/ between 0.1 MPa and 50 MPa, and 8% when comparing between 0.1 MPa and 100 MPa (see r 1IW in Table S2), while those of the second shell increased 5% and 4%, respectively (see r 2IW in Table S2). This result also showed that, in both the protein surface and protein-protein interface, the first shell water density tends to increase more than that of the second shell. High pressure increases the water density of the first hydration shell around the proteins and provokes water penetration into the interface. As shown above, high pressure did not change the timing of the salt bridge breakage, which implied that high pressure has a relatively weak effect on the salt bridge stability. Once salt bridge breakage occurred in the initial stage of dissociation, enhanced water penetration at high pressure accelerated the complete dissociation.
During the dissociation process of the complex, N 1IW and N 2IW gradually increased until the detachment of the C-terminal helix of aCheYp at the end of Phase 1, which is indicated by triangles in Fig. 4b. At this moment, the number of interface waters at 100 MPa was significantly greater than those at 0.1 and 50 MPa, which indicates the nonlinearity of the high-pressure effect on hydration at the protein-protein interface. As noted above, complete dissociation occurred at shorter d as pressure increased. This acceleration was caused by the increased hydration at higher pressure (Fig. 4c,d). At 0.1 MPa and d = 40 Å, the complex was still in Phase 2 ( Fig. 4c) but completely dissociated at the same distance at 100 MPa because more water molecules penetrated into the complex interface (Fig. 4d).
pressure impairs cheYp-fliM n binding. The free energy profile (potential of mean force as a function of d) from the bound to unbound states of CheYp and FliM N at 0.1, 50, and 100 MPa was calculated using the MSM analysis of the MD trajectories and is shown in Fig. 5a. This showed that the free energy difference between the bound and unbound states significantly decreases as the pressure increases, which clearly indicates weaker www.nature.com/scientificreports www.nature.com/scientificreports/ binding affinity at higher pressure. Although only five PaCS-MD trials were conducted, the free energy profiles showed small standard deviations, indicating that free energy converges to a certain value independent of dissociation direction if two molecules are sufficiently separated and interactions between the molecules become negligibly small, as shown in earlier work 29,31 . It should be noted that the free energy converged upon complete detachment of the complex, shown by red triangles in Fig. 5a. Since free energy differences among different pressures in Phase 1 (up to red open triangles) also were relatively small, the overall free energy differences originated mainly from the changes that occurred in Phase 2. From the point of view of the free energy profile, this difference may be interpreted as the result of weakened long-range interactions between the two proteins by more hydrated waters infiltrating between the proteins at higher pressure.
The standard binding free energy (ΔG°b ind ) of aCheYp and FliM N was calculated from the free energy profile (see Methods) and is shown in Fig. 5b. The ΔG°b ind value significantly decreased with an increase of hydrostatic pressure. This observation is qualitatively consistent with the model that protein-protein binding is suppressed by hydration of protein surfaces induced by high hydrostatic pressure 1 . The change in the ΔG°b ind value from 50 to 100 MPa (6.5 kcal/mol) was greater than that from 0.1 to 50 MPa (4.6 kcal/mol), which is consistent with the non-linear pressure response shown above. Overall, the simulation results suggested that the binding affinity of the CheYp protein to the flagellar motor protein FliM is significantly lowered by high pressure.

Pressure inhibits CW rotation of E. coli flagellar motors.
Finally, we experimentally examined the effect of hydrostatic pressure on the rotation of the flagellar motor in E. coli strain SYC12 (wild type for chemotaxis) by high-pressure microscopy. Rotation of a single motor was monitored as rotation of a cell body by attaching a single flagellar filament extended from the cell body directly to an inserted coverslip in the high-pressure chamber (Fig. 6a). Figure 6b shows typical time courses of the rotational speed of a single motor at 0.1 and 40 MPa, respectively. At 0.1 MPa, the flagellar motor rotated smoothly, and the rotational direction frequently switched from CCW to CW, and vice versa. In contrast, the motor rotated exclusively in a CCW direction at 40 MPa. The CW bias values, the time occupancy of CW rotation, at 0.1 and 40 MPa were 0.19 and zero, respectively (n = 25 cells) (Fig. 6c). On the other hand, application of 40 MPa showed no significant effect on the rotational speed (Fig. 6d). After the release of pressure, the frequent switching of the motor was observed to resume. These results indicated that application of 40 MPa of pressure suppresses motor switching, without affecting the overall structure and function of the motor. Therefore, the application of pressure is suggested to inhibit the binding of CheYp to FliM in E. coli cells, which is consistent with our observations using the simulations.
The bacterial switch complex in E. coli is a large ring comprising ~34 FliM proteins 54 . To understand how this large complex can instantaneously switch the conformations between CCW and CW rotations, the conformational spread model was proposed 35,55 . This cooperative model allows for a multistate switch in which each protomer switches its conformation non-instantaneously by biased random walk between conformational states with broadly distributed durations. Since the binding of CheYp to FliM N affects the inter-molecular interactions between FliM proteins 56 , the high-pressure inhibition of the binding is expected to reduce the probability to initiate the conformational switch to CW rotation, and also increase the probability to back to CCW rotation.
Compared to the typical time scale of the switching, the monomeric binding of CheYp to FliM N may be too strong (~30 k B T). For example, typical duration of CW rotation is ~ 0.1 s but the dissociation rate constant k off estimated with the transition state theory is 0.6 s −1 (The Eyring equation, k off = k B T/h exp(−ΔG°b ind /RT), was used with ΔG°b ind = 30 k B T. R and h are gas and Planck's constants, respectively). One possible explanation for the time scale gap is inaccurate estimation of this method. Since ΔG°b ind is related to k off , k on (association rate constant) and K d (dissociation constant) with K d = k off /k on = C 0 exp(ΔG°b ind /k B T) (C 0 = 1 molar), the value of k off can have significant variations depending on k on with a fixed value of ΔG°b ind . Another possible explanation may be also related to the results that protein-protein K d measured in vivo are one to three order greater than those in vitro 57 , suggesting that the protein-protein binding in vivo can be significantly weaker than that in vitro. The other possible reason is overestimation of ΔG°b ind in this study, however, the previous two applications of the same calculation procedure reproduced experimental binding free energy well within the error of 1 kcal/mol 29,31 , implying that the calculation error in this work is also not large. Finally, the aforementioned collective nature of the switch complex might be also related to the time gap.

conclusion
In this work, we first showed, by MD simulation, that the binding of chemotaxis signaling protein CheYp to motor rotor protein FliM N is inhibited at high pressure. Even at 100 MPa, the active form of CheYp was shown to retain its overall structure at ambient pressure, but high pressure increased the water density in the first hydration shell and caused conformational change of the C-terminal helix in the monomeric aCheYp case. Using dissociation www.nature.com/scientificreports www.nature.com/scientificreports/ simulation of the CheYp-FliM N complex by PaCS-MD and subsequent MSM analysis, we demonstrated that the binding affinity of the two proteins weakens as pressure rises, as shown by the calculated standard binding free energy. Consistent with these results, we showed, by high-pressure microscopy, that high pressure reversibly suppresses CW rotation of the bacterial flagellar motor. Altogether, our results provide a clear picture of high-pressure inhibition of CheYp-FliM N binding.
To the best of our knowledge, this is the first molecular simulation that observed the pressure-induced dissociation process of protein complexes in solution without applying bias. Since the calculation of binding free energy between proteins requires very long computational time, the pressure dependence of the binding free energy between proteins previously has been argued based on calculations for smaller methane-like molecules 58,59 . According to the heteropolymer collapse theory, hydrophobic interactions are destabilized by pressure far above 100 MPa, where the pressure-induced change in the interfacial free energy between protein and water lead to the penetration of water into a protein's hydrophobic core and result in protein denaturation 60 . In the present work, we demonstrated that the binding free energy between protein molecules decreases as pressure increases in the range ≤100 MPa. High pressure in this range induces a significant increase of water density around proteins; this increase is higher than that of bulk water density, especially in the first hydration shell, which apparently facilitates water penetration into the protein complex interface. We also observed a breakage of a key salt bridge accompanied by dissociation of the protein complex. Water molecules around such broken ionic bonds have been known to form a hydration shell that has a higher packing density than that of bulk water (i.e., electrostriction), which decreases the volume of a system and favors rupture of ionic bonds at high pressure 61 .
In conclusion, high pressure modulates a protein-protein interaction by changing its hydration state even at 100 MPa or below. This pressure range is equivalent to the pressure in the deep sea, although studies for pressure effects on proteins have focused primarily on the denaturation (or unfolding) of proteins caused by higher pressures 4,9,61 . Studies on the effects of such relatively high pressure on proteins in the biosphere are essential for understanding proteins' bioresponsive functions at the molecular level and also for designing functions governed by protein-protein interactions.

Materials and Methods
MD simulation of monomeric CheY, FliM n , and the complex. First, pressure effects on monomeric CheY, monomeric FilM N , and their complex were investigated by standard MD simulation. We prepared four different models of CheY (aCheYp, iCheYp, aCheY, and iCheY). The structure of aCheY was generated by removing FliM N from the aCheY-FliM N complex structure (PDB ID: 1F4V 39 ). The monomeric CheY crystal structure (PDB ID: 3CHY 40 ) was used for iCheY. The aCheYp and iCheYp structures were generated by modeling the phosphorylated Asp57 residue based on the crystal structure of an α-thiophosphonate derivative of the CheY D57C point mutant (PDB ID: 1C4W 42 ). The N-terminal fragment (residues 1-16) of the monomeric FliM N was obtained by removing aCheY from the crystal structure of the complex. The crystal structure of BeF 3 − -activated CheY in complex with FliM N comprising residues 1-16 shows that the first 7 residues of FliM N takes an extended conformation and residues 8-15 forms an α-helix 39 . The CheYp-FliM N complex structure was modeled based on the crystal structure (PDB ID: 1F4V 39 ) by replacing BeF 3 − with a phosphate group. After energy minimization, each of the above protein models was solvated into a cubic periodic box and gaps were filled with water molecules and 0.15 M KCl. The prepared systems contain a total of ~ 34,000, 23,000 and 38,000 atoms for monomeric CheY, monomeric FliM N , and the CheYp-FliM N complex, respectively. All of the MD simulations were performed using the GPU implementation 62 of the PMEMD module of the Amber14 package. The AMBER ff14SB force field 63 was used for the proteins and ions. The TIP3P model 50 was used for water molecules. The atomic partial charges of the phosphorylated aspartate were determined by quantum calculation using GAMESS 64 and the RESP 65 charge calculation. The other parameters (bond, angle, and dihedral) of the atoms were obtained using the general Amber force field (GAFF) 66 . These parameter files were generated by the module Antechamber of the AMBER package. The covalent bonds involving hydrogens were constrained using the SHAKE algorithm 67 , and the water molecules were kept rigid using the SETTLE algorithm 68 . After a 4000-step energy minimization, the system was brought to thermodynamic equilibrium at 300 K and 0.1, 50, or 100 MPa, using a Langevin thermostat 69 and a Monte Carlo barostat. Equations of motion were integrated with a time step of 2 fs. The long-range Coulomb energy was evaluated using the particle mesh Ewald method 70 . The MD simulations were conducted for 1.0 μs for monomeric CheY, monomeric FliM N , and the complex of CheYp and FliM N , respectively. Two pressure conditions of 0.1 and 100 MPa were used for the simulations of monomeric CheY, and three conditions (0.1, 50, and 100 MPa) were employed for monomeric FliM N and the CheYp-FliM N complex.
The solvent accessible surface area (SASA) and the excluded volume (V ex ; defined as the inside space of SASA) for CheY, FliM N , and the CheYp-FliM N complex were calculated using the CAVE software package 71 . The calculation of SASA for the interface residues of the CheYp-FliM N complex and the other analyses were performed using the cpptraj module of the AmberTools 14 Package. The representative MD structures were determined by k-means clustering into five clusters from the last half of each MD trajectory. Secondary structure was assigned by the DSSP program 72 . The isothermal compressibility, κ T , was obtained as follows 11 .
where V, p, k B , and T represent the system volume, pressure, the Boltzmann constant, and the absolute temperature, respectively. The angle bracket here denotes the average over the last half of the simulation. To calculate the density of bulk water, MD simulations with 10,056 water molecules were performed at 300 K and 0.1, 50, and 100 MPa for 100 ns. The average density during the last 50 ns of the MD trajectories were described in the text as the bulk water density. www.nature.com/scientificreports www.nature.com/scientificreports/ The free energy difference between two forms of CheY monomers was calculated from the population of the active form P active and inactive form P inactive as Δ = Δ −Δ =− ( ) The two forms were distinguished by the χ 1 angle of Y106; the inactive form corresponds to the range from −45° to 135°, the active form corresponds to the other range.

Dissociation simulation of the CheYp-FliM n complex by PaCS-MD.
To observe the dissociation process of the aCheYp and FliM N complex, we used a PaCS-MD simulation 27,28 in which mutually-overlapping conformational trajectories are generated by cycles of parallel short MD simulations and selection of MD snapshots closer to dissociation without applying additional bias to the system. In earlier works, PaCS-MD was successful in generating natural dissociation pathways of tri-N-acetyl-d-glucosamine from hen egg white lysozyme 31 and those of the transactivation domain of the p53 protein from the MDM2 protein 29 . In the present work, the dissociation process of the complex was investigated by a procedure similar to that used in the earlier reports.
As the initial structures for PaCS-MD, the structures of the CheYp-FliM N complex after 200 ns MD at 0.1, 50, and 100 MPa were adopted. To observe the dissociation up to d = 50 Å, the simulation box was expanded to a cubic box of 83 × 83 × 83 Å 3 and the gaps were filled with water and 0.15 M KCl. The total numbers of atoms in the systems were ~ 56,000 atoms including water and ions. After a short energy minimization, the simulation system was equilibrated at 300 K and 0.1 (or 50 or 100) MPa for 10 ns. The dissociation simulation was conducted five times with different initial conformations selected every 1 ns from the last snapshot of the 10-ns MD trajectory. For each PaCS-MD cycle, ten 100-ps MD simulations were performed in parallel and the trajectories were recorded every 10 fs. The snapshots sampled every 1 ps were rank-ordered according to the inter-COM distance, d, between CheYp and the helical segment of FliM N (residues [8][9][10][11][12][13][14][15][16], and the top ten snapshots were selected for the next cycle. The first 7 residues of FliM N were excluded from the calculation of COM because we noticed in preliminary trials that this region is very flexible and its fluctuation creates 'noise' in detecting the dissociation; the inclusion of the flexible region made detection of the dissociative movements more difficult. This restriction was applied only in the selection in PaCS-MD; all residues of FliM N were considered for the calculation of the inter-COM distances in other determinations. No positional restraints were applied, thereby allowing reorientation of CheYp. In the present work, we used a larger cubic box so that FliM N was permitted to dissociate in any direction. This approach differed from that used in the earlier work, in which reorientation of a protein was restrained within a rectangular box 12 . Calculation of binding free energy by MSM. Standard binding free energies of CheYp and FliM N at pressures of 0.1, 50, and 100 MPa were calculated using the MSM analysis of the MD trajectories 30 generated by each trial. For a reasonable estimation of the transition probabilities among microstates in MSM, the trajectories from 10 independent MD simulations typically were not sufficient to achieve meaningful statistics. Instead, 10 additional MD simulations were performed, starting from the initial structures of each cycle of the PaCS-MD but with different initial velocities. This addition was repeated until sufficient statistics was achieved according to the implied time scale test 30,73 . For the 0.1 MPa condition, one set of the additional MD simulations (10 MDs) was performed for two of the five PaCS-MD trials, and two sets (20 MDs) were performed for the other three trials. At 50 MPa, one additional set (10 MDs) was performed for two trials, three sets (30 MDs) were performed for one trial, and the other two trials were analyzed by MSM without additional MD simulations. At 100 MPa, one additional MD set was performed for four trials, and two sets were performed for the rest.
Clustering and MSM construction were performed for each dissociation simulation using MSMBuilder 3.5.0 74 . The inter-COM distance between CheYp and FliM N were calculated for all of the snapshots and were clustered into 30 clusters based on the k-means. The same procedure previously was used successfully in reproducing experimentally measured standard binding free energy of a protein-peptide complex 29 . A lag time of 30 ps was selected based on the implied time scales plot as a function of lag time. Free energy profiles (potential of mean force) were first calculated from the stationary probabilities of the microstates of each trial and then averaged over the multiple trials. The bound state was defined as a region before the free energy curve became flat, i.e., 40, 40, and 25 Å of the inter-COM distance for 0.1, 50 and 100 MPa, respectively. The other region (i.e., up to 50 Å after the energy curve became flat) was defined as the unbound state. The free energy difference (ΔG bind ) was calculated from the probabilities of the bound (P b ) and unbound (P u ) states: where k B is the Boltzmann constant and T is the absolute temperature. The values of P b and P u were obtained by summing the probabilities of the microstates included in the respective state. The standard free energy of binding (Δ°G bind ) was calculated using the following expression with a volume correction term 75 : where V b is the bound volume and V° is the standard-state volume (1661 Å 3 ). ∆W is the average free energy value of the unbound state compared to that of the bound state. The value of V b was calculated as the volume of the convex hull defined by the COM coordinates of FliM N in the bound state relative to the COM of CheYp at the origin. The convex hull calculation was performed for each dissociation simulation using Qhull 76 . Snapshots sampled every 2 ps were used for the V b calculation. The averaged value of V b for each pressure condition was used for the ΔG° calculation.