Activation and desensitization of ionotropic glutamate receptors by selectively triggering pre-existing motions.

Ionotropic glutamate receptors (iGluRs) are ligand-gated ion channels that are key players in synaptic transmission and plasticity. They are composed of four subunits, each containing four functional domains, the quaternary packing and collective structural dynamics of which are important determinants of their molecular mechanism of function. With the explosion of structural studies on different members of the family, including the structures of activated open channels, the mechanisms of action of these central signaling machines are now being elucidated. We review the current state of computational studies on two major members of the family, AMPA and NMDA receptors, with focus on molecular simulations and elastic network model analyses that have provided insights into the coupled movements of extracellular and transmembrane domains. We describe the newly emerging mechanisms of activation, allosteric signaling and desensitization, as mainly a selective triggering of pre-existing soft motions, as deduced from computational models and analyses that leverage structural data on intact AMPA and NMDA receptors in different states.

Introduction

In order to understand the brain, one needs to understand how its cells are connected (via synaptic junctions), the nature of the signals that are sent between them, and how these properties change over time for information processing and storage. These processes can only be understood in depth through a study of the molecular events involved in synaptic communication [1]. One key family of synaptic signaling molecules is the ionotropic glutamate receptors (iGluRs), which mostly reside on the post-synaptic side of a synapse. They open cation channels upon binding of pre-synaptically released glutamate to elicit various responses depending on their subunit composition, the pattern of pre-synaptic input, and their environmental conditions [[2], [3], [4], [5]]. There are four iGluR subfamilies two of which have been particularly well studied given their expression at most excitatory synapses and hence their predominant roles in synaptic transmission and plasticity: α-amino-3-hydroxy-5-methyl-4-isozazole propionic acid receptors (AMPARs), which respond on the millisecond timescale for moment-to-moment signaling [4] and initiation of synaptic plasticity [6], and N-methyl-d-aspartate receptors (NMDARs), which act as coincidence detectors that can pass much slower but larger currents enabling a calcium influx that drives synaptic plasticity in response to repetitive stimulation [3]. AMPARs also interact with a variety of other proteins that influence their trafficking, localization and function [[7], [8], [9]]. The two other iGluR subfamilies are the kainate receptors and delta iGluRs.

A vast number of functional and structural studies have helped elucidate the molecular origin of channel responses observed in experiments, such as agonist-dependent channel opening to multiple conductance states and a process called desensitization where the channel closes with ligand still bound [8,10]. Structural studies effectively started once iGluRs were found to contain four glutamate binding sites in specialized ligand-binding domains (LBDs) that are structurally related to the clamshell-like periplasmic binding proteins (PBPs) found in bacteria and could be genetically excised and studied as separate modules [[11], [12], [13], [14], [15]]. The first structures showed that AMPAR LBDs indeed formed clamshells that close around ligands to drive channel opening akin to bacterial PBPs [16] (Fig. 1A).

Fig. 1

iGluRs are made of periplasmic binding protein (PBP) domains that close and undergo dimeric rearrangements. (A) The ligand-binding domain (LBD) is a clamshell domain that closes around the agonist glutamate. Shown are single subunits from LBD crystal structures in the ligand-free (apo) and glutamate-bound states together with their protein databank (PDB) codes. (B) LBDs form dimers with the two clamshells back to back. Closure of the two clefts pulls apart either the lower lobes to pull open the channel (activation) or the upper lobes to decouple from channel opening and remain closed (desensitization). (C) The NTD could also be produced as a separate module whose structure could be solved by X-ray crystallography. Two structures are shown representing the two main conformations and subfamilies. (D) A structure of the whole receptor is shown that is based on the resting state structure from cryo-EM (PDB 4UQJ), along with a schematic representation (on the right). Some linkers and the C-termini were modelled in with MODELLER.

