Benzene Probes in Molecular Dynamics Simulations Reveal Novel Binding Sites for Ligand Design

Protein flexibility poses a major challenge in binding site identification. Several computational pocket detection methods that utilize small-molecule probes in molecular dynamics (MD) simulations have been developed to address this issue. Although they have proven hugely successful at reproducing experimental structural data, their ability to predict new binding sites that are yet to be identified and characterized has not been demonstrated. Here, we report the use of benzenes as probe molecules in ligand-mapping MD (LMMD) simulations to predict the existence of two novel binding sites on the surface of the oncoprotein MDM2. One of them was serendipitously confirmed by biophysical assays and X-ray crystallography to be important for the binding of a new family of hydrocarbon stapled peptides that were specifically designed to target the other putative site. These results highlight the predictive power of LMMD and suggest that predictions derived from LMMD simulations can serve as a reliable basis for the identification of novel ligand binding sites in structure-based drug design.

T he identification and characterization of binding pockets is an important component of structure-based drug design. It is also often the first step in evaluating the druggability of a protein target. 1,2 In recent decades, various computational algorithms and methods that rely on the use of static protein structures have been developed for rapid identification of binding pockets for ligand design. 1 They are, however, severely limited by their dependence on available protein structures. Proteins are intrinsically flexible and frequently undergo conformational changes on ligand binding. 3−6 A major concern is that cryptic binding pockets that are absent in the input structures and present themselves only in the presence of a suitable interacting ligand will be missed. This is often the case for hydrophobic pockets, which tend to remain occluded in polar solvents and open up only in the presence of less polar ligands. 7 To address this issue, there have been recent efforts to develop molecular dynamics (MD)-based methods that incorporate small molecules into the protein's solvent box for pocket detection. 8−12 In these simulations, the probes interact dynamically with the protein surface, allowing for ligandinduced conformational changes. The use of hydrophobic probes is of particular interest because it reduces the solvent polarity, thus facilitating the opening and enlargement of hydrophobic pockets that may otherwise remain undetected in pure water simulations of the protein. 7 Ligand-mapping MD (LMMD) 13,14 is one of two probe-based MD simulation methods that employ hydrophobic probes for pocket detection. In contrast to the related site identification by ligand competitive saturation (SILCS) method, 9 LMMD does not require the addition of artificial interligand repulsive energy terms because of the use of relatively low concentrations of hydrophobic probes to avoid ligand aggregation. LMMD simulations have been shown to be especially useful at revealing cryptic binding sites 14 and were previously used to guide the design of a ligand to target a cryptic pocket. 13 Recently, LMMD has also been established as a reliable method for the identification of hydrophobic peptide binding sites. 15 To date, probe-based MD simulations have mostly been limited to the reproduction of known structural data. Unlike the non MDbased pocket detection methods, 16 there have been no previous reports of the successful prediction of a previously unknown binding site by these simulations, although a recent study suggests that SILCS has the potential to propose alternative binding sites. 17 A demonstration of the predictive power of probe-based MD simulations will provide confidence for and encourage their application in structure-based drug design projects. Here, we focus on the promising anticancer therapeutic target MDM2 as a prototypical example for the detection of novel ligand binding sites by LMMD.
The E3 ubiquitin ligase MDM2 is a potent inhibitor of the tumor suppressor protein p53, 18 which plays an essential role in coordinating cellular responses, including cell cycle arrest, apoptosis, and senescence, to a variety of stress signals. 19 MDM2 binds to the transactivation domain of p53 to block p53-mediated transactivation 20 and targets it for ubiquitinmediated proteolysis. 21 It is overexpressed in many cancers and is thought to be one of the primary causes of p53 network inactivation in p53 wild-type (WT) tumors. 22 Antagonists of the MDM2−p53 interaction can reactivate the p53 response, leading to cell cycle arrest and apoptosis in tumor cells. 23,24 Several small-molecule inhibitors of the MDM2−p53 interaction have been developed, and some of them have reached clinical trials. 25,26 These molecules mimic the three key binding residues (Phe19, Trp23, and Leu26) in the p53 transactivation domain, which binds as an amphipathic α-helix to a deep hydrophobic cleft in the N-terminal domain of MDM2. 27 Besides small molecules, peptides derived from the transactivation domain of p53 have also been used to inhibit the MDM2−p53 interaction. Unlike small molecules, however, linear peptides are susceptible to proteolytic cleavage, lack a well-defined conformation prior to target engagement, and are poorly cell-permeable. 28 These shortcomings can potentially be overcome by hydrocarbon stapling, in which two unnatural residues bearing olefin side chains of varying lengths are introduced into the α-helix of the peptide, followed by a ruthenium-catalyzed ring-closing metathesis reaction to form a covalent staple across one or two α-helical turns. 29 Hydrocarbon stapling generally leads to enhanced helicity, protease resistance, and cell permeability. 30,31 Three families of MDM2binding hydrocarbon stapled peptides have been identified to date: one derived from the WT p53 sequence 32 and two others from phage display studies. 33,34 Two of these peptides have been cocrystallized with MDM2, revealing that the hydrocarbon staple also interacts with the protein surface, thus augmenting the hydrophobic interactions at the binding interface. 35,36 In this Letter, we describe the use of LMMD to identify two new binding sites close to the consensus p53-binding cleft in MDM2. Subsequent biophysical and crystallographic studies showed that one of these sites is important for the binding of a novel set of hydrocarbon stapled peptides.
The crystal structure of human MDM2 in complex with the p53 transactivation domain peptide (PDB code 1YCR) 27 was used as the initial structure for the LMMD simulations. As cryptic binding sites tend to be hydrophobic in nature, 7 benzene was chosen as the mapping ligand to better explore such sites. An initial set of 10 independent 5 ns LMMD simulations with different starting benzene distributions was performed on apo MDM2 after removing p53 from the complex structure. A relatively low benzene concentration (0.2 M) was used to prevent phase separation and ligand aggregation on the protein surface, which may cause denaturation. 37 A second set of LMMD simulations was performed on p53-bound MDM2 to determine whether the binding sites that are discrete from the p53-binding cleft could be reproduced in the presence of the peptide.
Analysis of the benzene occupancy maps generated from the apo and holo LMMD simulations revealed two putative discrete binding sites close to the p53-binding cleft, which have never been exploited by any known inhibitors of MDM2 and may be utilized to improve the affinity and selectivity of existing

