Revealing Atomic-Level Mechanisms of Protein Allostery with Molecular Dynamics Simulations.

Molecular dynamics (MD) simulations have become a powerful and popular method for the study of protein allostery, the widespread phenomenon in which a stimulus at one site on a protein influences the properties of another site on the protein. By capturing the motions of a protein's constituent atoms, simulations can enable the discovery of allosteric binding sites and the determination of the mechanistic basis for allostery. These results can provide a foundation for applications including rational drug design and protein engineering. Here, we provide an introduction to the investigation of protein allostery using molecular dynamics simulation. We emphasize the importance of designing simulations that include appropriate perturbations to the molecular system, such as the addition or removal of ligands or the application of mechanical force. We also demonstrate how the bidirectional nature of allostery-the fact that the two sites involved influence one another in a symmetrical manner-can facilitate such investigations. Through a series of case studies, we illustrate how these concepts have been used to reveal the structural basis for allostery in several proteins and protein complexes of biological and pharmaceutical interest.

Introduction

Allostery—the process by which biological macromolecules transmit the effect of a stimulus (typically ligand binding) at one site to a physically distinct site—is among the most important and heavily studied properties of proteins [1,2]. Many proteins require allostery to carry out their basic functions; for example, the binding of a hormone to a receptor may cause a structural change in that receptor, allowing the receptor to bind another protein and trigger cellular signaling. Other proteins rely on allostery for regulation; for example, small molecules that bind an enzyme at sites far from the catalytic active site often tune the enzyme’s activity. Identifying the atomic-level mechanisms of allostery is critical not only for understanding protein function but also as a foundation for structure-based drug design and protein engineering [3,4].

