Next Article in Journal
Nematicidal Coumarins from Cnidium monnieri Fruits and Angelica dahurica Roots and Their Physiological Effect on Pine Wood Nematode (Bursaphelenchus xylophilus)
Next Article in Special Issue
DeepBindGCN: Integrating Molecular Vector Representation with Graph Convolutional Neural Networks for Protein–Ligand Interaction Prediction
Previous Article in Journal
Chitosan, Chitosan/IgG-Loaded, and N-Trimethyl Chitosan Chloride Nanoparticles as Potential Adjuvant and Carrier-Delivery Systems
Previous Article in Special Issue
In Silico Methods for Identification of Potential Active Sites of Therapeutic Targets
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Molecules
Volume 28
Issue 10
10.3390/molecules28104108
molecules-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Binding and Dynamics Demonstrate the Destabilization of Ligand Binding for the S688Y Mutation in the NMDA Receptor GluN1 Subunit

by
Jake Zheng Chen
1,2,
William Bret Church
1,
Karine Bastard
1,
Anthony P. Duff
3 and
Thomas Balle
1,2,*
1
Sydney Pharmacy School, Faculty of Medicine and Health, The University of Sydney, Camperdown, NSW 2006, Australia
2
Brain and Mind Centre, The University of Sydney, Camperdown, NSW 2050, Australia
3
National Deuteration Facility, Australian Nuclear Science and Technology Organization, New Illawarra Road, Lucas Heights, NSW 2234, Australia
*
Author to whom correspondence should be addressed.
Molecules 2023, 28(10), 4108; https://doi.org/10.3390/molecules28104108
Submission received: 5 April 2023 / Revised: 10 May 2023 / Accepted: 11 May 2023 / Published: 15 May 2023
(This article belongs to the Special Issue Role of Computer Aided Drug Design in Drug Development)
Browse Figures
Versions Notes

Abstract

:
Encephalopathies are brain dysfunctions that lead to cognitive, sensory, and motor development impairments. Recently, the identification of several mutations within the N-methyl-D-aspartate receptor (NMDAR) have been identified as significant in the etiology of this group of conditions. However, a complete understanding of the underlying molecular mechanism and changes to the receptor due to these mutations has been elusive. We studied the molecular mechanisms by which one of the first mutations within the NMDAR GluN1 ligand binding domain, Ser688Tyr, causes encephalopathies. We performed molecular docking, randomly seeded molecular dynamics simulations, and binding free energy calculations to determine the behavior of the two major co-agonists: glycine and D-serine, in both the wild-type and S688Y receptors. We observed that the Ser688Tyr mutation leads to the instability of both ligands within the ligand binding site due to structural changes associated with the mutation. The binding free energy for both ligands was significantly more unfavorable in the mutated receptor. These results explain previously observed in vitro electrophysiological data and provide detailed aspects of ligand association and its effects on receptor activity. Our study provides valuable insight into the consequences of mutations within the NMDAR GluN1 ligand binding domain.

1. Introduction

Early infantile epileptic encephalopathies are some of the most devastating early-onset epilepsies and can lead to cognitive, sensory, and motor development impairments. The etiology of these syndromes is wide-ranging and can be structural, metabolic, neurological, or genetic in nature [1]. Recently, due to advances in high throughput sequencing technologies, there has been an increase in the identification of genetic mutations associated with early infantile epileptic encephalopathies, including mutations in the N-methyl-D-aspartate receptor (NMDAR) [2].
NMDARs are ionotropic glutamate receptors (iGluR) and consist of four subunits, which are typically two GluN1 subunits and two GluN2 subunits (Figure 1A); however, the existence of other triheterotetrameric arrangements has been reported [3,4,5]. They are unique in their mode of activation, requiring simultaneous binding of glutamate at the GluN2 subunit and a co-agonist such as a glycine [6], D-serine [7], D-alanine [6,8], or D-cycloserine [7] at the GluN1 subunits [9,10,11], along with membrane depolarisation to displace a Mg2+ ion block present in the channel [12,13]. Together, these events activate the receptor and enable Ca2+ influx and further downstream effects within the cell [14].
Structurally, all GluN subunits consist of an amino terminal domain (ATD), a ligand binding domain (LBD), a transmembrane domain (TMD), and an intracellular carboxyl-terminal domain (CTD) (Figure 1A). X-ray structures of the NMDAR GluN1 LBD were originally reported in 2003 [7], and advances in cryogenic electron microscopy have since led to the membrane and extracellular receptor structures being solved [15,16,17,18,19]. Significantly aided by crystallographic work [20], the individual domains appear largely physically discreet. The NMDAR is normally thought to be involved in learning and memory formation [21,22]. However, dysfunction of the NMDAR has been associated with complex, multifactorial conditions such as schizophrenia [23,24,25,26] and Alzheimer’s disease [27,28]. There have been several studies conducted which utilize molecular dynamics simulations to study the function and effects of the NMDAR [29,30,31,32,33].
Generally, mutations within the GluN1 subunit are associated with conditions such as hypotonia [34], hyperkinetic and stereotyped movements, as well as other conditions [35,36]. Predominantly, mutations in the GluN1 subunit have been identified within the TMD [33,34,35,36,37,38,39,40,41,42,43]; however, mutations within the ATD have also been reported [36]. The first mutation within the LBD linked to infantile epileptic encephalopathies, Ser688Tyr (S688Y), was described by Zehavi et al., in 2017 [34]. Skrenkova et al. [33] established the effects of the S688Y mutation in terms of a significant decrease in functional potency for both glycine and D-serine along with a differential regulation in surface delivery of the mutant NMDAR in HEK293 cells.
However, the molecular mechanisms by which the S688Y mutation leads to decreased potency are yet to be established. We utilize computational methods, including a combination of docking and molecular dynamics simulations, to assist in understanding the ligand dynamics within the binding site. The binding thermodynamics relevant to protein–ligand binding was also estimated in the study. Computer-aided drug design (CADD) utilizes these techniques to determine the function and structure of the protein as a part of structure-based drug design (SBDD) and alchemical methods to determine ligand binding affinity as a part of ligand-based drug design (LBDD).
Our results reveal how the receptor functions and, specifically, how the S688Y mutation leads to impaired receptor function. These calculations demonstrate changes in the receptor–ligand complex structure with important concomitant alterations of the electrostatic and water-bridge network present within the binding site. This gives a deeper insight into the effects of mutations within the NMDAR, specifically those crucial within the LBD and for ligand binding.
Molecules 28 04108 g001 550
Figure 1. GluN1/GluN2B receptor (PDB: 6CNA) and GluN1 LBD (PDB: 4NF8). (A) Full NMDA receptor GluN1/GluN2B assembly in a phospholipid membrane. GluN1 is shown in salmon, and GluN2B is shown in green. (B,C) LBD of GluN1 shown in cartoon representation. Glycine is shown in space fill representation in cyan, and Ser688 is shown in green. Note that the CTD is not visible in the electron microscopy structure shown. Figure generated using PyMOL [44].
Figure 1. GluN1/GluN2B receptor (PDB: 6CNA) and GluN1 LBD (PDB: 4NF8). (A) Full NMDA receptor GluN1/GluN2B assembly in a phospholipid membrane. GluN1 is shown in salmon, and GluN2B is shown in green. (B,C) LBD of GluN1 shown in cartoon representation. Glycine is shown in space fill representation in cyan, and Ser688 is shown in green. Note that the CTD is not visible in the electron microscopy structure shown. Figure generated using PyMOL [44].
Molecules 28 04108 g001

2. Results

2.1. Ligand Docking to GluN1 LBD