The Journal of Physical Chemistry Letters
Letter ligands. The first binding site lies between Tyr100 and Tyr104 and was mapped by the benzene probes in both apo and holo LMMD simulations. It was found to correspond to a crystal packing contact for a second nutlin molecule in the crystal structure of its complex with MDM2 ( Figure 1A). 38 This secondary nutlin interaction site is defined by seven residues: Val28, Met50, Tyr100, Ile103, Tyr104, Leu107, and Val109 ( Figure 1B). These residues potentially constitute a cryptic binding pocket that is occluded by Tyr104 in the majority of MDM2 structures, including the one used for the LMMD simulations. There is experimental evidence to suggest that the region around this cryptic site may be relevant for nutlin binding at the p53-binding cleft; 39 however, its ability to bind ligands has not been verified.
The second putative binding site was identified in the apo LMMD simulations within a region contiguous to where Pro27 of p53 binds ( Figure 1C); however, it was not detected in the corresponding holo simulations because it was occluded by the C-terminal region of the p53 peptide. This "proximal P27 site" is essentially a relatively shallow hydrophobic basin constituted by Thr26, Met50, Lys51, Leu54, Tyr100, and Tyr104 ( Figure  1D).
The proximity of the second nutlin interaction site to the p53-binding cleft raises the possibility of accessing it by extension of the p53 peptide. Suitable MDM2 structures, in which the identified cryptic pocket was present, were selected from the apo and holo LMMD simulations. The p53 peptide was modeled onto the apo MDM2 structure by extraction from the 1YCR crystal structure. Either a phenylalanine or a tyrosine residue was then appended to the peptide's C-terminus, such that there was optimal overlap of the phenyl ring with the benzene molecule bound at the second nutlin interaction site ( Figure 2A). The peptide in the holo MDM2 structure was similarly extended. We then replaced the WT p53 sequence in these extended peptides with that of a previously reported potent MDM2-binding and p53-activating stapled peptide, sMTide-02. 33 The i, i + 7 staple used in sMTide-02, however, could not be transferred to these peptides because it would induce an additional helical turn that prevents the appended residue from accessing the second nutlin interaction site. Instead, a shorter i, i + 4 staple formed by two (R)-2-(4′pentenyl)alanine residues (R 5 ) was used to maintain the extended structure of the C-terminal residues. The R,R configuration of the hydrocarbon staple was preferred to the more common S,S configuration because it allows the staple to point toward and form intimate hydrophobic contacts with the MDM2 surface. These novel stapled peptides were termed YS-1 (TSFR 5 EYWR 5 LLPENF) and YS-2 (TSFR 5 EYWR 5 LLPENY), respectively.
MD simulations of the complexes of MDM2 with WT p53, sMTide-02, YS-1, and YS-2 were performed for 50 ns each to evaluate their stability. The structures of YS-1 and YS-2 did not deviate significantly from their initial minimized structures, with root-mean-square deviation (RMSD) values remaining below 1.6 Å ( Figure S1). The Phe ( Figure 2B) and Tyr extensions (not shown) remained bound at their targeted binding sites at the end of the simulations, indicating the stability of the predicted binding mode. Binding free energies of the various MDM2−peptide complexes were then estimated using the molecular mechanics/generalized Born surface area (MM/ GBSA) method. 40 The sMTide-02 peptide was predicted to have a higher binding affinity for MDM2 compared to WT p53 peptide (Table S5), in agreement with previous experimental results. 33 Notably, the calculations suggest that both YS-1 and YS-2 are significantly more potent binders of MDM2 than sMTide-02 (Table S5). To validate the computational predictions, the stapled peptides were subsequently characterized in biophysical and X-ray crystallographic experiments.
Binding affinities of the stapled peptides for MDM2 were determined in competitive fluorescence polarization (FP) assays ( Table 1). The designed stapled peptides exhibited approximately up to 5-fold improvement in binding potency relative to sMTide-02, in qualitative agreement with the computational predictions. Incorporation of a tyrosine instead of phenylalanine at the C-terminus of the stapled peptides made little difference to their binding affinities.
Crystal structures of the N-terminal domain of MDM2 in complex with YS-1 and YS-2 were solved to elucidate their binding modes. The two peptides bind to MDM2 very similarly, with only minor variations for some solvent-exposed residues ( Figure 3A and Figure S3). Both stapled peptides adopt an α-helical structure with Phe19, Trp23, and Leu26  The Journal of Physical Chemistry Letters Letter bound in their consensus pockets. The i, i + 4 hydrocarbon staple packs snugly against the MDM2 glycine shelf in a manner broadly similar to that seen for longer i, i + 7 staples of other MDM2-binding stapled peptides. 35,36 Surprisingly, the last three residues of YS-1 and YS-2 do not adopt the linear conformation observed at the C-terminus of WT p53 when bound to MDM2; instead they extend the helical fold to create a kinked helix encompassing residues 19− 30 ( Figure 4A, residue numbers follow that in WT p53 sequence). While Pro27 acts as a helix breaker in the p53 peptide because of steric constraints and its inability to maintain the backbone hydrogen bonding network, it is unexpectedly incorporated into the helices of YS-1 and YS-2, allowing it to form close hydrophobic contacts with the MDM2 surface. It is likely that this is driven by the presence of the hydrocarbon staple, which has the ability to stabilize α-helices beyond the stapled α-helical turn. 41 This results in the Ctermini of the peptides packing into the cleft between helices α3 and α5 ( Figure 3B), with Phe/Tyr30 bound at the proximal P27 site instead of the second nutlin binding site ( Figure S5), contrary to the computational models. The presence of a proline in an α-helix is strongly linked to helix kinking, 42 similar to what is observed here. A recent survey of the Protein Data Bank (PDB) 43 reveals that this phenomenon is commonly found in long α-helices, occurring in 30% of membrane protein helices and 20% of length-matched (≥20 residues) globular protein helices. 44 Tyr100 and Tyr104, which line the proximal P27 site, are both rotated toward the peptide, relative to the WT p53-MDM2 complex structure ( Figure 3B). Tyr100 forms an edge-to-face interaction with Phe/Tyr30 and a hydrogen bond with the backbone carbonyl of Leu26 from the peptides, while rotation of Tyr104 allows it to form a hydrogen bond with the hydroxyl of Tyr30 and slightly occlude the putative second nutlin binding site.
The binding of the N-termini of YS-1 and YS-2 (residues 17−23) is very similar to that of two other MDM2-binding stapled peptides, SAH-p53-8 and M06 ( Figure 4B). However, the different stapling strategies resulted in residues 24−28 being bound in a slightly different manner. The introduction of the second stapling point at residue 24 of YS-1 and YS-2 causes Leu25 and Leu26 to be displaced further along the helical axis. This results in a more expanded helical structure with a 1.7 Å displacement at the Cα position of Leu26 (compared to SAH-p53-8), which forces Leu26 into a position and conformation that is midway between that of SAH-p53-8/M06 and WT p53. Residues 27−29 then adopt a tighter helical conformation than in SAH-p53-8, bringing Asn29 of SAH-p53-8 and Gln29 in YS-1 and YS-2 closer together.
MD simulations show that this unexpected binding mode of the stapled peptides is stable ( Figure S6). Although MM/GBSA analysis suggests that adoption of the crystallographic binding mode results in a loss of binding free energy compared to the predicted binding mode, YS-1 and YS-2 were predicted to be tighter MDM2 binders than sMTide-02 (Table S5), in agreement with the FP assay results. The unanticipated binding modes of the stapled peptides observed in the crystal structures highlight the challenges in characterizing multiple conformational states with their distinct thermodynamic profiles. 45 Fine balances between enthalpic and entropic contributions often make it difficult to predict the behavior of the protein and its ligand upon complex formation. 46 The effect of Phe/Tyr30 on the interaction with MDM2 was further investigated by synthesizing two control stapled peptides, one lacking the C-terminal aromatic residue (YS-3) and the other having alanine as the C-terminal residue (YS-4), and then evaluating their binding affinities in competitive FP binding assays. The dissociation constants (K d ) of YS-3 and YS-4 were 5-fold and 4-fold higher than that of YS-2, respectively.