The dimeric arrangement of the LBDs within a receptor, with two clamshells facing away from each other, gave rise to an initial model where cleft closure could either drive channel opening by pulling apart the lower lobes (Fig. 1B, left) or lead to desensitization by destabilizing the upper lobe dimer interface and decoupling cleft closure from channel opening (Fig. 1B, right) [16,17]. Structures of whole receptors revealed a complex multi-layered architecture with three structured domain layers (Fig. 1D) – the NTD, LBD and pore-forming transmembrane (TM) domain (TMD) – in addition to an unstructured C-terminal domain (CTD, which has not been resolved structurally to date) [[18], [19], [20], [21]]. The N-terminal domain (NTD) also belongs to the PBP superfamily and occupies different dimeric conformations (Fig. 1C) [[22], [23], [24], [25], [26], [27]], which enable it to allosterically regulate the channel in NMDARs [28,29].

In recent years, a number of agonist-bound crystal structures of whole AMPARs revealed some conformational changes without channel opening [[30], [31], [32]] and a cryo-electron microscopy study provided insights into activation and desensitization but could not clearly resolve whether and how the channel had opened [33], and the gating mechanism remained unclear. The situation has improved with recent structural characterization of auxiliary subunit complexes in resting [34,35], desensitized [36,37] and active states [36,38], as well as NMDAR structures in allosterically inhibited [39,40] and antagonist- and agonist-bound [27,41] states, but there are still a number of unanswered questions [42].

Structural data alone are not sufficient to make inferences on the mechanisms of function. Yet, the multitude of structures in different states have helped elucidate the conformational landscape accessible to these receptors. These, along with single-molecule experiments such as fluorescence resonance energy transfer (FRET) and related methods, nuclear magnetic resonance (NMR) spectroscopy, and double electron–electron resonance (DEER) spectroscopy have been applied to the study of iGluR dynamics [[43], [44], [45], [46], [47], [48], [49], [50], [51], [52], [53], [54], [55]]. Structural data have also opened the way to several structure-based computations, at various levels of resolution. Many questions regarding the activation and desensitization mechanism of iGluRs, as well as the allosteric events triggered upon ligand binding can now be explored by computational modelling and simulations of structural dynamics [[56], [57], [58]]. We provide below an overview of two major group of computations, molecular dynamics (MD) simulations and elastic network model (ENM) analyses for predicting structure-encoded dynamics, and the emerging mechanisms of function.

Computational methods for modelling protein dynamics

Molecular dynamics (MD) simulations use a detailed force field to describe all the covalent and non-covalent interactions between all atoms in the system including the protein as well as surrounding water and ions (Fig. 2A). This method is based on a large number of time steps involving cycles of calculating all forces on all atoms and numerically integrating Newton’s equations of motion to update their positions and velocities, generating trajectories with full atomic detail [59]. The requirement of calculating the forces on all atoms every few femtoseconds [59] means that such methods are rather slow and can barely access time scales relevant for biology. They are indeed usually limited to hundreds of nanoseconds (ns) or less, especially for large systems such as multimeric receptors in membranes (785,175 atoms for the AMPAR-membrane-water system shown in Fig. 2A, made using MODELLER [60,61] and CHARMM-GUI [62]). The capabilities of MD simulations have somewhat improved with the development of specialized supercomputers [63] and with enhanced sampling schemes [64] that made it possible to calculate free energy landscapes associated with particular conformational changes and binding processes for small systems such as iGluR LBDs [65,66].

Fig. 2

Two main methods for computational analysis of protein dynamics are illustrated. (A) All-atom molecular dynamics (MD) simulations include all atoms from the protein as well as surrounding water, salt and membrane. Shown is another model of a whole GluA2 AMPAR homotetramer including the C-terminal tails, generated using MODELLER and the CHARMM-GUI (785,175 atoms, of which 54,176 belong to the protein). The inset shows a zoomed view of three amino acids that are interacting and the myriad interactions between all their atoms. Covalent bonds are indicated by thick black dashed lines, Coulomb electrostatic interactions are yellow and van der Waals contacts are grey. (B) An elastic network model (ENM) representation is shown for the same AMPAR model with the C-tails removed. It has a node for each residue (a total of 3168 nodes). Uniform springs connect residues within 15 Å (thick black lines). The three residues shown in panel A are now represented by three nodes with springs describing their interactions, which oscillate around the starting distances shown.