Determining the structural basis for allostery, however, often proves challenging. Allostery generally involves structural changes that are caused by a stimulus at one site and that also affect another site. These changes may be very subtle, and their investigation is complicated by the fact that even in the absence of a stimulus, proteins fluctuate constantly between ensembles of related conformations; a stimulus generally induces allostery by shifting this ensemble [5-PLoS Comput Biol. 2014 ">7]. Such effects are often not evident from the static snapshots provided by crystallographic protein structures.

Molecular dynamics (MD) simulations have emerged as a valuable complement to experimental methods [8] in the study of allostery because they capture the motion of proteins in full atomic detail. These simulations—whose initial development several decades ago was among the achievements recognized by the 2013 Nobel Prize in Chemistry—predict the position of every atom in a biomolecular system as a function of time using Newtonian mechanics [9]. They have become substantially more powerful in recent years, thanks to improvements in computer power, improved simulation algorithms, and improved potential energy functions (known as force fields) that better approximate the underlying physics [10].

MD simulations have proven particularly useful in the investigation of allostery for two reasons. First, they provide an all-atom description of protein behavior with very high resolution in both space (sub-angstrom) and time (femtoseconds). Second, they allow a researcher to introduce carefully controlled perturbations to the biomolecular system—for example, through the addition of ligands or the application of external forces to localized regions of the protein—and to determine their effect on the protein’s structure and motion.

Here, we demonstrate how to use MD to study allostery. Several reviews have surveyed specialized simulation techniques designed to probe allosteric systems [11-13]. Our focus is instead on the fundamental concepts that guide the design of a simulation-based investigation of allostery in a biomolecular system. We hope that this article will prove useful both to structural biologists beginning to use MD and to MD practitioners turning to the study of allostery.

Determining Binding Modes of Allosteric Ligands

Our primary focus in this review is on using MD to investigate the nature of allosteric coupling within a protein or protein complex (next section). MD has also proven useful, however, in determining the binding sites and poses of ligands previously shown pharmacologically to bind a given protein and in determining putative binding sites for which new ligands might be designed. In particular, these techniques can help determine binding sites for allosteric ligands and thus contribute to the study of allostery.

Thanks to the availability of faster and cheaper computing hardware, it has recently become practical to perform MD simulations in which ligands bind spontaneously to proteins without any prior knowledge of the binding site [14-16]. The most conceptually straightforward approach is to simply perform a long simulation in which one or more ligands are positioned randomly in the water surrounding the receptor. The ligands then diffuse around the entire protein before “discovering” their binding site. Such methods have been used to uncover binding sites and poses for allosteric ligands ([17] and Case Study 1). They may capture ligand binding events even when the protein needs to change conformation locally to allow binding.

The binding of certain ligands, however, takes place on timescales too slow to be captured by simulation, even on supercomputers. In some cases, one can circumvent this problem by using large numbers of short simulations together with statistical modeling techniques that describe long-timescale events as sequences of short-timescale events [14,18]. Alternatively, if the binding site is known, one can position the ligand within the binding site at the beginning of the simulation, allowing its pose and the local protein conformation to adapt during the simulation [19-21].

MD can also be used to find pockets on a protein that might represent targets for ligand design. In particular, some drugs bind to “cryptic” pockets that are not visible in a ligand-free crystal structure [22]. Such pockets may form during simulation, particularly if drug-like ligands are included in the simulation to help induce pocket formation [16,23-26]. These techniques may be used to discover binding sites different from the native binding site. Such alternative sites are often referred to as “allosteric,” sometimes even in cases when binding of ligands at these alternative sites does not have an allosteric effect on the native binding site.

Determining the Nature of Allosteric Coupling

A more widespread application of MD simulations involves determining whether, why, and how a stimulus at one site on a protein or protein complex causes effects at other sites. In many cases, two particular sites have been shown experimentally to be allosterically coupled—for example, binding of a ligand at one site increases or decreases affinity or activity at another site—and one may wish to explain the molecular mechanism underlying this coupling, thus enabling the design of additional ligands with similar or different effects. In other cases, the goal is to discover what other locations on a protein might be allosterically coupled to a particular site of interest. Despite the diversity of allosteric proteins and allosteric effects, several common principles guide their investigation by MD simulations.

Observing correlated motions

The most straightforward—and perhaps the most common—approach to studying allostery with MD simulations is simply to initiate a simulation from an experimental protein structure and then look for interdependent motions of residues at and between the sites of interest. These simulations are most often analyzed by visual inspection (using molecular rendering software such as VMD [27]) and calculation of simple measurements such as distances between pairs of atoms or dihedral angles of amino acid side chains.

Because MD simulations produce large amounts of data, and because interdependent motions can be subtle, these coupled motions can be difficult to pick out by eye. Several analysis methods have been developed specifically to find residues that are coupled allosterically (reviewed in [11,13]), using metrics such as mutual information [28,29], interaction energies [30], or correlated motion [31-34]. General MD analysis methods that automatically detect structural changes or identify structural substates may also be useful for identifying interdependent motions of interest [35-37].

Researchers sometimes speak of allosteric “pathways,” suggesting a picture in which a single continuous chain of residues between the two sites accounts for allosteric coupling. In reality, though, allostery between two sites is often mediated by difficult-to-define networks involving many loosely coupled residues or by global conformational changes involving the entire protein [38,39]. These different possibilities should be kept in mind when examining simulation results with an eye toward deciphering an allosteric mechanism.

Introducing in silico perturbations

In some cases, one can identify allostery and explain its mechanism based simply on an “equilibrium” simulation—a simulation initiated from an experimental structure, with the same bound ligands as in the experimental structure and with no external forces applied. To characterize allostery in the physical world, however, one generally needs to apply a stimulus, such as adding a ligand, and the same is often true in simulation.

Indeed, one of the strengths of MD simulation for the study of allostery is that one can readily introduce a wide variety of well-defined perturbations to a biomolecular system in silico. Most importantly, one can add or remove ligands. Removing a ligand present in the experimental structure is particularly straightforward; adding a ligand that was absent in the experimental structure is more involved and requires knowledge of the ligand binding mode, but often this can be determined by MD (see previous section) or other computational or experimental methods. Other naturally occurring stimuli have easily implemented computational equivalents. For example, one can apply a mechanical force ([40] and Case Study 4) or transmembrane voltage [41], or one can modify protein protonation states to capture the effects of local changes in acidity [42].

The introduction of “artificial” perturbations in silico may also help to identify allosteric mechanisms. For example, one might force a site or residue in a protein into a particular conformation using techniques known as steered MD ([43] and Case Study 4), biased MD [44], or targeted MD [45] in order to see how a different site or residue responds. A number of specialized methods have been developed to characterize allostery by introducing small local alterations in position or motion and then scanning for resulting changes in the behavior of neighboring regions of the protein (for example, the perturbation-response–scanning method [46,47] or pump-probe MD [48]). Several other reviews provide more detail on such methods [11,13].

Finally, mutation of residues is straightforward to implement in silico. Mutagenesis differs from the other perturbations discussed above in that it involves changing the internals of the protein under study rather than perturbing it in the way that a change in its external environment would. Connecting mutagenesis results to the allosteric mechanism of the wild-type protein is not always straightforward. Mutagenesis is a widely used experimental technique, however, and replicating or predicting the results computationally can be of value.

Comparing simulation results across multiple conditions, including perturbed and unperturbed systems, requires careful analysis. One must first identify differences between simulations under different conditions and then determine the statistical significance of those differences. Quantitative analyses, including statistical tests, are often required to identify changes that result specifically from the application of a perturbation [49].

Deciding which perturbations to make, and to which starting structures, is perhaps the most critical aspect of designing an MD-based study of allostery. We reflect on this point in the following subsection and in the case studies.

Bidirectionality and the principle of microscopic reversibility

Two physical principles are particularly useful when designing MD studies of allostery. The first is the bidirectional nature of allostery—the fact that the effects of two sites on one another are symmetrical (Fig 1). For example, if binding of ligand A at site A increases the binding affinity of ligand B at site B (a situation known as “positive cooperativity”), then the binding of ligand B at site B will increase the binding affinity of ligand A at site A. (In fact, the two affinity increases will be exactly equal, if measured as changes in the free energy of binding.) This means that if one wishes to understand the allosteric effect of site A on site B, one could instead study the effect of site B on site A. Depending on the locations of the sites and the nature of the ligands involved, the latter simulations might prove much more feasible; for example, the required simulations might be much shorter (see Case Studies 1 and 3).

Fig 1

Bidirectional symmetry in allosteric systems.

Suppose that a protein can bind two different ligands, A and B, each at their respective sites. The difference in free energy ΔG0→AB between a state with both ligands bound and one with neither ligand bound is independent of the order in which the ligands bind, so we can write it both as ΔG0→AB = ΔG0→B + ΔGB→AB and as ΔG0→AB = ΔG0→A + ΔGA→AB, implying that ΔG0→B + ΔGB→AB = ΔG0→A + ΔGA→AB. Rearranging this equation yields ΔGA→AB − ΔG0→B = ΔGB→AB − ΔG0→A. That is, the difference in the binding energy of ligand A with or without ligand B bound is the same as the difference in the binding energy of ligand B with or without ligand A bound. If one side of this equation is negative, the other must also be negative, and if one side is positive, the other must also be positive. Thus, if binding of ligand B is more favorable in the presence of A, then binding of ligand A is more favorable in the presence of ligand B (positive cooperativity). Likewise, if the presence of ligand A disfavors the binding of ligand B, then the presence of ligand B disfavors the binding of ligand A (negative cooperativity).

The second principle is microscopic reversibility, which states that, at equilibrium, a chemical system (e.g., a protein) that transitions between two structural states will follow exactly the same pathways in the forward and reverse directions [50,51]. This means that if we have two different experimental structures for a protein, we can choose to simulate the transition between them in either direction. Although the forward and reverse transitions will follow the same pathways (in the absence of an external energy input), one of them might occur much more quickly and thus be more amenable to simulation (see Case Study 2). Likewise, we can simulate a ligand-binding process to learn about the unbinding process, or vice versa [15]. Properly exploiting microscopic reversibility does require careful consideration of simulation initial conditions, although transitions initiated from well-defined physical states will generally be less sensitive to these initial conditions (for a thorough discussion, see [52]).

Limitations of MD Simulations in the Context of Allostery

Both the design of MD simulations and the interpretation of their results should take into account the limitations of these simulations, several of which we highlight here. First, the force fields used to compute the forces acting on each atom during an MD simulation represent approximations to the underlying quantum mechanics. These force fields have improved substantially in recent years, but they remain imperfect [53]. Some types of simulation results are more likely than others to be affected by force field error, as discussed elsewhere [54]; such uncertainty should be taken into account in analyzing simulation results, and these results should be interpreted in the context of all available sources of information, particularly relevant experimental data.

Second, important biomolecular processes, including ligand binding and conformational change, sometimes take place on timescales longer than the timescales accessible by simulation, which range from nanoseconds to milliseconds depending on the simulation software, methodology, and hardware employed [10,55,56]. As we show in the case studies below, careful study design can help overcome these limitations in some cases.

Third, an accurate simulation typically depends on the availability of an accurate experimental protein structure, or a good homology model, for use as an initial condition. Design of simulation studies is thus heavily influenced by the availability of experimental structures.

Case Studies

In order to illustrate the concepts described in the previous sections, we next present several case studies that used MD simulations to investigate protein allostery. We have chosen examples from our own previous work simply because we are most familiar with the details of these studies. In each example, we describe the rationale underlying simulation design, with a focus on the types of questions asked and the perturbations involved.

Case Study 1: Structural basis for allosteric modulation of a GPCR

In our first sample study, we investigated the binding modes and molecular mechanisms of allosteric modulators of a G protein-coupled receptor (GPCR) (Fig 2A) [17]. Approximately one-third of all drugs act by binding to GPCRs, generally at the same site as the native ligands of these receptors—a site known as the orthosteric site. There is a great deal of pharmaceutical interest, however, in developing allosteric modulators, which bind at other sites on GPCRs and modulate the affinity of orthosteric ligands [57,58]. A number of allosteric modulators have been identified experimentally, but before we performed this study, it was not clear how they bound to their receptors or why they affected ligand affinity at the orthosteric site.

Fig 2

Simulation schematics for case studies.

Key perturbations to the experimental structure are indicated in red text. (A) Binding and cooperativity of a negative allosteric modulator (NAM) in the M2 muscarinic receptor (M2R). The binding pose of the NAM and its cooperativity with the orthosteric ligand are probed by performing unguided binding simulations with an unliganded receptor and a receptor with bound orthosteric ligand. Negative cooperativity leads to stronger binding of the NAM to the unliganded M2R. Simulations indicate that cooperativity in this system is due both to direct electrostatic repulsion between cationic ligands and coupled conformational changes of the two binding sites. (B) Activation mechanism of the β2-adrenergic receptor (β2AR). Removal of the bound G protein from the active-state, agonist-bound crystal structure of β2AR leads to a spontaneous transition to the inactive state in simulation, capturing the activation process in reverse. (C) Mechanism of GPCR-catalyzed nucleotide release from a heterotrimeric G protein. Simulations of a G protein with and without bound guanosine diphosphate (GDP) suggested that nucleotide release from a G protein—which leads to G protein activation—takes place via a previously unexpected mechanism. (D) Force-induced uncoupling of a bacterial adhesin from a fibronectin fragment. MD simulations of a fibronectin–adhesin complex led to the discovery that application of stretching forces to fibronectin reduces its affinity to adhesin.

We worked with the M2 muscarinic acetylcholine receptor (M2R), a prototypical GPCR for which a substantial number of allosteric modulators have been identified experimentally and for which a crystal structure was available [59]. Because the crystal structure had an antagonist bound in the orthosteric site, we focused on several allosteric modulators known to affect the binding affinity of orthosteric antagonists, particularly the widely used antagonist N-methyl scopolamine (NMS).

To identify the binding sites and poses of the allosteric modulators, we performed a set of MD simulations in which each modulator was allowed to diffuse around the receptor until it bound spontaneously. Repeated binding simulations with several modulators showed not only that the ligands consistently bound an extracellular-facing vestibule on the receptor but also that the ligand binding mode was governed in all cases by a set of cation-π interactions between protein and ligand, rather than by previously proposed interactions involving electrostatic attraction or aromatic–aromatic contacts. To confirm our predicted ligand binding poses, we computationally designed a set of M2R mutants intended to increase or decrease allosteric modulator affinity. The affinities of two allosteric modulators for each of these mutants were then determined through radioligand binding experiments. In every case, the affinity changed in the direction predicted.

Next, we sought to determine the mechanism by which the allosteric modulators affect the affinity of the orthosteric ligand NMS. Our study included both modulators that decrease the affinity of the orthosteric ligand, known as negative allosteric modulators (NAMs), and modulators that increase the affinity of the orthosteric ligand, known as positive allosteric modulators (PAMs).

The bidirectional symmetry of allosteric systems dictates that a bound orthosteric ligand will also affect the binding of allosteric modulators: the orthosteric ligand will decrease the affinity of NAMs and increase the affinity of PAMs. We thus performed simulations of each allosteric modulator binding in both the presence and absence of the orthosteric ligand (in essence, the “perturbation” required to uncover the allosteric mechanism). Because the allosteric modulators bind and dissociate far more quickly than the orthosteric ligand, this was much more practical than performing simulations of orthosteric ligand binding in the presence and absence of allosteric modulators.

Analysis of these simulation results revealed two mechanisms of cooperativity between the allosteric modulators and the orthosteric ligand. The first, interestingly, did not require any structural change: electrostatic repulsion between charged allosteric and orthosteric ligands destabilized the binding of one in the presence of the other. The second mechanism did involve structural change: the orthosteric ligand increased the width of both the allosteric and orthosteric binding sites, and the PAMs had a similar effect; thus, the presence of either the orthosteric ligand or a PAM increased the affinity of the other.

We used our understanding of this mechanism to design an allosteric modulator with altered cooperativity: we added bulk to a modulator such that its binding was predicted computationally to require a wider geometry of the allosteric site. We then showed experimentally that this chemical modification shifted the modulator’s cooperativity with the orthosteric ligand in the positive direction. Our results may thus help guide the rational design of allosteric GPCR drugs.

Several studies have also used MD simulations to examine the allosteric effects of sodium ions at various GPCRs [60-62] and to show that various orthosteric ligands bind transiently at allosteric sites on their way to and from the orthosteric site [15, 63–65].

Case Study 2: Mechanism of GPCR activation

We next turn to a study that investigates allostery between two distant sites on a protein, mediated by a conformational change of the entire protein (Fig 2B) [66]. Here, we sought to determine the mechanism of GPCR activation—the process by which binding of certain ligands (known as agonists) to the extracellular-facing orthosteric site favors a structural change that allows the intracellular side of the GPCR to bind to a G protein. This allosteric process underlies the basic function of GPCRs, which is to stimulate intracellular signaling when agonists are detected in the extracellular environment [67]. A number of other MD studies have addressed similar questions in a variety of GPCRs [68-72].

We aimed to discover how the agonist stabilizes the active state—that is, to identify the allosteric pathway connecting the orthosteric site to the intracellular G protein–coupling interface, as well as the nature of the coupling along that pathway. We focused on the β2-adrenergic receptor (β2AR), an archetypal GPCR for which crystal structures were available in both an inactive state (bound to various orthosteric ligands) [73] and an active state (bound simultaneously to an orthosteric agonist and a G protein or G protein–mimetic nanobody) [74,75].

A straightforward approach to studying the β2AR activation mechanism by MD would be to initiate simulations from the receptor’s inactive state, with an agonist placed in the orthosteric site. Unfortunately, the estimated timescales of receptor activation are many milliseconds, substantially beyond the timescales that were accessible to MD simulations.

We found, however, that simulations of an agonist-bound receptor initiated from the active state, with the G protein removed, transitioned spontaneously to the inactive state within a few microseconds. Indeed, a crystallographic structure of an agonist-bound β2AR with no intracellular binding partner was also inactive [76]; the intracellular binding partner plays a key role in stabilizing the fully active β2AR conformation.

Exploiting the principle of microscopic reversibility, we reasoned that simulations of receptor deactivation—the transition of an agonist-bound receptor from the active state to the inactive state—could be used to infer the mechanism of the reverse process, activation, in which an agonist-bound receptor transitions to the active state. We thus performed and analyzed several dozen such simulations and identified residues that acted as connectors between the extracellular-facing orthosteric site and the receptor’s intracellular face. Interestingly, structural changes to these different regions did not occur simultaneously or in a well-defined sequence. Instead, the various regions were loosely coupled, with each one more likely to assume an active-like conformation when its neighbors did so. We also found that the transition between the receptor’s active and inactive states generally proceeded through a well-defined, metastable intermediate state.

One caveat is that we did not simulate G protein dissociation, an event that would take place on far longer timescales, and our simulations may thus not reflect activation and deactivation at equilibrium. The receptor might follow a different pathway between inactive and active states in the presence of a G protein, particularly if such a binding partner is necessary to induce early conformational changes on the activation pathway. Experimental data suggests, however, that activation-like conformational fluctuations of the agonist-bound receptor take place even in the absence of an intracellular binding partner [74,77]. In addition, a subsequent nuclear magnetic resonance (NMR) study confirmed several of our computational predictions, including the loose nature of the allosteric coupling involved in activation and the existence of an intermediate state with properties similar to those suggested by our simulations [77].

Case Study 3: Mechanism of nucleotide release in heterotrimeric G proteins

We next shift our focus to a study of allostery in a receptor–G protein complex, in which we asked how the interaction of a receptor (GPCR) with a heterotrimeric G protein accelerates GDP release from the G protein (Fig 2C) [78]. This process—which leads directly to nucleotide exchange and G protein activation—has long been viewed as a cornerstone of cell signaling, but its mechanism had remained unclear [79]. Simulations implicated a mechanism substantially different from those hypothesized previously.

The existence of crystal structures representing both a GDP-bound G protein [80] and a nucleotide free receptor–G protein complex [75], with distinct structural differences in the conformation of the G protein, made MD a particularly attractive method for studying the mechanistic basis for this process. The former crystal structure shows two domains, the Ras domain and the helical domain, pinched together around GDP. In the latter structure, the two domains appear to have swung apart in a hinge-like fashion. Based on this structural data, one might hypothesize that binding of a receptor to a G protein forces domain separation, thus enabling nucleotide release.

Surprisingly, when we initiated simulations from the closed, GDP-bound G protein structure, the helical and Ras domains separated spontaneously and frequently. This separation was more than sufficient to clear a pathway for GDP exit, but GDP remained bound. (Smaller interdomain motions had previously been observed in shorter MD simulations of G proteins [81-85].) The domains separated further in simulation when GDP was removed, but the hinge-like domain opening motion did not show any clear evidence of allosteric coupling to the receptor-binding interface of the G protein. These observations suggest that the receptor does not act on the G protein by increasing domain separation.

So what does the receptor do to accelerate nucleotide release from the G protein? Directly simulating the process of receptor–G protein association was not feasible, as it occurs on timescales longer than those we could access in simulation. By exploiting the principle of bidirectionality, however, we were able to address this question without directly simulating the association process.

In particular, physical principles imply that removal of GDP will favor G protein conformations that exhibit lower affinity for GDP over those that exhibit higher affinity for GDP. We thus compared simulations of a G protein with and without GDP bound, all initiated from the same structure solved in the absence of a receptor. We focused on the behavior of structural features on the G protein at the receptor-binding interface. We found that nucleotide removal favors the spontaneous transition of a Ras domain helix (α5) to a conformation seen in the receptor–G protein complex structure; this conformation presumably allows the G protein to couple more tightly to the receptor.

These simulations suggest that receptor binding causes a subtle conformational change within the Ras domain, reducing its affinity for GDP. This allows the bound GDP to escape when the domains separate spontaneously, in a manner that is independent of receptor binding. Simulations thus implicated a mechanism substantially different from those previously hypothesized on the basis of available experimental data.

Biophysical experiments, including several designed specifically to validate our computational results, supported this novel mechanism. Double electron–electron resonance spectroscopy experiments indicate that, in the presence of GDP, the distribution of distances between the helical domain and Ras domain has two distinct peaks, suggestive of an equilibrium between closed and open conformations [78,86]. To test the prediction that spontaneous domain separation plays an important role in nucleotide release, our collaborators engineered a G protein in which the Ras and helical domains are tethered together [78]. Kinetics experiments show that, in the absence of a GPCR, nucleotide release takes place about 20 times more slowly in this engineered G protein than in the wild type.

Case Study 4: Mechanosensitive interaction between bacterial adhesins and fibronectin

Our final case study [87] differs from the previous ones in two ways (Fig 2D). First, the allosteric phenomenon itself was discovered by simulation before being verified experimentally. Second, the allostery involves regulation of a protein–protein interaction by a mechanical force rather than by binding of another ligand.

The large multidomain protein fibronectin plays important roles in cell adhesion, embryogenesis, and wound healing and exists in fibrillar form in the extracellular matrix of vertebrates [88]. Cells adhering to fibronectin fibers can generate mechanical forces that induce strain in individual fibronectin molecules [89], which convert tensile forces into biochemical cues [90]. Interestingly, several fibronectin domains are binding targets for pathogenic bacteria, which exploit the abundance of fibronectin in host tissue and have evolved fibronectin-binding adhesins to initiate host invasion [88].

Since fibronectin exhibits a wide range of stretched conformations in vivo, we wondered if the application of force to adhesin-bound fibrils altered fibronectin–adhesin interactions. If force applied to the ends of these fibers perturbed the distant adhesin–fibronectin binding interface, this would constitute a form of allostery. This study, then, aimed to discover whether an allosteric mechanism could govern a protein–protein interaction. Our goal was to investigate whether and how fibronectin domains that interact with bacterial adhesins undergo conformational changes in response to force and whether these putative changes lead to altered binding affinities between bacterial adhesins and fibronectin.

We attempted to answer these questions by simulating an NMR-resolved structure of a bacterial adhesin bound to a fibronectin fragment containing two domains [91], with stretching forces applied to the two ends of the fibronectin fragment. The forces used in these steered MD simulations were orders of magnitude greater than those observed biologically, allowing a force-based transition to take place on shorter, MD-accessible timescales. Nonetheless, over the course of our simulations, the fibronectin fragment stretched by an amount matching that observed in vivo as a result of cell-generated tension.

In our simulations, the fibronectin fragment transitioned into a stretched state in which the interdomain separation between the two fibronectin domains significantly increased. This was accompanied by partial unbinding of the bacterial adhesin: initially bound to both domains, the bacterial adhesin detached from one of the two domains in simulations but remained bound to the other. This transition is allosteric in nature because the forces applied at the two ends of the fibronectin fragment were distant from the binding motif from which the adhesin disengaged. In simulations of the same complex with no applied force, adhesin did not separate from fibronectin. Our conclusions from this computational analysis are in agreement with an independently conducted simulation study [92].

We used this data to make predictions about affinities within the fibronectin–adhesin complex: because bacterial adhesins bind to only a single domain when fibronectin is in a stretched state, reduction of the binding mode from bivalent to monovalent should result in lower affinities of bacterial adhesins for stretched fibronectin. Binding assays using fluorescently labeled bacterial adhesins and fibronectin fibers on stretchable silicone sheets subsequently confirmed that affinity of bacterial adhesins to fibronectin fibers varied with mechanical strain [87]. Furthermore, insertion of a few extra residues into the linker region of a bacterial adhesin abolished its mechanosensitive binding properties, because the engineered adhesin was able to maintain binding to two fibronectin domains even when their interdomain distance increased upon stretching [93].

Application of steered MD to probe force-induced allostery in other systems has demonstrated that results obtained through this technique agree well with experimental data: studied systems include integrins [94], large-conductance mechanosensitive channels [95], bacterial pili [96], and the focal adhesion proteins talin and vinculin [97,98].

Conclusion

The conceptual approaches described here for studying allostery are general and powerful enough for application to a wide variety of biomolecular systems and problems. Yet much remains to be done. We lack a simple “recipe” for determining mechanisms of allostery computationally. We also need better approaches for mining large quantities of simulation data, particularly in cases in which the changes involved are subtle. As analysis methods continue to improve and as MD simulations become faster and more accurate, we envision that insights from simulation will play an ever greater role in discovering the structural basis for allostery, designing allosteric drugs, and engineering allosteric proteins.