Docking of both glycine and D-serine in the wild-type GluN1 LBD (Figure 2A,B and Figure 3A,B) utilizing a flexible receptor-docking algorithm revealed key salt bridge interactions between the carboxylate group of both ligands and the guanidium group of Arg523 of the protein as well as between the ammonium groups and the side chain carboxylate group of Asp732. Further interactions involved π-cation interactions with the phenyl group of Phe484, and D-serine displayed π-cation interactions with Trp731. Hydrogen bonding was also present in the docked conformations with the backbone of residues Pro516, Thr518, Ser688 and side chains of residues Arg523 and Trp732. These interactions are comparable with those observed in the X-ray structures of the NMDAR GluN1 LBD co-crystallized with glycine (PDB ID: 1PB7) and D-serine (PDB ID: 1PB8) by Furukawa and Gouaux [7]. This demonstrates that our docking protocol is valid and appropriate for studying the S688Y mutant. The root mean square deviation (RMSD) of all ligand-heavy atoms for the docked pose when compared to the X-ray structure was 0.8 Å (PDB: 1PB7) for glycine and 1.7 Å for D-serine when superimposed on protein heavy atoms.
In the docked poses, the carboxylate and ammonium groups of glycine and D-serine were present in a polar pocket of the binding site, with the hydroxymethyl group of D-serine present in a hydrophobic sub-pocket. However, in the S688Y mutant receptor, the environment is changed with fewer polar residues present in the environment of both glycine and D-serine due to the alterations of the binding pocket induced because of the change to tyrosine. Furthermore, the key Arg523 residue is not present in the binding pocket due to steric hindrance from the bulky Tyr688 hydroxyphenyl side chain (Figure 3C,D), which upon docking leads to an associated loss of a salt bridge and hydrogen bond in comparison with the wild-type configuration with serine (Figure 2C,D). Other interactions, however, remained intact. Of interest is the main chain nitrogen atom of Tyr688, which continues to facilitate hydrogen bonding with the carboxylate group of glycine and D-serine (Figure 3C,D) with the formation of an additional aromatic hydrogen bond between the side chain of Tyr688 and glycine. Positioning of the carboxylate group of both ligands is similar in the wild-type protein; however, due to the absence of the guanidium group of Arg523 inside the binding pocket in the mutant receptor, there is a conformational difference between glycine and D-serine whereby the carboxylate group on glycine has a “downward” rotated orientation relative to the carboxylate group in D-serine. Meanwhile, the ammonium group is held in place by Asp732 and is comparable between both the wild-type protein and the S688Y mutant (Figure 2 and Figure 3).
The analysis demonstrated that there was only one binding pose for glycine with a docking score of −5.28 kcal/mol, whereas there were two poses present for D-serine with a top Glide XP score of −7.37 kcal/mol. However, for the S688Y GluN1 LBD, there were four possible docked conformations for glycine, with a top docking score of −4.10 kcal/mol. There were nine output poses for D-serine with a top docking score of −3.91 kcal/mol. In both ligands, the top S688Y docking score was less negative compared to the wild-type, which can be indicative of poor ligand affinity in the receptor site.

2.2. Microsecond MD Analysis

To assess the stability of the systems and their behavior, we conducted 1 μs MD simulations on the apo state, as well as on top, binding poses of glycine and D-serine in wild-type and S688Y mutant receptors. Overall, the 1 μs MD simulations demonstrated that the simulation systems consisting of the solvated GluN1 LBD with glycine and D-serine ligands were stable.
The apo simulations performed demonstrated no significant fluctuations with protein RMSD fluctuating and stabilizing at approximately 4.0 Å from the input structure (Supplementary Figure S2). Overall protein RMSD in the wild-type glycine and D-serine docked poses stabilized at 3.7 Å for glycine and around 4.0 Å for D-serine as compared to the input structure (Figure 4A,B). RMSD of residues within 5 Å of the binding pocket stabilized at 2.5 Å for glycine and 3.2 Å for D-serine in the wild-type receptor pocket. The binding pocket RMSD for the S688Y receptor increases to 4.1 Å for glycine and 4.3 Å for D-serine. Furthermore, flexible loop regions within the protein are similar compared to those observed in the crystal structure [46] (Supplementary Figure S2).
Glycine ligand RMSD values, when compared against the wild-type protein, hover approximately at 1 Å while D-serine ligand RMSD in the analogous comparison is approximately 2.6 Å (Figure 4A,B). In the case of the S688Y mutation, glycine RMSD, when fitted against the protein, is approximately 3 Å, while D-serine has a ligand RMSD fitted against the protein of approximately 2.4 Å (Figure 4C,D). D-serine demonstrated two binding modes with the wild-type receptor, with contact being present with Gln686 in some frames of the trajectory (Supplementary Figure S3). Ligand–receptor interactions during the entire simulation for all four liganded simulation systems remained comparable with our observations from the docking experiments.
Individual residue root–mean–square fluctuation (RMSF) comparison between the wild-type and the S688Y mutation LBD shows that per-residue fluctuation increased by 27% for the glycine-bound GluN1 receptor and 10% for the D-serine-bound GluN1 receptor, respectively (Figure 5).
Solvent accessibility of the ligand in the wild-type LBD during the entire course of the simulation was primarily limited to the ammonium group of glycine and D-serine. In the case of D-serine, the side chain hydroxyl group was also solvent accessible for 5% of the simulation time. Observations of a higher SASA of the ligand are indicative of the receptor providing a complex state with less solvent exposure. The solvent-accessible surface area (SASA) (Figure 6) for glycine bound to wild-type protein reached a peak of ~12 Å2 and stabilized in the range of 0–4 Å2 for most of the simulation while the SASA for D-serine bound to the wild-type protein fluctuated between ~6 Å2 and a more solvent accessible conformation of ~12 Å2.
Analysis of the trajectories of the S688Y receptor showed that there was extensive water accessibility of the ligand compared to the wild-type protein, with solvent contact of the ligand reaching as high as 71% of the simulation time. The SASA for glycine bound to the S688Y LBD demonstrated stability around 4 Å2 until ~450 ns, where the loss of interactions with Thr518 and Asp732 (Figure 6 and Figure 7) led to a large increase in SASA reaching as high as 45 Å2. For D-serine bound to the S688Y LBD, there was an initial period of high solvent accessibility; however, SASA stabilized around 10–15 Å2. There were no large changes in the SASA due to changes in the Thr518 and Asp732 interactions, as they remained stable for the duration of the simulation.
Ligand contact with the guanidinium side chain of the key residue Arg523 was present in 99% of frames in the wild-type receptor. Ligand contact with Ser688 revealed a preference for the main chain nitrogen but not the side chain as seen in the docking poses with hydrogen bonds with the main chain nitrogen present for 76% of simulation time for glycine and 94% in D-serine while hydrogen bonds for the side chain were at 15% and 27% of simulation time for glycine and D-serine, respectively. Ligand interactions with the guanidium side chain of Arg523 were lost upon mutation for glycine; however, D-serine was able to maintain hydrogen bonds with the guanidium side chain of Arg523 in the S688Y GluN1 LBD.

2.3. Molecular Dynamics Simulations

To determine the reproducibility of the simulations and to potentially observe any other plausible protein–ligand interactions, we performed five 200 ns simulations that were randomly seeded to have alternative starting velocities and one that used the identical seed as in the microsecond production runs for each of the four liganded simulation systems (Supplementary Table S1). When using the same seed used in the microsecond simulations, the 200 ns simulation gave results identical to the first 200 ns of the longest production run. Compared to the microsecond simulations, the simulations with randomized initial velocities for both glycine and D-serine bound to the wild-type protein were comparable and had similar protein-ligand contacts for all simulations. Solvent contacts were present in three out of the five simulations for glycine and present in all simulations for D-serine. In the wild-type protein, all interactions observed were similar to those observed in the microsecond simulations. SASA for both glycine and D-serine stabilized around 1 Å2 with a maximum of 10 Å2 (Figure 8 and Figure 9).
In the S688Y mutated protein, solvent exposure was observed in all simulations with both glycine and D-serine. Both glycine and D-serine bound to the S688Y mutated protein experienced an average SASA in the range of 15–30 Å2 for all seeds, contrasting with observations of the wild-type protein, which showed SASA in the range of 4 Å for glycine and 7 Å for D-serine. Interestingly, seed 1293 of glycine with S688Y mutated protein demonstrated a complete dissociation of glycine from the receptor binding site.

2.4. Binding Free Energy Calculations

To estimate the binding free energy of glycine and D-serine, we performed MM-GBSA free energy calculations to approximate the ΔGbind to both the wild-type and S688Y mutant receptors for an estimation of the relative overall binding affinity of the ligand to the receptor and the individual components (Table 1).
Compared to the wild-type receptor, the calculated binding free energy in the S688Y mutant receptor was increased from −30.0 kcal/mol to −26.3 kcal/mol for glycine and −34.6 kcal/mol to −22.7 kcal/mol for D-serine. This result is directly comparable to previous electrophysiology studies on the S688Y mutant, where the EC50 of the ligands was increased significantly compared to the native receptor [33]. The contribution from electrostatic interactions (ΔGbind coulomb) increased from −5.62 kcal/mol to 11.8 kcal/mol for glycine and −6.62 kcal/mol to 7.58 kcal/mol for D-serine. This strongly supports that it is no longer energetically favorable for the system to form salt bridges which are consistent with our results from the MD trajectories where electrostatic interactions between the ligand and Arg523 were present in far fewer frames in the S688Y protein compared to the wild-type, especially for the glycine ligand where Arg523 interactions were almost completely absent. There were minor changes to other binding affinities, which included an increased contribution from hydrogen bonding (ΔGbind hydrogen bonding) in both glycine and D-serine post-mutation, which is reflective of our docking results. Lipophilic interactions (ΔGbind lipophilicity) decreased slightly for glycine, which reflects our previous observations where the ligand environment is slightly more hydrophobic (Figure 2A,C); however, this did not change significantly for D-serine.
Solvent interactions (ΔGbind solvent) became more favorable, changing from −7.04 kcal/mol to −18.9 kcal/mol for glycine and similarly a drop from −6.18 kcal/mol to −9.65 kcal/mol for D-serine, which reflects our observations from the simulations where there was a significant increase in solvent exposure in the mutated protein. This greater extent of ΔGbind solvent decrease for glycine is also consistent with our observations where complexes containing glycine had a higher number of water bridges in the S688Y receptor compared to D-serine, which had a modest increase. Van der Waals interactions (ΔGbind van der Waals) also decreased for glycine from −9.97 kcal/mol to −13.2 kcal/mol and −13.8 kcal/mol to −14.9 kcal/mol for D-serine. This may be due to the change in the environment around the glycine to one which is more hydrophobic.