Letter
The slightly lower K d of YS-4 compared to YS-3 could be due to the weak hydrophobic interaction of the methyl side chain of Ala30 with the proximal P27 site. Together with the crystal structures of MDM2 bound to YS-1 and YS-2, these results suggest that the proximal P27 site is a functional binding site that plays an important role in ligand binding to MDM2.
The secondary structures of unbound YS-1 and YS-2 were characterized by circular dichroism spectroscopy. Both peptides show very low overall α-helicity as compared to sMTide-02, based on their molar circular dichroism value at 222 nm ( Figure  S7). To increase the α-helicity of the designed peptides, we replaced Pro27 in YS-1 and YS-2 with a serine residue, generating the respective analogs YS-5 and YS-6. As predicted, YS-5 and YS-6 exhibited a marked increase in α-helicity relative to YS-1 and YS-2 ( Figure S7). However, enhancing peptide αhelicity did not improve the binding affinity (Table 1). This could be due to the loss of favorable hydrophobic interactions of Pro27 with the MDM2 surface when it is replaced by serine, which negates the reduced entropic cost of increasing αhelicity.
In this proof-of-concept study, we have used a probe-based MD method called LMMD to detect novel binding sites on the surface of the anticancer protein target MDM2. Two adjacent putative binding sites on the N-terminal domain of MDM2 close to the p53 binding pocket were identified. Through biophysical binding assays and X-ray crystallography, we serendipitously confirmed the proximal P27 site as a functional ligand binding site by using hydrocarbon stapled peptides that were designed to target the downstream putative site. We measured a 5-fold improvement in binding affinity for stapled peptides that interact with the proximal P27 binding site (YS-2 versus YS-3, Table 1). Because of its proximity to the p53binding cleft, it may be exploited to enhance the binding potencies of current MDM2 ligands. Further optimization of inhibitor structure, such as the use of alternative staple architectures 47 and replacement of Pro27 with helix-stabilizing aliphatic residues that retain hydrophobic interactions with the MDM2 surface, may also result in binding affinity improvements greater than those observed here.
Here, we also present the first ever crystal structures of an i, i + 4 (R,R) stapled peptide bound to its protein target. The structures of these complexes reveal the rare and unusual incorporation of a proline residue into the α-helix of the stapled peptides. This implies that the hydrocarbon staple is able to exert a strong helix stabilization effect that extends to residues beyond the staple, including those with low α-helix propensities such as proline. We also report four modified stapled peptides with MDM2 binding affinities superior to that of the parent peptide sMTide-02, indicating their potential as templates for the development of a new family of potent p53-activating stapled peptides. Studies are now underway to evaluate their biological activities and develop variants with enhanced potencies.
We have used LMMD simulations to identify a novel binding site on the N-terminal domain of MDM2 close to the p53binding cleft. The work described here is a significant step forward for the emerging field of probe-based MD because to the best of our knowledge, this is the first study to utilize such simulations to successfully predict and subsequently validate with definitive biophysical and structural data a previously unknown binding site. The results obtained corroborate the use of probe-based MD techniques to identify novel binding sites and inform structure-based drug design.

* S Supporting Information
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpclett.6b01525.
Computational and experimental details; additional MD simulation details; supplementary tables and figures (PDF)

Accession Codes
The crystal structures of MDM2 in complex with YS-1 and YS-2 were deposited in the PDB under the accession codes 4UE1 and 4UD7, respectively.

Notes
The authors declare no competing financial interest.

■ ACKNOWLEDGMENTS
We are grateful for funding from the A*STAR Joint Council Office. J.R. and C.E.J. are supported by grants from Astex Pharmaceuticals Ltd and Cancer Research UK (grant reference C240/A15751). We thank Diamond Light Source for access to beamlines I04 and I04-1 and Dr. Arnaud Basle for data collection.