An alternative and complementary approach is to use a coarse-grained model, such as elastic network models (ENMs), which enables one to efficiently visualize the large collective modes of motions of the system at the expense of atomic details [67]. The adoption of ENMs in normal mode analysis (NMA) greatly speeds up the calculations, permitting one to evaluate the dynamics of whole iGluRs within minutes (real time) even in the presence of membrane. The most commonly used ENM is the anisotropic network model (ANM) [[68], [69], [70]]. Therein each residue is treated as a bead/node located at the position of the α-carbon and all interactions within a cutoff distance (of 12–15 Å) are modelled as uniform springs (Fig. 2B). The ANM provides a simple analytical expression for the Hessian matrix H of the second derivatives of the potential with respect to deformations. The decomposition of H yields a spectrum of normal modes of motions. The lowest frequency modes among them describe the most accessible collective motions that are usually functional, and occur at biologically relevant (microsecond to millisecond) timescales [71]. ENMs can also be used to explore allosteric signaling [72,73] and to identify important residues such as hinge sites and kinetic hotspots. The latter can be identified using the earlier and even faster Gaussian network model (GNM), which is more accurate at the expense of 3D information [[74], [75], [76], [77]].

Dynamics of iGluR domains revealed by MD simulations

The arrival of the first LBD structures from AMPAR subunit GluA2 in 1998 [78] and 2000 [11] was immediately followed by a number of MD simulation studies on this domain that continue to this day. Early studies were limited by computational power and therefore focused on small systems and short timescales. The main highlights of this period include: 1) simulations lasting 5 ns or less of LBD monomers in the presence and absence of agonists, which found that the open apo form showed larger fluctuations in line with the idea that agonists stabilized closed cleft conformations [79,80]; 2) MD simulations where kinetic energy was applied to specific binding site residues based on quantum mechanical calculations of vibrational energy transfer from agonist [81,82] and this energy was found to propagate through the LBD to dimer interface helix J where it would impact activation and desensitization [83]; and 3) the first enhanced sampling simulations of the LBD using a “partial template forcing” that drove large scale cleft opening and enabled analysis of associated atomic rearrangements and ligand dissociation [84,85].

Since then, studies have increasingly accessed longer timescales and carried out more detailed characterizations of conformational changes. This phase started around 2006 with simulations of the GluA2 LBD [86] and the GluN1 NMDAR LBD [87] that captured cleft closure in 20 ns of conventional MD. These were immediately followed by studies that used umbrella sampling for exploring cleft closure energetics, which found that the apo LBD was preferentially in the open state with a free energy difference of ∼4 kcal/mol with respect to the closed state [65,66]; the antagonist DNQX stabilized the open state; whereas the agonist glutamate stabilizes the closed state [65]. Similar analyses carried out on a larger panel of ligands showed that both binding and cleft closure energetics were important for efficacy [88] and that NMDA receptor LBDs showed a conformational selection mechanism that was not seen for the GluA2 AMPAR LBD [89,90].

Recent increases in computing power have further enabled the simulations of larger systems such as LBD dimers, which provided insights into rearrangements linked to desensitization [91,92], and NTD dimers, which revealed clamshell motions and interprotomer rearrangements akin to LBDs and metabotropic glutamate receptor ligand-binding cores, which may explain their allosteric signaling function in NMDARs [93,94].

Once whole receptor structures became available, simulations were carried out on tetramers including both the LBD and TMD for the GluA2 AMPAR [95] and the GluN1/2B NMDAR heteromer [96,97]. These were performed using targeted MD protocols to drive LBD conformational changes related to activation and desensitization, and the response of the TMD and the connecting linkers was studied. These studies indicated potential mechanisms whereby the linkers enable the LBDs to pull open the TMD channel in a subunit-specific manner.