3. Discussion

In the microsecond simulations of the NMDAR ligand binding domain and ligand, there was a demonstrable increase in the presence of water bridges interacting with both glycine and D-serine in the mutated receptor containing the S688Y mutation compared to the wild-type. A group of five water molecules, normally in contact with the ligand [7], is involved in binding for ligands at the GluN1 LBD. This involvement is observed in our microsecond simulations of wild-type protein with both glycine and D-serine (Figure 7A,B). However, these water molecules do not appear to contribute significantly to the wild-type binding network compared to the S688Y mutant. D-serine further demonstrates two possible poses it could use to bind to the wild-type receptor (Figure 3), including an extra contact with Gln686 in the second pose. In the mutated receptor, the water in contact with the ligand increased (Figure 7C,D) relative to that in the wild-type receptor. There was also instability observed with the SASA of the ligands, where values hovered around 5–10 Å for both glycine and D-serine in the wild-type receptor but increased significantly in the mutated receptor. This is important to note as the bulk water in contact water in contact with the ligand can behave differently from structured, isolated water within the binding site [49] and instead contribute to the instability of the system. Furthermore, large protein RMSD differences were observed especially in the case of mutated receptors, reaching 8 Å (Figure 4C,D) compared to those in the wild-type receptor, which fluctuated around 4.5 Å (Figure 4A,B). The present study focuses on the simulation of the LBD in isolation and the interactions of the ligand with the LBD. This approach is appropriate for the present implementation, as shown by Wang et al., (2021) [47], where the LBD and transmembrane domain (TMD) are linked, while each of the subunits remains independent and is able to maintain their structure throughout the transition process.
Loss and weakening of key interactions are also observed in the simulations, which contribute to an increase of ΔGbind and especially the key interaction of the guanidium group Arg523 with the ligand. This interaction is absent in the S688Y receptor with glycine in the microsecond simulation and present in fewer than 20% of frames in the random-seeded simulations. The loss of this key interaction, combined with the change in the binding pocket, contributed to the increase in ΔGbind for both glycine and D-serine. In the current study, we performed the MM-GBSA calculations on the final 50 ns of the stable 200 ns trajectories; however, it would be of interest for future studies to determine whether a calculation over a longer trajectory yields similar results to the findings presented here.
Results from the ΔGbind calculations show that there is a significant decrease in binding affinity for the S688Y mutant LBD compared to that of the wild-type protein with an increase of 11.9 kcal/mol in ΔGbind for D-serine and an increase of 3.75 kcal/mol in ΔGbind for glycine. This is of particular importance due to GluN1/GluN2A’s higher affinity for D-serine compared to glycine in the forebrain regions [48] and that synaptic NMDARs are gated by D-serine [46,50,51,52]. ΔGbind coulomb values increase to 11.83 kcal/mol and 7.58 kcal/mol for glycine and D-serine, respectively, upon mutation, showing that the interference of Tyr688 between the guanidium group of Arg523 and the ligand was able to disrupt this key interaction present in the wild-type NMDA GluN1 LBD. However, the rotamer state for Tyr688 in this study differs from those previously reported [33]. This difference leads to a change in conformation for other residues in the receptor binding site and is significant due to the interference with Arg523. Though we have run substantive molecular dynamics calculations, which should sample favorable conformations, a further study to confirm the extent of any biases from the initial side chain conformation of the receptor model should also be performed. Especially for Tyr688, due to the multiple different conformations, the residue side chain could adapt and impact the surrounding receptor environment. Due to the small size of the glycine and D-serine ligands, this observed difference is potentially very significant for understanding the binding of the ligand. MM-GBSA calculations are suitable for the present scoping study as an initial starting point for the understanding of ligand binding affinity of the endogenous glycine and D-serine ligands to the NMDA GluN1 LBD. Implicit solvent-based approximation methods have been shown to perform well [53] compared to explicit solvent models such as FEP. Future work will employ the implementation of free energy perturbation (FEP) [54,55], which is significantly more complex but will result in ligand binding calculations with much greater confidence.
Traditionally, MD simulations are carried out on the nanosecond scale; however, we have performed microsecond simulations which gives us confidence that the possible ligand poses are well sampled, and this, combined with randomly seeded 200 ns simulations, helps to overcome issues of interpretation from reproducibility of MD simulations and has ensured that the initial system is stable under multiple conditions. Therefore, we are confident in our approach of using the refined docked structure as a starting point for the simulation. The use of the OPLS3e force field is widespread in docking [56,57]; however, due to limitations regarding local hardware and software implementations, the OPLS_2005 force field was used for the current protocol. Future protocols will take advantage of the latest OPLS4 force field. Our choice of the SPC water model is appropriate, despite the potential advantages of other models [58,59]. Furthermore, future studies should integrate simulations of the full NMDAR combined with in vitro analysis.
Our docked poses agree with previously published X-ray structures and show that this approach is suitable for studying the detailed changes which may occur between the wild-type and mutated receptor. The use of the XP scoring function allows the elimination of false positive docking poses to ensure confidence in the generated poses. Furthermore, differences in docking scores between the wild-type receptor and the S688Y receptor reported by the XP scoring function correlate with the calculated free energy of binding from the MD simulation trajectories. The current study looks at one of the currently known mutations within the NMDAR GluN1 LBD; however, other possible mutations within the LBD and in regions in proximity could influence the binding of ligands to this receptor in other ways. Therefore, future studies should aim to ascertain the effects these other mutations have on the LBD and ligand affinity to the binding site. Docking was selected over insertion for ligand placement as binding may not be identical in the S688Y mutant LBD compared to the wild-type, and the use of docking for both parts ensures that there is consistency and comparison between the two systems. As shown by our docking results, the output from glycine and D-serine docking to the wild-type receptor is comparable to those observed in the X-ray crystallography structures as determined by Furukawa and Gouaux [7].
From our results, the S688Y mutation leads to poor binding of both glycine and D-serine to the NMDAR GluN1 subunit. This may result in partial or non-activation of the receptor, which results in issues during early development, which may contribute to an explanation on the molecular level of the encephalopathies as described by Zehavi et al. [34] and further characterized by Skrenkova et al. [33] in vitro where the mutation of the receptor resulted in a drop in EC50 of both D-serine and glycine and also reduced cell surface trafficking of the NMDAR. This manuscript outlines the two primary endogenous ligands for the NMDAR GluN1 subunit, D-serine, and glycine, as characterized by Skrenkova et al. [33]. However, other endogenous ligands, such as L-aspartate, and antagonists, such as 5,7-dichlorokynurenic acid, could also be studied in future studies.

4. Materials and Methods

Figures in this paper were prepared using PyMOL [44] and the Maestro visualization interface [45]. Calculations were performed on Intel Xeon E5-2680 V3 CPUs and Nvidia V100 SXM2 GPUs. Graphs were generated using the Simulation Interactions Diagram module in Desmond 2018-4 [45,60] and GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.

4.1. Receptor and Ligand Preparation

Criteria for model selection include the absence of allosteric modulators, the absence of stabilizing mutations or linkers, high-resolution structures, and the absence of large inhibitors, which may influence the structure of the LBD directly. Therefore, the crystal structure of the GluN1/GluN2A ligand-binding domain in complex with glycine and glutamate (PDB: 4NF8) [61] was selected due to its high resolution of 1.86 Å and absence of allosteric modulators or large inhibitors.
The structure was downloaded and prepared using the protein preparation module in Schrödinger Maestro [45,62]. In brief, this involved assignment of bond orders, adding missing hydrogens, and connection of disulfide bonds. Protonation states were assigned corresponding to pH 7.0 ± 2.0 using PROPKA [63,64]. All water molecules in the PDB entry were removed, and the protein structure was minimized with positional restraints as per default settings using the OPLS3e force field [65].
The S688Y mutation was introduced at the GluN1 binding site using the Maestro residue mutation interface and the rotamer with the highest standard population and probability based on their respective dihedral angles determined using the rotamer selection tool. This corresponds to a chi1 of −66° and chi2 of 98°. Minimization was subsequently carried out with position restraints on Cα atoms using the VSGB solvation model with the OPLS3e force field [65].

4.2. Molecular Docking

Glycine and D-serine were docked into the wild-type and mutant GluN1 binding pockets in their fully ionized form using the Glide docking module [45,66,67,68] with flexible ligand sampling and the extra precision (XP) scoring function. Van der Waals radii were scaled by a factor of 0.80, and a partial charge cut-off of 0.15 eV for non-polar ligand atoms was applied. The receptor grid to guide the docking (Supplementary Figure S1) was defined using the co-crystallized glycine from the original PDB structure by Furukawa and Gouaux [7] as the centroid of the receptor site and therefore, the receptor grid for docking with a boundary of 20 Å × 20 Å × 20 Å. The van der Waals, radius scale factor, was set to 1.0, and a partial charge cut-off to 0.25 eV. The maximum number of conformations and output per ligand was set to 10. Glide utilizes a relative docking energy for each docking pose, with lower energies being more favorable. Docking was carried out using a flexible receptor sampling algorithm to ensure receptor residue conformations were sampled.

4.3. Molecular Dynamics (MD) Simulation

The top binding pose for each of the ligand/receptor combinations was subjected to molecular dynamics (MD) simulations for further analysis. Simulation systems were built using the System Builder tool in Desmond 2018-4 [45,60] using the OPLS_2005 force field, suitable for local hardware [69], and a cubic system with a 10 Å × 10 Å × 10 Å buffer around the protein. Each system was solvated using the explicit SPC water model [59,70]. The overall charge of the system was neutralized by substituting water molecules with Cl− and further NaCl was added to a concentration of 0.15 M by random replacement of water molecules. Each system was subject to the default relaxation protocol, which involved the relaxation of the simulation system using default NVT and NPT ensembles; for more detail regarding the minimization and equilibration code, please see Supplementary Code S1.
Simulations were conducted for each system using identical starting parameters, differing only in the seed values of the randomized initial starting velocities. Simulations were performed for a period of 1 μs with a sampling interval of 20 ps and a default seed of “2007”. Subsequent simulations were 200 ns in duration with a sampling interval of 5 ps using the default seed of “2007” and five randomized seeds for the calculation of initial velocities for the MD simulations (Supplementary Table S1). Simulation trajectories were subsequently analyzed using the simulation interaction analysis interface of Desmond [47,50] to determine the RMSD, RMSF, and SASA values for the receptor and ligand of interest.

4.4. Binding Free Energy Calculations

The final 50 ns from the 200 ns MD trajectories (fixed seed “2007”) were used to estimate the affinity of glycine and D-serine for the GluN1 and mutant receptor using the generalized Born and surface area (MM-GBSA) [71,72] method. The Gibbs binding free energy (ΔGbind) was computed based on the equation shown below:
Δ G b i n d = G c o m p l e x − G p r o t e i n − G l i g a n d
ΔGbind is a measure of estimating the binding affinity of a ligand to a protein, with negative values corresponding to energetically favorable reactions and a positive value being energetically unfavorable.
The MM-GBSA calculations were performed using the thermal_mmgbsa.py plugin of Prime [73,74] as supplied [45]. Calculations were performed using the OPLS3e [65] force field and the polarisable variable-dielectric Generalised Born (VSGB) 2.1 implicit solvation model, which has been thoroughly validated [75]. In brief, the script extracts the final 50 ns and performs geometry optimization to extract energies of the complex, the protein, and the ligand to calculate ΔGbind. MM-GBSA is a method for predicting the binding energy of protein–ligand complexes. This method is based on several approximations, including the molecular mechanics (MM) approximation for describing the atomic interactions, the Generalised Born (GB) approximation for modeling solvent effects, and the surface area (SA) approximation for calculating the non-polar contribution to the binding energy. Despite these approximations, MM-GBSA is a valuable tool for rational drug design that provides a good balance between computational cost and accuracy.

5. Conclusions

Early infantile epileptic encephalopathies are a debilitating group of conditions that have multiple causes. The existence of mutations within the NMDAR GluN1 subunit is well documented in the literature; however, mutations within the LBD have only been discovered recently, and a complete understanding of the impacts of these mutations on receptor function is yet to be achieved. Our current research demonstrates how computational drug discovery tools, including ligand docking and extended molecular dynamics simulations, can be used to bridge this gap. Our results show that the S688Y mutation of the NMDA GluN1 receptor is detrimental to ligand binding and affinity, particularly for the binding of D-serine. This may provide an explanation for some of the clinical observations [2]. In future studies, simulation of the full-length receptor similar to those performed by Černý et al. [30] and Sinitskiy and Pande [32], as well as experimental determination of the mutated S688Y receptor structure would also be welcome, as would other mutant NMDAR structures, to assist the study of the available ligand binding modes, and ultimately the design of agents that can modulate the activity of mutant receptors to aid personalized treatments.
The present study presents an insight into the impact of the S688Y mutation on NMDAR GluN1 subunit function and associated endogenous ligand affinity, particularly glycine and D-serine. Clinically, this study contributes to one aspect regarding the underlying molecular mechanisms, which may lead to symptoms associated with severe early infantile encephalopathy and associated symptoms such as intellectual disability [34,38], epilepsy [35], and postnatal microcephaly [76]. Further understanding of the implications of mutations within the LBD of the NMDA GluN1 subunit may aid in the development of personalized treatments based on patient genotyping and may be picked up during genetic screening during pregnancy [77].

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/molecules28104108/s1, Table S1: Seed values of randomly seeded MD simulations. Figure S1: Position of the ligand docking grid; Figure S2: Trajectory analysis for apo state NMDA ligand binding domain during 1 μs molecular dynamics simulation; Figure S3: Ligand–protein interaction timeline during the 1000 ns (1 μs) simulation; Code S1: Relaxation protocol for Desmond.

Author Contributions

Conceptualization, J.Z.C., T.B. and W.B.C.; methodology, J.Z.C., T.B. and W.B.C.; validation, J.Z.C.; formal analysis, J.Z.C.; investigation, J.Z.C.; data curation, J.Z.C.; writing—original draft preparation, J.Z.C.; writing—review and editing, J.Z.C., T.B., K.B., A.P.D. and W.B.C.; visualization, J.Z.C.; supervision, W.B.C. and T.B; project administration, T.B.; funding acquisition, T.B. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the National Computational Merit Allocation Scheme (NCMAS) 2022 (Grant number NCMAS 2022-154) and the Sydney Informatics Hub HPC Access Scheme 2022 and 2023. J.Z.C is in receipt of an Australian Government Research Training Program (RTP) Scholarship.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Simulation data and other materials will be provided upon request.

Acknowledgments

This project was undertaken with the assistance of resources and services from the National Computational Infrastructure (NCI), which is supported by the Australian Government. The authors would also like to acknowledge the technical assistance and computing resources provided by the Sydney Informatics Hub, a Core Research Facility at The University of Sydney. We would also like to acknowledge Ali Kusay for his maintenance of the Balle MD cluster at the Brain and Mind Centre at The University of Sydney. The National Deuteration Facility is partly supported by the National Collaborative Research Infrastructure Strategy—an initiative of the Australian Government.

Conflicts of Interest

W.B.C. declares he is a shareholder of Filamon Pty. Ltd., which owns intellectual property relating to the use of small peptides and other compounds in the diagnosis and treatment of human disease. Other authors declare no conflict of interest.

Sample Availability

Samples of the compounds are not available from the authors.