Most recently, a set of MD simulations each on the order of five microseconds carried out with the Anton supercomputer for AMPAR LBD monomers and dimers revealed the detailed mechanism, energetics and kinetics of glutamate binding to individual LBD clamshells [98].

Dynamics of intact AMPAR and NMDAR revealed by elastic network models

Elastic network model NMA has been applied much more extensively to large assemblies than all-atom MD simulations where the system size becomes an issue. Large conformational changes such as the splaying apart of the NTD dimers of iGluRs observed in low resolution electron microscopy studies [18,33] are beyond the scope of MD simulations; whereas they can be readily modelled by ENMs [94,99] (see Fig. 3, Fig. 4).

Fig. 3

Global dynamics of iGluRs related to allosteric inhibition and desensitization. (A) A conformation generated by traversing ANM mode 4 of the original AMPAR crystal structure (PDB 3KG2) in one direction shows a compact packing of the NTD and LBD, which is highlighted in red and blue for the two pairs of NTD dimers and LBD dimers. The lower diagram displays the relative positions of the NTD dimers, which reveal a compact O shape when viewed from the extracellular region. (B) A similar compaction of the NTD and LBD is observed in the allosterically inhibited NMDAR crystal structure (PDB 4PE5) and a similarly compacted NTD arrangement is seen in the top view. (C-D) Traversing this ANM mode back in the other direction, we pass through the starting structure (C) where the NTD forms the classical N shape, and then reach an alternative conformation with the NTD dimers lifted away from the LBD and tilted away from each other and from the central symmetry axis similar to low resolution desensitized structures (D). See Movie 1 for more details, including cross-sectional views of the NTD, LBD and TMD layers.

Fig. 4

Global dynamics of iGluRs related to activation. (A-B) Comparison of the original AMPAR structure (PDB 3KG2) (A) and ANM intermediate along mode 6 (B) shows an iris-like rotation and opening of the TMD and NTD inter-dimer pivoting coupled to complex LBD dynamics, involving a rolling of the LBD dimers towards and away from each other together with a sliding rotation past each other. Top views of the LBD and TMD are displayed in each panel, along with color-coded arrows to more clearly show the collective motions of the four subunits and three layers (NTD, LBD and TMD) in mode 6. The latter can be viewed more clearly in Movie 2. (C) The TARP-associated open channel state (PBD 5WEO) reveals similar motions of the LBD and TMD but minimal NTD dynamics as a result of the modified NTD-LBD linker present in the crystallography construct used. The TARPs are shown in surface representation to highlight their optimized positioning for modulation of LBD and TMD dynamics. Black surfaces on the TARPs in the top view of the TMD layer indicate contacts with the LBD.

ANM analysis was initially applied to NTD monomers and dimers where clamshell closure and inter-subunit rearrangements were linked to allosteric signaling [26,93,94,100], and to an LBD tetramer [101]. The latter study revealed that the softest mode intrinsically accessible to the LBD tetramer is a cooperative rolling of the two dimers towards and away from each other that occurs during activation [101]. This cooperative motion has since been characterized further in the context of whole receptors both in structural studies [31] and using the ANM [99] as described next.

ANM studies of whole receptors provided a first glimpse of the dynamic coupling between the NTD, LBD and TMD layers [94,99,102]. These analyses further demonstrated the existence of cooperative modes of motions that are shared between AMPARs (exemplified by those sampled by the resting state GluA2 homomer) and NMDARs (exemplified by the ifenprodil-inhibited GluN1/N2B heteromer), as well as subtype-specific differences [94,99]. For example, the most energetically favorable mode (mode 1) is a global bending in both cases, with the NTD bending towards the TMD (or towards the plasma membrane, in a cellular context, into proximity of auxiliary subunits [36,38,103]), facilitated by the unstructured (∼16 aa) NTD-LBD linker [103]. In this mode, the loosely packed NTD-LBD interface acts as the main hinge region in the AMPAR whereas this region is tightly packed in the NMDAR so the LBD-TMD linker region is its main hinge site. Higher frequency modes showed more local motions such as a pivoting at the interface between NTD dimers coupled to rotations of the other domains and the aforementioned rolling apart of the LBD dimers [99] as well as clamshell motions of individual NTDs [94]. Iris-like rotations of the TMD were also seen reminiscent of channel opening. These motions were also seen in the presence of a lipid bilayer (also modelled using an ENM) which suppressed some modes and favored others [99].