References

  1. Hwang, S.-K.; Kwon, S. Early-onset epileptic encephalopathies and the diagnostic approach to underlying causes. Korean J. Pediatr. 2015, 58, 407–414. [Google Scholar] [CrossRef] [PubMed]
  2. Allen, A.S.; Berkovic, S.F.; Cossette, P.; Delanty, N.; Dlugos, D.; Eichler, E.E.; Epstein, M.P.; Glauser, T.; Goldstein, D.B.; Han, Y.; et al. De novo mutations in epileptic encephalopathies. Nature 2013, 501, 217–221. [Google Scholar] [CrossRef] [PubMed]
  3. Stroebel, D.; Casado, M.; Paoletti, P. Triheteromeric NMDA receptors: From structure to synaptic physiology. Curr. Opin. Physiol. 2018, 2, 1–12. [Google Scholar] [CrossRef] [PubMed]
  4. Traynelis, S.F.; Wollmuth, L.P.; Mcbain, C.J.; Menniti, F.S.; Vance, K.M.; Ogden, K.K.; Hansen, K.B.; Yuan, H.; Myers, S.J.; Dingledine, R. Glutamate Receptor Ion Channels: Structure, Regulation, and Function. Pharmacol. Rev. 2010, 62, 405. [Google Scholar] [CrossRef]
  5. Paoletti, P.; Bellone, C.; Zhou, Q. NMDA receptor subunit diversity: Impact on receptor properties, synaptic plasticity and disease. Nat. Rev. Neurosci. 2013, 14, 383–400. [Google Scholar] [CrossRef]
  6. Kleckner, N.W.; Dingledine, R. Requirement for glycine in activation of NMDA-receptors expressed in Xenopus oocytes. Science 1988, 241, 835. [Google Scholar] [CrossRef]
  7. Furukawa, H.; Gouaux, E. Mechanisms of activation, inhibition and specificity: Crystal structures of the NMDA receptor NR1 ligand-binding core. EMBO J. 2003, 22, 2873–2885. [Google Scholar] [CrossRef]
  8. Tsai, G.E.; Yang, P.; Chang, Y.-C.; Chong, M.-Y. D-Alanine Added to Antipsychotics for the Treatment of Schizophrenia. Biol. Psychiatry 2006, 59, 230–234. [Google Scholar] [CrossRef]
  9. Johnson, J.W.; Ascher, P. Glycine potentiates the NMDA response in cultured mouse brain neurons. Nature 1987, 325, 529–531. [Google Scholar] [CrossRef]
  10. Benveniste, M.; Mayer, M.L. Structure-activity analysis of binding kinetics for NMDA receptor competitive antagonists: The influence of conformational restriction. Br. J. Pharmacol. 1991, 104, 207–221. [Google Scholar] [CrossRef]
  11. Clements, J.D.; Westbrook, G.L. Activation kinetics reveal the number of glutamate and glycine binding sites on the N-methyl-d-aspartate receptor. Neuron 1991, 7, 605–613. [Google Scholar] [CrossRef] [PubMed]
  12. Mayer, M.L.; Westbrook, G.L.; Guthrie, P.B. Voltage-dependent block by Mg2+ of NMDA responses in spinal cord neurones. Nature 1984, 309, 261–263. [Google Scholar] [CrossRef] [PubMed]
  13. Nowak, L.; Bregestovski, P.; Ascher, P.; Herbet, A.; Prochiantz, A. Magnesium gates glutamate-activated channels in mouse central neurones. Nature 1984, 307, 462–465. [Google Scholar] [CrossRef] [PubMed]
  14. Pagano, J.; Giona, F.; Beretta, S.; Verpelli, C.; Sala, C. N-methyl-d-aspartate receptor function in neuronal and synaptic development and signaling. Curr. Opin. Pharmacol. 2021, 56, 93–101. [Google Scholar] [CrossRef] [PubMed]
  15. Chou, T.-H.; Tajima, N.; Romero-Hernandez, A.; Furukawa, H. Structural Basis of Functional Transitions in Mammalian NMDA Receptors. Cell 2020, 182, 357–371.e13. [Google Scholar] [CrossRef]
  16. Jalali-Yazdi, F.; Chowdhury, S.; Yoshioka, C.; Gouaux, E. Mechanisms for Zinc and Proton Inhibition of the GluN1/GluN2A NMDA Receptor. Cell 2018, 175, 1520–1532.e15. [Google Scholar] [CrossRef]
  17. Regan, M.C.; Grant, T.; Mcdaniel, M.J.; Karakas, E.; Zhang, J.; Traynelis, S.F.; Grigorieff, N.; Furukawa, H. Structural Mechanism of Functional Modulation by Gene Splicing in NMDA Receptors. Neuron 2018, 98, 521–529.e3. [Google Scholar] [CrossRef]
  18. Tajima, N.; Karakas, E.; Grant, T.; Simorowski, N.; Diaz-Avalos, R.; Grigorieff, N.; Furukawa, H. Activation of NMDA receptors and the mechanism of inhibition by ifenprodil. Nature 2016, 534, 63–68. [Google Scholar] [CrossRef]
  19. Meyerson, J.R.; Kumar, J.; Chittori, S.; Rao, P.; Pierson, J.; Bartesaghi, A.; Mayer, M.L.; Subramaniam, S. Structural mechanism of glutamate receptor activation and desensitization. Nature 2014, 514, 328–334. [Google Scholar] [CrossRef]
  20. Karakas, E.; Furukawa, H. Crystal structure of a heterotetrameric NMDA receptor ion channel. Science 2014, 344, 992. [Google Scholar] [CrossRef]
  21. Papouin, T.; Oliet, S.H.R. Synaptic and Extra-Synaptic NMDA Receptors in the CNS. In The NMDA Receptors; Hashimoto, K., Ed.; Springer International Publishing: Cham, Switzerland, 2017; pp. 19–49. [Google Scholar]
  22. Furukawa, H.; Singh, S.K.; Mancusso, R.; Gouaux, E. Subunit arrangement and function in NMDA receptors. Nature 2005, 438, 185–192. [Google Scholar] [CrossRef] [PubMed]
  23. Balu, D.T. The NMDA Receptor and Schizophrenia: From Pathophysiology to Treatment. Adv. Pharmacol. 2016, 76, 351–382. [Google Scholar] [CrossRef] [PubMed]
  24. Jayawickrama, G.S.; Sadig, R.R.; Sun, G.; Nematollahi, A.; Nadvi, N.A.; Hanrahan, J.R.; Gorrell, M.D.; Church, W.B. Kynurenine Aminotransferases and the Prospects of Inhibitors for the Treatment of Schizophrenia. Curr. Med. Chem. 2015, 22, 2902–2918. [Google Scholar] [CrossRef]
  25. Coyle, J.T. NMDA receptor and schizophrenia: A brief history. Schizophr. Bull. 2012, 38, 920–926. [Google Scholar] [CrossRef] [PubMed]
  26. Guochuan, E.T.; Pao-Yen, L. Strategies to Enhance N-Methyl-D-Aspartate Receptor-Mediated Neurotransmission in Schizophrenia, a Critical Review and Meta-Analysis. Curr. Pharm. Des. 2010, 16, 522–537. [Google Scholar] [CrossRef]
  27. Majláth, Z.; Török, N.; Toldi, J.; Vécsei, L. Memantine and Kynurenic Acid: Current Neuropharmacological Aspects. Curr. Neuropharmacol. 2016, 14, 200–209. [Google Scholar] [CrossRef]
  28. Lau, C.G.; Zukin, R.S. NMDA receptor trafficking in synaptic plasticity and neuropsychiatric disorders. Nat. Rev. Neurosci. 2007, 8, 413–426. [Google Scholar] [CrossRef]
  29. Durham, R.J.; Paudyal, N.; Carrillo, E.; Bhatia, N.K.; Maclean, D.M.; Berka, V.; Dolino, D.M.; Gorfe, A.A.; Jayaraman, V. Conformational spread and dynamics in allostery of NMDA receptors. Proc. Natl. Acad. Sci USA. 2020, 117, 3839–3847. [Google Scholar] [CrossRef]
  30. Černý, J.; Božíková, P.; Balík, A.; Marques, S.M.; Vyklický, L. NMDA Receptor Opening and Closing—Transitions of a Molecular Machine Revealed by Molecular Dynamics. Biomolecules 2019, 9, 546. [Google Scholar] [CrossRef]
  31. Palmai, Z.; Houenoussi, K.; Cohen-Kaminsky, S.; Tchertanov, L. How does binding of agonist ligands control intrinsic molecular dynamics in human NMDA receptors? PLoS ONE 2018, 13, e0201234. [Google Scholar] [CrossRef]
  32. Sinitskiy, A.V.; Pande, V.S. Computer Simulations Predict High Structural Heterogeneity of Functional State of NMDA Receptors. Biophys. J. 2018, 115, 841–852. [Google Scholar] [CrossRef] [PubMed]
  33. Skrenkova, K.; Song, J.-M.; Kortus, S.; Kolcheva, M.; Netolicky, J.; Hemelikova, K.; Kaniakova, M.; Krausova, B.H.; Kucera, T.; Korabecny, J.; et al. The pathogenic S688Y mutation in the ligand-binding domain of the GluN1 subunit regulates the properties of NMDA receptors. Sci. Rep. 2020, 10, 18576. [Google Scholar] [CrossRef] [PubMed]
  34. Zehavi, Y.; Mandel, H.; Zehavi, A.; Rashid, M.A.; Straussberg, R.; Jabur, B.; Shaag, A.; Elpeleg, O.; Spiegel, R. De novo GRIN1 mutations: An emerging cause of severe early infantile encephalopathy. Eur. J. Med. Genet. 2017, 60, 317–320. [Google Scholar] [CrossRef] [PubMed]
  35. Ohba, C.; Shiina, M.; Tohyama, J.; Haginoya, K.; Lerman-Sagie, T.; Okamoto, N.; Blumkin, L.; Lev, D.; Mukaida, S.; Nozaki, F.; et al. GRIN1 mutations cause encephalopathy with infantile-onset epilepsy, and hyperkinetic and stereotyped movement disorders. Epilepsia 2015, 56, 841–848. [Google Scholar] [CrossRef] [PubMed]
  36. Lemke, J.R.; Geider, K.; Helbig, K.L.; Heyne, H.O.; Schütz, H.; Hentschel, J.; Courage, C.; Depienne, C.; Nava, C.; Heron, D.; et al. Delineating the GRIN1 phenotypic spectrum. Neurology 2016, 86, 2171. [Google Scholar] [CrossRef] [PubMed]
  37. Hamdan, F.F.; Gauthier, J.; Araki, Y.; Lin, D.T.; Yoshizawa, Y.; Higashi, K.; Park, A.R.; Spiegelman, D.; Dobrzeniecka, S.; Piton, A.; et al. Excess of De Novo Deleterious Mutations in Genes Associated with Glutamatergic Systems in Nonsyndromic Intellectual Disability. Am. J. Hum. Genet. 2011, 88, 306–316. [Google Scholar] [CrossRef] [PubMed]
  38. Redin, C.; Gérard, B.; Lauer, J.; Herenger, Y.; Muller, J.; Quartier, A.; Masurel-Paulet, A.; Willems, M.; Lesca, G.; El-Chehadeh, S.; et al. Efficient strategy for the molecular diagnosis of intellectual disability using targeted high-throughput sequencing. J. Med. Genet. 2014, 51, 724. [Google Scholar] [CrossRef]
  39. Chen, W.; Shieh, C.; Swanger, S.A.; Tankovic, A.; Au, M.; Mcguire, M.; Tagliati, M.; Graham, J.M.; Madan-Khetarpal, S.; Traynelis, S.F.; et al. GRIN1 mutation associated with intellectual disability alters NMDA receptor trafficking and function. J. Hum. Genet. 2017, 62, 589–597. [Google Scholar] [CrossRef]
  40. Yu, W.; MacKerell, A.D. Computer-Aided Drug Design Methods. In Antibiotics: Methods and Protocols; Sass, P., Ed.; Springer: New York, NY, USA, 2017; pp. 85–106. [Google Scholar]
  41. Chen, W.; Tankovic, A.; Burger, P.B.; Kusumoto, H.; Traynelis, S.F.; Yuan, H. Functional Evaluation of a De Novo GRIN2A Mutation Identified in a Patient with Profound Global Developmental Delay and Refractory Epilepsy. Mol. Pharmacol. 2017, 91, 317. [Google Scholar] [CrossRef]
  42. Ogden, K.K.; Chen, W.; Swanger, S.A.; Mcdaniel, M.J.; Fan, L.Z.; Hu, C.; Tankovic, A.; Kusumoto, H.; Kosobucki, G.J.; Schulien, A.J.; et al. Molecular Mechanism of Disease-Associated Mutations in the Pre-M1 Helix of NMDA Receptors and Potential Rescue Pharmacology. PLoS Genet. 2017, 13, e1006536. [Google Scholar] [CrossRef]
  43. Rossi, M.; Chatron, N.; Labalme, A.; Ville, D.; Carneiro, M.; Edery, P.; des Portes, V.; Lemke, J.R.; Sanlaville, D.; Lesca, G. Novel homozygous missense variant of GRIN1 in two sibs with intellectual disability and autistic features without epilepsy. Eur. J. Hum. Genet. EJHG 2017, 25, 376–380. [Google Scholar] [CrossRef] [PubMed]
  44. Delano, W.L. The PyMOL Molecular Graphics System. 2002. Available online: http://www.pymol.org (accessed on 12 December 2022).
  45. Schrödinger. Schrödinger Release 2020–2023, Schrödinger LLC: New York, NY, USA, 2020.
  46. Papouin, T.; Ladépêche, L.; Ruel, J.; Sacchi, S.; Labasque, M.; Hanini, M.; Groc, L.; Pollegioni, L.; Mothet, J.P.; Oliet, S.H. Synaptic and Extrasynaptic NMDA Receptors Are Gated by Different Endogenous Coagonists. Cell 2012, 150, 633–646. [Google Scholar] [CrossRef] [PubMed]
  47. Wang, H.; Lv, S.; Stroebel, D.; Zhang, J.; Pan, Y.; Huang, X.; Zhang, X.; Paoletti, P.; Zhu, S. Gating mechanism and a modulatory niche of human GluN1-GluN2A NMDA receptors. Neuron 2021, 109, 2443–2456.e5. [Google Scholar] [CrossRef] [PubMed]
  48. Le Bail, M.; Martineau, M.; Sacchi, S.; Yatsenko, N.; Radzishevsky, I.; Conrod, S.; Ouares, K.A.; Wolosker, H.; Pollegioni, L.; Billard, J.M.; et al. Identity of the NMDA receptor coagonist is synapse specific and developmentally regulated in the hippocampus. Proc. Natl. Acad. Sci. USA 2015, 112, E204. [Google Scholar] [CrossRef] [PubMed]
  49. Maurer, M.; Oostenbrink, C. Water in protein hydration and ligand recognition. J. Mol. Recognit. JMR 2019, 32, e2810. [Google Scholar] [CrossRef] [PubMed]
  50. Fossat, P.; Turpin, F.R.; Sacchi, S.; Dulong, J.; Shi, T.; Rivet, J.-M.; Sweedler, J.V.; Pollegioni, L.; Millan, M.J.; Oliet, S.H.; et al. Glial D-Serine Gates NMDA Receptors at Excitatory Synapses in Prefrontal Cortex. Cereb. Cortex 2012, 22, 595–606. [Google Scholar] [CrossRef]
  51. Curcio, L.; Podda, M.V.; Leone, L.; Piacentini, R.; Mastrodonato, A.; Cappelletti, P.; Sacchi, S.; Pollegioni, L.; Grassi, C.; D’Ascenzo, M. Reduced d-serine levels in the nucleus accumbens of cocaine-treated rats hinder the induction of NMDA receptor-dependent synaptic plasticity. Brain 2013, 136, 1216–1230. [Google Scholar] [CrossRef]
  52. Basu, A.C.; Tsai, G.E.; Ma, C.L.; Ehmsen, J.T.; Mustafa, A.K.; Han, L.; Jiang, Z.I.; Benneyworth, M.A.; Froimowitz, M.P.; Lange, N.; et al. Targeted disruption of serine racemase affects glutamatergic neurotransmission and behavior. Mol. Psychiatry 2009, 14, 719–727. [Google Scholar] [CrossRef]
  53. Pu, C.; Yan, G.; Shi, J.; Li, R. Assessing the performance of docking scoring function, FEP, MM-GBSA, and QM/MM-GBSA approaches on a series of PLK1 inhibitors. MedChemComm 2017, 8, 1452–1458. [Google Scholar] [CrossRef]
  54. Fratev, F.; Sirimulla, S. An Improved Free Energy Perturbation FEP+ Sampling Protocol for Flexible Ligand-Binding Domains. Sci. Rep. 2019, 9, 16829. [Google Scholar] [CrossRef]
  55. Jorgensen, W.L.; Thomas, L.L. Perspective on Free-Energy Perturbation Calculations for Chemical Equilibria. J. Chem. Theory Comput. 2008, 4, 869–876. [Google Scholar] [CrossRef] [PubMed]
  56. Pekel, H.; Ilter, M.; Sensoy, O. Inhibition of SARS-CoV-2 main protease: A repurposing study that targets the dimer interface of the protein. J. Biomol. Struct. Dyn. 2021, 40, 7167–7182. [Google Scholar] [CrossRef]
  57. Chakraborti, S.; Chakraborty, M.; Bose, A.; Srinivasan, N.; Visweswariah, S.S. Identification of Potential Binders of Mtb Universal Stress Protein (Rv1636) through an in silico Approach and Insights into Compound Selection for Experimental Validation. Front. Mol. Biosci. 2021, 8, 177. [Google Scholar] [CrossRef] [PubMed]
  58. Zielkiewicz, J. Structural properties of water: Comparison of the SPC, SPCE, TIP4P, and TIP5P models of water. J. Chem. Phys. 2005, 123, 104501. [Google Scholar] [CrossRef] [PubMed]
  59. Mark, P.; Nilsson, L. Structure and Dynamics of the TIP3P, SPC, and SPC/E Water Models at 298 K. J. Phys. Chem. A 2001, 105, 9954–9960. [Google Scholar] [CrossRef]
  60. Bowers, K.J.; Chow, E.; Xu, H.; Dror, R.O.; Eastwood, M.P.; Gregersen, B.A.; Klepeis, J.L.; Kolossvary, I.; Moraes, M.A.; Sacerdoti, F.D.; et al. Scalable algorithms for molecular dynamics simulations on commodity clusters. In Proceedings of the 2006 ACM/IEEE conference on Supercomputing, Association for Computing Machinery, Tampa, FL, USA, 11–17 November 2006; pp. 84–es. [Google Scholar]
  61. Jespersen, A.; Tajima, N.; Fernandez-Cuervo, G.; Garnier-Amblard, E.C.; Furukawa, H. Structural Insights into Competitive Antagonism in NMDA Receptors. Neuron 2014, 81, 366–378. [Google Scholar] [CrossRef]
  62. Sastry, G.M.; Adzhigirey, M.; Day, T.; Annabhimoju, R.; Sherman, W. Protein and ligand preparation: Parameters, protocols, and influence on virtual screening enrichments. J. Comput. Aided Mol. Des. 2013, 27, 221–234. [Google Scholar] [CrossRef] [PubMed]
  63. Olsson, M.H.M.; Søndergaard, C.R.; Rostkowski, M.; Jensen, J.H. PROPKA3: Consistent Treatment of Internal and Surface Residues in Empirical pKa Predictions. J. Chem. Theory Comput. 2011, 7, 525–537. [Google Scholar] [CrossRef]
  64. Søndergaard, C.R.; Olsson, M.H.M.; Rostkowski, M.; Jensen, J.H. Improved Treatment of Ligands and Coupling Effects in Empirical Calculation and Rationalization of pKa Values. J. Chem. Theory Comput. 2011, 7, 2284–2295. [Google Scholar] [CrossRef]
  65. Roos, K.; Wu, C.; Damm, W.; Reboul, M.; Stevenson, J.M.; Lu, C.; Dahlgren, M.K.; Mondal, S.; Chen, W.; Wang, L.; et al. OPLS3e: Extending Force Field Coverage for Drug-Like Small Molecules. J. Chem. Theory Comput. 2019, 15, 1863–1874. [Google Scholar] [CrossRef]
  66. Friesner, R.A.; Murphy, R.B.; Repasky, M.P.; Frye, L.L.; Greenwood, J.R.; Halgren, T.A.; Sanschagrin, P.C.; Mainz, D.T. Extra Precision Glide:  Docking and Scoring Incorporating a Model of Hydrophobic Enclosure for Protein−Ligand Complexes. J. Med. Chem. 2006, 49, 6177–6196. [Google Scholar] [CrossRef] [PubMed]
  67. Halgren, T.A.; Murphy, R.B.; Friesner, R.A.; Beard, H.S.; Frye, L.L.; Pollard, W.T.; Banks, J.L. Glide:  A New Approach for Rapid, Accurate Docking and Scoring. 2. Enrichment Factors in Database Screening. J. Med. Chem. 2004, 47, 1750–1759. [Google Scholar] [CrossRef] [PubMed]
  68. Friesner, R.A.; Banks, J.L.; Murphy, R.B.; Halgren, T.A.; Klicic, J.J.; Mainz, D.T.; Repasky, M.P.; Knoll, E.H.; Shelley, M.; Perry, J.K.; et al. Glide:  A New Approach for Rapid, Accurate Docking and Scoring. 1. Method and Assessment of Docking Accuracy. J. Med. Chem. 2004, 47, 1739–1749. [Google Scholar] [CrossRef] [PubMed]
  69. Banks, J.L.; Beard, H.S.; Cao, Y.; Cho, A.E.; Damm, W.; Farid, R.; Halgren, T.A.; Mainz, D.T.; Maple, J.R.; Murphy, R. Integrated Modeling Program, Applied Chemical Theory (IMPACT). J. Comput. Chem. 2005, 26, 1752–1780. [Google Scholar] [CrossRef]
  70. Berendsen, H.J.C.; Postma, J.P.M.; Van Gunsteren, W.F.; Hermans, J. Interaction Models for Water in Relation to Protein Hydration. In Intermolecular Forces, Proceedings of the Fourteenth Jerusalem Symposium on Quantum Chemistry and Biochemistry, Jerusalem, Israel, 13–16 April 1981; Pullman, B., Ed.; Springer: Dordrecht, The Netherlands, 1981; pp. 331–342. [Google Scholar]
  71. Genheden, S.; Ryde, U. The MM/PBSA and MM/GBSA methods to estimate ligand-binding affinities. Expert Opin. Drug Discov. 2015, 10, 449–461. [Google Scholar] [CrossRef]
  72. Hou, T.; Wang, J.; Li, Y.; Wang, W. Assessing the Performance of the MM/PBSA and MM/GBSA Methods. 1. The Accuracy of Binding Free Energy Calculations Based on Molecular Dynamics Simulations. J. Chem. Inf. Model. 2011, 51, 69–82. [Google Scholar] [CrossRef]
  73. Jacobson, M.P.; Pincus, D.L.; Rapp, C.S.; Day TJ, F.; Honig, B.; Shaw, D.E.; Friesner, R.A. A hierarchical approach to all-atom protein loop prediction. Proteins Struct. Funct. Bioinform. 2004, 55, 351–367. [Google Scholar] [CrossRef]
  74. Jacobson, M.P.; Friesner, R.A.; Xiang, Z.; Honig, B. On the Role of the Crystal Environment in Determining Protein Side-chain Conformations. J. Mol. Biol. 2002, 320, 597–608. [Google Scholar] [CrossRef]
  75. Li, J.; Abel, R.; Zhu, K.; Cao, Y.; Zhao, S.; Friesner, R.A. The VSGB 2.0 model: A next generation energy model for high resolution protein structure modeling. Proteins 2011, 79, 2794–2812. [Google Scholar] [CrossRef]
  76. Fry, A.E.; Fawcett, K.A.; Zelnik, N.; Yuan, H.; Thompson, B.A.N.; Shemer-Meiri, L.; Cushion, T.D.; Mugalaasi, H.; Sims, D.; Stoodley, N.; et al. De novo mutations in GRIN1 cause extensive bilateral polymicrogyria. Brain 2018, 141, 698–712. [Google Scholar] [CrossRef]
  77. Samura, O. Update on noninvasive prenatal testing: A review based on current worldwide research. J. Obstet. Gynaecol. Res. 2020, 46, 1246–1254. [Google Scholar] [CrossRef] [PubMed]