A recent study by Popescu and coworkers used ENMs for flexible fitting of electron microscopy data and transition pathway modelling, yielding insights into NMDAR activation [104]. In addition to agreeing with conclusions from the original rigid-body fitting of isolated domains and portions thereof such as an opening of the NTD clefts and separation of the NTD dimers [27], the NMA-based flexible fitting revealed rearrangements of the LBD-TMD linkers consistent with pulling open the pore [104], but the data for the TMD was insufficient to enable resolution of the channel itself. This study further confirmed that two individual normal modes overlap well with the transition from the allosterically inhibited state to the active state [104]: both modes feature a separation of the two NTD dimers coupled to rotations of the other two domain layers akin to modes of motion observed in the earlier ANM study of whole AMPARs and NMDARs [99]. More detailed analysis of the transition with the interpolated ENM approach revealed sequential motion of the domains starting with the NTD relieving inhibition, followed by the NTD-LBD interface then the LBD itself, passed on to the LBD-TMD linkers and finally to the TMD.

Emerging mechanism of function of iGluRs from integration of ENM and MD data

iGluRs have a large extracellular region that facilitates an extensive spectrum of dynamic rearrangements, which are being extensively explored in structural and computational studies. All these observations translate into a selective excitation of pre-existing motions i.e. glutamate binding to the LBD selectively triggers cooperative motions from amongst those soft modes already encoded by the receptor structure. Such intrinsic modes of motion of whole receptors and isolated domains have been visualized using MD simulations [86,87,93,98] and ENM NMA [94,99,104], providing a deeper understanding of the relevance of particular conformational changes to activation, desensitization and allosteric modulation than that provided by the structures alone. Our ANM analysis of the intact AMPAR and NMDAR structures [94,99,102] revealed two cooperative global motions of particular interest, referred to as mode 4 and mode 6. These rearrangements appear to directly drive channel operation and likely influence binding to synaptic interaction partners to indirectly modulate receptor localization and signaling [7,8].

Mode 4 of the AMPAR (Fig. 3 and Movie 1) shows a good example of how the structural changes are triggered between domains. Coupled dynamics of the NTD and LBD induces contraction of the receptor with these domains coming towards each other and the LBD approaching the TMD [99], forming a compact structure that has been trapped in a cryo-EM study and also appears to be an inactive conformation [102]. Similar movements have been observed upon comparison of kainate/(R,R)-2b and ZK200775-bound states of AMPA receptors using ANM modes in a dynamic alignment scheme [105]. Interestingly, the compact structure (Fig. 3A) resembles the allosterically inhibited NMDAR (Fig. 3B), suggesting that this mode enables the structural transition related to entry into and release from an allosterically inhibited state as seen by following the mode in the other direction (see Fig. 3D). Associated with this is a change in the relative positioning of the NTD pairs of dimers (bottom diagrams and Movie 1), which vary from an O-shape organization in the compact state (panel A) through the classical N-shape in the starting structure (panel C) to a splayed apart, stretched shape at the other end of the mode (panel D).

The rearrangement of the NTDs in the compact AMPAR (panel A) approximates that of the inhibited NMDAR (panel B). The compact, inhibited-like state exhibits a close packing at the LBD-TMD interface where closure of the LBD clamshells could have less influence through the coiled LBD-TMD linkers not receiving as much tension. The release transition along this AMPAR ANM mode (Fig. 3C–D) shows the two NTD dimers separating from each other and from the LBD akin to transitions seen by comparing cryo-EM structures of NMDARs [27]. This collective motion also resembles conformational changes linked to desensitization in AMPARs [18,33]. Hence, it provides a possible explanation for the relationship between allosteric inhibition and desensitization in NMDARs [106] that may also hold in AMPARs.