Molecules 28 04108 g002 550
Figure 2. 2D interactions diagram of ligands bound to GluN1 LBD. Hydrophobic residues are shown in green, polar residues are shown in cyan, positively charged residues are shown in purple, and negatively charged residues are shown in orange. Hydrogen bonds are shown as purple arrows, π-cation interactions are shown as red lines, and salt bridges are shown as red-blue lines. The tip of the tear-drop shape points towards the side chain of the residue; dots on a connection indicate a residue not making contact with the ligand. (A) Glycine bound to the wild-type GluN1 LBD (B) D-serine bound to the wild-type GluN1 LBD. (C) Glycine bound to the S688Y GluN1 LBD (D) D-serine bound to the S688Y GluN1 LBD. Figure generated using the ligand interaction visualization module in Maestro 2021.4 [44,45].
Figure 2. 2D interactions diagram of ligands bound to GluN1 LBD. Hydrophobic residues are shown in green, polar residues are shown in cyan, positively charged residues are shown in purple, and negatively charged residues are shown in orange. Hydrogen bonds are shown as purple arrows, π-cation interactions are shown as red lines, and salt bridges are shown as red-blue lines. The tip of the tear-drop shape points towards the side chain of the residue; dots on a connection indicate a residue not making contact with the ligand. (A) Glycine bound to the wild-type GluN1 LBD (B) D-serine bound to the wild-type GluN1 LBD. (C) Glycine bound to the S688Y GluN1 LBD (D) D-serine bound to the S688Y GluN1 LBD. Figure generated using the ligand interaction visualization module in Maestro 2021.4 [44,45].
Molecules 28 04108 g002
Molecules 28 04108 g003 550
Figure 3. Top binding poses of glycine and D-serine at the GluN1 LBD (PDB: 4NF8). Carbon atoms of ligands are shown in magenta, with interacting residues shown in light blue and the ribbon of the receptor shown in salmon. Hydrogen bonds are shown as yellow dashed lines; ionic interactions are shown as purple dashed lines; π-cation interactions are shown as green dashed lines and aromatic hydrogen bonds are shown as cyan dashed lines. (A) Glycine bound to wild-type GluN1. (B) D-serine bound to wild-type GluN1. (C) Glycine bound to S688Y GluN1. (D) D-serine bound to S688Y GluN1. Figure generated using PyMOL [44].
Figure 3. Top binding poses of glycine and D-serine at the GluN1 LBD (PDB: 4NF8). Carbon atoms of ligands are shown in magenta, with interacting residues shown in light blue and the ribbon of the receptor shown in salmon. Hydrogen bonds are shown as yellow dashed lines; ionic interactions are shown as purple dashed lines; π-cation interactions are shown as green dashed lines and aromatic hydrogen bonds are shown as cyan dashed lines. (A) Glycine bound to wild-type GluN1. (B) D-serine bound to wild-type GluN1. (C) Glycine bound to S688Y GluN1. (D) D-serine bound to S688Y GluN1. Figure generated using PyMOL [44].
Molecules 28 04108 g003
Molecules 28 04108 g004 550
Figure 4. Averaged triplicate protein and ligand RMSD for GluN1 LBD with glycine and D-serine bound over the course of 1000 ns (1 μs) molecular dynamics simulation. RMSD values calculated for Cα are shown in blue, RMSD values calculated for residues within 5 Å of the ligand are shown in black, and ligand RMSD fit to protein is shown in red. (A) Wild-type with glycine. (B) Wild-type with D-serine. (C) S688Y with glycine. (D) S688Y with D-serine. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Figure 4. Averaged triplicate protein and ligand RMSD for GluN1 LBD with glycine and D-serine bound over the course of 1000 ns (1 μs) molecular dynamics simulation. RMSD values calculated for Cα are shown in blue, RMSD values calculated for residues within 5 Å of the ligand are shown in black, and ligand RMSD fit to protein is shown in red. (A) Wild-type with glycine. (B) Wild-type with D-serine. (C) S688Y with glycine. (D) S688Y with D-serine. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Molecules 28 04108 g004
Molecules 28 04108 g005 550
Figure 5. Averaged triplicate root–mean–square fluctuation (RMSF) of the GluN1 receptor with bound ligands. (A) Wild-type GluN1 with glycine bound. (B) Wild-type GluN1 with D-serine bound. (C) S688Y GluN1 mutant with glycine bound. (D) S688Y GluN1 mutant with D-serine bound. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Figure 5. Averaged triplicate root–mean–square fluctuation (RMSF) of the GluN1 receptor with bound ligands. (A) Wild-type GluN1 with glycine bound. (B) Wild-type GluN1 with D-serine bound. (C) S688Y GluN1 mutant with glycine bound. (D) S688Y GluN1 mutant with D-serine bound. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Molecules 28 04108 g005
Molecules 28 04108 g006 550
Figure 6. Averaged triplicate Solvent Accessible Surface Area (SASA) of the ligand during the microsecond simulations. Wild-type is shown in purple, and S688Y mutant receptor is shown in teal. (A) Glycine bound to GluN1 LBD; (B) D-serine bound to GluN1 LBD. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Figure 6. Averaged triplicate Solvent Accessible Surface Area (SASA) of the ligand during the microsecond simulations. Wild-type is shown in purple, and S688Y mutant receptor is shown in teal. (A) Glycine bound to GluN1 LBD; (B) D-serine bound to GluN1 LBD. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Molecules 28 04108 g006
Molecules 28 04108 g007 550
Figure 7. Averaged triplicate GluN1 LBD-ligand interactions of the 1 µs simulations. Note that figure scales are different to improve the visibility of the data. Hydrogen bonds are defined as an H-Acceptor distance less than 2.8 Å and a Donor-H-Acceptor angle greater than 120. Hydrophobic interactions include pi-pi stacking, pi-cation interactions, and van der Waals interactions within 3.6 Å of ligand. Ionic interactions are defined as charged interactions within 3.6 Å of the ligand. Water bridges are defined as H-Acceptor distances less than 2.7 Å. Hydrogen bonds are shown in green, hydrophobic interactions are shown in lilac, ionic interactions are shown in magenta, and water bridges are shown in blue. For residues with more than one type of interaction, the interaction fraction can exceed 1. (A) Glycine bound to wild-type GluN1 LBD, (B) D-serine bound to wild-type GluN1 LBD, (C) Glycine bound to S688Y GluN1 LBD, (D) D-serine bound to S688Y GluN1 LBD. Graphs were generated using the Simulation Interactions Diagram module of Maestro 2021.4 [47,48].
Figure 7. Averaged triplicate GluN1 LBD-ligand interactions of the 1 µs simulations. Note that figure scales are different to improve the visibility of the data. Hydrogen bonds are defined as an H-Acceptor distance less than 2.8 Å and a Donor-H-Acceptor angle greater than 120. Hydrophobic interactions include pi-pi stacking, pi-cation interactions, and van der Waals interactions within 3.6 Å of ligand. Ionic interactions are defined as charged interactions within 3.6 Å of the ligand. Water bridges are defined as H-Acceptor distances less than 2.7 Å. Hydrogen bonds are shown in green, hydrophobic interactions are shown in lilac, ionic interactions are shown in magenta, and water bridges are shown in blue. For residues with more than one type of interaction, the interaction fraction can exceed 1. (A) Glycine bound to wild-type GluN1 LBD, (B) D-serine bound to wild-type GluN1 LBD, (C) Glycine bound to S688Y GluN1 LBD, (D) D-serine bound to S688Y GluN1 LBD. Graphs were generated using the Simulation Interactions Diagram module of Maestro 2021.4 [47,48].
Molecules 28 04108 g007
Molecules 28 04108 g008 550
Figure 8. Averaged protein and ligand RMSD for GluN1 LBD with glycine and D-serine bound over the course of quintuplet 200 ns molecular dynamics simulation. The results shown are the average of five randomly seeded replicates. RMSD values calculated for Cα are shown in blue, and ligand RMSD fit to protein is shown in red. (A) Wild-type with glycine. (B) Wild-type with D-serine. (C) S688Y with glycine. (D) S688Y with D-serine. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Figure 8. Averaged protein and ligand RMSD for GluN1 LBD with glycine and D-serine bound over the course of quintuplet 200 ns molecular dynamics simulation. The results shown are the average of five randomly seeded replicates. RMSD values calculated for Cα are shown in blue, and ligand RMSD fit to protein is shown in red. (A) Wild-type with glycine. (B) Wild-type with D-serine. (C) S688Y with glycine. (D) S688Y with D-serine. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Molecules 28 04108 g008
Molecules 28 04108 g009 550
Figure 9. Averaged root–mean–square fluctuation (RMSF) of the GluN1 receptor with bound ligands. The results shown are the average of five randomly seeded replicates. (A) Wild-type GluN1 with glycine bound. (B) Wild-type GluN1 with D-serine bound. (C) S688Y GluN1 mutant with glycine bound. (D) S688Y GluN1 mutant with D-serine bound. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Figure 9. Averaged root–mean–square fluctuation (RMSF) of the GluN1 receptor with bound ligands. The results shown are the average of five randomly seeded replicates. (A) Wild-type GluN1 with glycine bound. (B) Wild-type GluN1 with D-serine bound. (C) S688Y GluN1 mutant with glycine bound. (D) S688Y GluN1 mutant with D-serine bound. Graphs were generated using GraphPad Prism version 9.5.0 for MacOS, GraphPad Software, San Diego, CA, USA.
Molecules 28 04108 g009
Table
Table 1. Binding free energy, ΔGbind, values for glycine and D-serine bound to the wild-type and S688Y mutant receptor. The values are averaged values of the five randomly seeded replicates.
Table 1. Binding free energy, ΔGbind, values for glycine and D-serine bound to the wild-type and S688Y mutant receptor. The values are averaged values of the five randomly seeded replicates.
System ΔGbindΔGbind CoulombΔGbind CovalentΔGbind H-BondΔGbind LipophilicityΔGbind SolventΔGbind Van der Waals
wt
Glycine		Mean (kcal/mol)−30.0−5.620.09−5.05−2.43−7.04−9.97
Standard Deviation3.084.270.420.540.153.562.64
wt
D-serine		Mean−34.6−6.621.09−5.00−4.03−6.18−13.8
Standard Deviation3.214.900.540.360.183.472.45
S688Y
Glycine		Mean−26.311.8−0.11−3.08−2.80−18.9−13.2
Standard Deviation1.943.440.590.260.103.331.78
S688Y
D-serine		Mean−22.77.580.84−3.32−3.18−9.65−14.9
Standard Deviation2.033.420.630.290.312.921.32
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Chen, J.Z.; Church, W.B.; Bastard, K.; Duff, A.P.; Balle, T. Binding and Dynamics Demonstrate the Destabilization of Ligand Binding for the S688Y Mutation in the NMDA Receptor GluN1 Subunit. Molecules 2023, 28, 4108. https://doi.org/10.3390/molecules28104108

AMA Style

Chen JZ, Church WB, Bastard K, Duff AP, Balle T. Binding and Dynamics Demonstrate the Destabilization of Ligand Binding for the S688Y Mutation in the NMDA Receptor GluN1 Subunit. Molecules. 2023; 28(10):4108. https://doi.org/10.3390/molecules28104108

Chicago/Turabian Style

Chen, Jake Zheng, William Bret Church, Karine Bastard, Anthony P. Duff, and Thomas Balle. 2023. "Binding and Dynamics Demonstrate the Destabilization of Ligand Binding for the S688Y Mutation in the NMDA Receptor GluN1 Subunit" Molecules 28, no. 10: 4108. https://doi.org/10.3390/molecules28104108

APA Style

Chen, J. Z., Church, W. B., Bastard, K., Duff, A. P., & Balle, T. (2023). Binding and Dynamics Demonstrate the Destabilization of Ligand Binding for the S688Y Mutation in the NMDA Receptor GluN1 Subunit. Molecules, 28(10), 4108. https://doi.org/10.3390/molecules28104108

Article Metrics

Back to TopTop