The second global mode, mode 6 of the AMPAR, involves a set of cooperative dynamics, which is likely important for activation (Fig. 4A–B and Movie 2). In this case, we observe coupled inter-dimer pivoting motions in the LBD and NTD layers that are coupled to an iris-like rotation and opening of the TMD channel. These motions of the LBD and NTD have been confirmed by disulfide crosslinking [99,101] and are seen in structural studies that could not resolve channel opening [30,31,33,36,38,101,102]. The NTD motion was not captured in recent cryo-EM structures that resolved an open state of AMPAR in the presence of TARPs, possibly because one study used a truncated NTD-LBD linker (that aided crystallization by restricting NTD dynamics) [38], and the other masked out the NTD to improve resolution of the LBD and TMD sectors [36]. The rest of the motion closely approximates the structural change observed in the transition to the active structure (considering only the LBD and TMD). Thus, the ANM analysis allows us to infer a mechanism for activation that incorporates motions of all the domains that have been seen in separate structural studies: glutamate binding favors closure of individual LBDs [65,88], which pulls on inter-LBD interfaces, resulting in rearrangements such as a rolling of the dimers towards and away from each other pivoted about interfacial helix G, and a sliding of these helices past each other, that increase the pulling force on the other domains resulting in closure of the NTD layer and opening of the TMD channel pore. For a visualization of these cooperative events see Movie 2, where the motions of the intact AMPAR (side view) and those of the individual domains (cross-sectional views) can be clearly seen. The LBD D1 dimer interface is fairly stable on the microsecond timescale [92], allowing the above activation-related inter-dimer rearrangements to occur before the LBD dimer interface breaks and desensitization sets in.

Auxiliary subunits are optimally placed in close proximity to the LBD and TMD (Fig. 4C) to modulate this behavior and hence channel gating and ion conduction. For example, TARPs are known to favor channel opening and stabilize the activated state [8,36,38] whereas GSG1L disfavors channel opening and stabilizes inactive/desensitized states [36,37]. These functions may be achieved through the auxiliary subunits acting to constrain or amplify the intrinsically accessible cooperative motions related to TMD channel opening and desensitization.

There are still a number of unanswered questions that are likely to be addressed over the next few years through a combination of experimental and computational approaches. While key states and transitions can now be visualized, we are still lacking the energy landscape and a complete structural view of the gating cycle. It is also still unclear how auxiliary subunits and allosteric ligands modulate iGluRs. Recent and continuing developments in cryo-EM as well as hybrid simulation methods that benefit from the advantages of both NMA and MD are expected to shed light to these pending questions.

Conclusion

iGluRs undergo complex dynamics related to their functions, which are starting to be understood using experimental and computational studies. Two main types of computational methods have been applied that have their own strengths and limitations. All-atom MD simulations use detailed force fields that allow for analyzing and generating free energy landscapes but this comes at great computational cost that prevents the examination of cooperative events in large structures such as the multimeric iGluRs and the attainment of biologically relevant time scales. Coarse-grained representations such as ENMs, on the other hand, enable identification of biologically relevant global (soft) modes of motion that often agree with and build upon experimental data from structural studies, albeit at low resolution. In particular, these intrinsically accessible and energetically favorable modes of motion, analyzed in the context of existing experimental (structural) data and atomic-simulations performed for substructures, reveal well-defined couplings between the motions of AMPAR domains, and provide better glimpses of channel activation, desensitization and inhibition mechanisms than static structures alone. An emerging picture is an allosteric communication between the NTD, LBD and TMD modulated by glutamate binding to the LBD, which selectively triggers pre-existing robust modes uniquely encoded by the overall tetrameric architecture.

Declarations of interest

None.