Deciphering the Allosteric Effect of Antagonist Vismodegib on Smoothened Receptor Deactivation Using Metadynamics Simulation

Introduction

The Hedgehog (Hh) signaling pathway plays central roles in the animal development and stem-cell function, linking to cell growth, and differentiation, with normal roles in embryonic pattern formation and adult tissue homeostasis and pathological roles in tumor initiation and growth (Pasca di Magliano and Hebrok, 2003; Lum and Beachy, 2004; Rohatgi and Scott, 2007). The smoothened receptor (Smo) serves as one of the key proteins in this signaling pathway, the phosphorylation, stabilization, and accumulate of Smo cause signal to intracellular effectors (Denef et al., 2000; Kalderon, 2000).

Smo has been regarded as an efficacious therapeutic target for basal cell carcinoma (BCC) and medulloblastoma (MB) (Kim et al., 2013). Nevertheless, therapeutic challenges remain in the tumors that acquire resistance to Smo antagonists (Rudin et al., 2009; Yauch et al., 2009; Dijkgraaf et al., 2011; Sharpe et al., 2015a), also in cases of active mutants of Smo that were reported as oncogenic drivers (Reifenberger et al., 1998; Sweeney et al., 2014). So the emergence of resistant mutations and active mutations has raised the need for the discovery of novel Smo inhibitors or novel target sites of Smo.

As a member of the F family of GPCRs, the structure of Smo, so as frizzled (FZD) receptors (Huang and Klein, 2004), contains an extracellular domain (ECD) composed of an extracellular cysteine-rich domain (CRD) and a linker domain (192-222), a seven-transmembrane helical domain (TMD) and an intracellular carboxy-terminal domain (ICD). The CRD and TMD are connected by linker domain (Schulte, 2010; Byrne et al., 2016). Zhang et al. published the structures of multiple domains of Smo and the linker domain called hinge domain (HD) (191-220) (Zhang et al., 2017). Smo possesses two separable ligand-binding sites, one in the TMD and another in the CRD (Sharpe et al., 2015b). The linker domain and TM helices form the deep hydrophobic TMD binding pocket, targeted by allosteric inhibitors and activators, including LY-2940680 (Bender et al., 2011; Wang et al., 2013), cyclopamine (Chen et al., 2002a; Weierstall et al., 2014), SANT-1 and SAG1.5 (Chen et al., 2002b; Miller-Moslin et al., 2009; Wang et al., 2014), and the anticancer drug vismodegib (Robarge et al., 2009; Byrne et al., 2016). Vismodegib has been used to treat advanced basal cell cancer in clinical practice. Another ligand-binding site is located in the hydrophobic groove of the CRD (Nachtergaele et al., 2013; Nedelcu et al., 2013), targeted by cholesterol (Cooper et al., 2003; Huang et al., 2016), 20(S)-hydroxycholesterol (OHC) and other cholesterol analogs (Dwyer et al., 2007; Nachtergaele et al., 2012), that activate Smo in Hh signaling pathway (Nedelcu et al., 2013). Cholesterol was regarded as the endogenous activator of Smo (Huang et al., 2016).

Previous biochemical studies have shown that sterols binding to CRD groove site activate Hh signaling pathway. Binding of an agonist or antagonist to the TMD-site can activate or inactivate Hh signaling pathway. In addition, the vismodegib binding to the TMD-site results in loss of cholesterol from the CRD–linker domain–TMD interface (Byrne et al., 2016). Zhang et al. have demonstrated that TM6, extracellular loop 3 (ECL3) and the HD play a central role in signal transmission, and their structures reveal a precise arrangement of TMD, HD, and CRD. This structure enables allosteric interactions between the three domains that are important to ligand recognition and receptor activation (Zhang et al., 2017). Yet, the detailed mechanisms of how vismodegib allosterically inhibits the activation of Smo and binding of cholesterol remain unknown.

The crystal structures of Smo binding with vismodegib and cholesterol have similar conformation, both in inactive states. We have little known from the static crystal structures about the conformational variation in Smo upon binding different ligands. The computational methods such as molecular dynamics (MD) simulations can provide the information about the dynamic process of conformational changes at the atomic level upon agonist and antagonist binding to GPCRs (McRobb et al., 2016; Miao and McCammon, 2016; Latorraca et al., 2017). The enhanced sampling method, metadynamics (Laio and Parrinello, 2002; Micheletti et al., 2004), supplies an effective and reliable way to explore the binding and unbinding of ligand from GPCRs, and conformational dynamics of GPCRs binding different ligands (Li et al., 2013; Schneider et al., 2015; Saleh et al., 2017). Bai et al performed metadynamics simulation on Smo-vismodegib complex to explore the binding mechanism between the vismodegib and Smo (Bai et al., 2014).

Herein we performed metadynamics simulation to shed light on the mechanism that vismodegib allosterically inhibits the activation of human Smo and binding of cholesterol, by analyzing the synergistic interaction of multiple domains of Smo. Our results revealed that the movements of TM5, TM6, and TM7 induced by vismodegib binding are crucial in deactivation of Smo. In the extracellular side, upon vismodegib binding, the hydrophobic pocket accommodating cholesterol, forming by hydrophobic residues of CRD groove, TM6, ECL3, and HD, rearrange results from the movement of TM6. Hence, the CRD groove takes the open conformational state in the inactivation state, accommodating the TM6 and ECL3. In the intracellular side of inactivation of Smo by vismodegib binding, the movements of TM3, TM5, TM6, and TM7 lead the ICL2 showing closed, ICL3 being open, and the ICL1 and W535 being away from each other. And the HD and ECLs play crucial roles in coordinating the synergistic interaction between the multiple domains of Smo in deactivation. The revealed deactivation mechanism of Smo and conformational changes will be helpful to develop of more effective modulators of Smo or detect potential active site.

Materials and Methods

Preparation of Simulation Systems

As a starting point in the simulations, and as a reference conformation to analyze the results, we used previously determined structures of the human Smo bound to vismodegib and cholesterol (Byrne et al., 2016). The crystal structures were obtained from the PDB database (PDB ID: 5L7I, 5L7D). The chain A of the Smo-cholesterol (hereafter called Smo-CLR) complex was selected, while the chain B of Smo-vismodegib (hereafter called Smo-Vismod) complex was selected. The cytochrome B-562 and solvate molecules of the crystal construct were omitted except water molecules. The missing residues of Smo-Vismod and Smo-CLR were built in the Protein Preparation Wizard in Schrodinger 2015 (Madhavi Sastry et al., 2013). Side chain ionization states were modeled with the PROPKA tool (Søndergaard et al., 2011). The membrane around the transmembrane domain of Smo was built by 85 × 85 Å POPC: cholesterols with 9:1 using CHARMM-GUI webserver (Lee et al., 2016), the receptor crystal structure pre-aligned in the OPM (Orientations of Proteins in Membranes) database (Lomize et al., 2006). Each system was solvated by 12 Å with a truncated rectangular box of TIP3P waters (Jorgensen et al., 1983) and neutralized to a concentration of 0.15 M NaCl.

The proteins were modeled using the AMBER FF99SB force field (Hornak et al., 2006), the ligands were modeled using the generalized AMBER force field (GAFF) (Wang et al., 2004). Geometry optimization and the electrostatic potential calculations on the ligands were performed at the HF/6-31G^* level in the Gaussian09 software (Frisch et al., 2009), and the partial charges were calculated with the RESP (Fox and Kollman, 1998). The force field parameters for the ligands were created by the Antechamber package.

Before metadynamics simulation, the energy minimization and equilibration were conducted by NAMD 2.9 simulation package (Phillips et al., 2005) in order to equilibrate the systems. Firstly, to remove bad contacts in the initial structures, steepest descent was carried out. After energy minimization, each system was gradually heated in NVT ensemble from 0 to 300 K in 300 ps. Subsequently, constant temperature equilibration at 300 K for a total of 5 ns was performed to adjusting the solvent density. Finally, 20 ns conventional molecular dynamic simulations were carried out for each system in NPT ensemble with periodic boundary conditions; an integration step of 2 fs was used. The particle mesh Ewald (PME) algorithm (Darden et al., 1993) was employed to treat long-range electrostatic interactions, while the non-bonded interactions were calculated based on a cutoff of 12 Å. The SHAKE algorithm (Ryckaert et al., 1977) was applied to constrain all covalent bonds involving hydrogen atoms.

Metadynamics Simulations

Metadynamics is an efficient enhanced sampling method, allows the system to escape from local minima in the free energy surface (FES) to explore the conformational space by filling the minima with an external history-dependent bias potential, and permits an accurate determination of the FES (Laio and Parrinello, 2002; Laio and Gervasio, 2008). This bias potential is built as a sum of Gaussians deposited along the trajectory in the pre-defined collective coordinates (CVs) space.

$\begin{array}{ll}V(\vec{s},t)={\displaystyle \sum _{k\tau <t}}W(k\tau )exp\left(-{\displaystyle \sum _{i=1}^d}\frac{{(s_i-s_i(q(k\tau )))}^2}{2{{\sigma }_i}^2}\right)& (1)\end{array}$

where $V(\vec{s},t)$ is the bias potential added to the system, τ is the Gaussian deposition stride, σ_i is the width of the Gaussian for the ith CV, and W(kτ) is the height of the Gaussian. The effect of the metadynamics biased potential is to push the system away from local minima into visiting new regions of the phase space. Furthermore, in the long time limit, the bias potential converges to minus the free energy as the function of the CVs:

$\begin{array}{ll}V(\vec{s},t\to \infty )=-F(\vec{s})+C& (2)\end{array}$

In standard metadynamics, Gaussians of constant height is added for the entire course of a simulation. As a result, the system is eventually pushed to explore high free-energy regions and the estimate of the free energy calculated from the bias potential oscillates around the real value.

In this work, we carried out 100 ns metadynamics simulations for each system, biasing the potential along the following two CVs: the distance between the center mass of W109 and R161, as well as the distance between the center mass of P263 and W535, ignoring hydrogen atoms and labeled as d1, d2, respectively. A Gaussian width of 0.15 Å was used for both CVs, and a Gaussian deposition rate of 0.1 kcalmolps was used.

Results and Discussions

Determination of Conformational States From Free Energy Surface

W109 and R161 of CRD groove located in the opposite place are used to characterize the conformational dynamics of the groove site, and they are involved in binding of sterols (Rana et al., 2013). P263 located in the C-terminal of ICL1 and W535 located in the intercellular tip of TM7. The communication between ICL1 and W535 may trigger Smo activation (Arensdorf et al., 2016), since the activating Smo mutation W535L has been already known in BCCs (Xie et al., 1998; Lam et al., 1999). The free energy surfaces (FES) are shown in Figure 1. The FES along d1 and d2 achieves convergence in last 10 ns for each system (Supplementary Figure 1). As seen in Figure 1A, the FES of Smo-Vismod in the spatial distribution is separated obviously by four energy basins, marked as 1, 2, 3, and 4. While the FES of Smo-CLR owns two energy basins and a minor energy basin, marked as 1, 2, and 3 (Figure 1B). That is, the Smo-Vismod undergoes larger conformational changes compared with Smo-CLR. As a comparison, we marked the location of crystal structures in the FES diagram according to the values of d1 and d2 (Figure 1). Obviously, the crystal structures of Smo-Vismod and Smo-CLR sit among the intermediates and finally stabilize, respectively, in state 1 with the simulations.

FIGURE 1

Figure 1. The two-dimensional maps of the free energy surfaces along d1 and d2 of Smo-Vismod (A) and Smo-CLR (B). The marked black triangles represent the position of d1 and d2 of the crystal structures.

Furthermore, d1 varies widely of Smo-Vismod while d2 varies in a narrow range both in Smo-Vismod and Smo-CLR. We can still differentiate the transformation among states (Figure 1). By combining the FES distribution of Smo-Vismod and Smo-CLR (Figure 1), as well as the crystal structures of Smo bound to different agonist or antagonist (Nachtergaele et al., 2013; Rana et al., 2013; Wang et al., 2014; Weierstall et al., 2014; Byrne et al., 2016; Huang et al., 2016; Zhang et al., 2017) (Supplementary Table 1), we speculated that the CRD groove of the inhibited state is open. ICL1 and TM7 are away from each other, with the corresponding d1 is 23 ± 3 Å and d2 is 10 ± 2 Å. The d1 and d2 shifting to 12 ± 2 Å and 5 ± 1 Å are considered as activated state of Smo, the CRD groove is in a closed conformation where can accommodate cholesterol, and ICL1 and TM7 remain communicating. The d1 varying between 17 ± 3 or 7 ± 2 Å indicates the intermediate states of Smo, meanwhile, the d2 varies between 7 ± 1 Å (Table 1).

TABLE 1

Table 1. The conformational states of Smo were determined by d1 and d2.

The Conformational Difference Between Smo-Vismod and Smo-CLR

As marked in Figure 1, Smo-Vismod and Smo-CLR transform from similar conformation to their respective stable conformations. In order to observe the conformational difference of stable states between Smo-Vismod and Smo-CLR, we extracted the representative structures from the global minimum of FES, respectively. As shown in Figure 2a, the most significant difference is that the CRD tilts to the membrane plane of Smo-CLR, almost forms 60° angle between CRD and the membrane plane. Contrarily, the CRD of Smo-Vismod further moves away from the membrane plane and almost perpendicular to the membrane plane. The distance of the center mass between CRDs of Smo-Vismod and Smo-CLR is 15.6 Å by aligning their TMDs. Compared with their crystal structures (Supplementary Figure 2), the CRD of Smo-CLR slightly tilts to membrane plane related to Smo-Vismod with respect to the aligned TMDs, and the distance of the center mass between CRDs of Smo-Vismod and Smo-CLR is 6.0 Å. This means that their crystal structures are not the stable states, both of the structures undergo substantial conformational rearrangement.

FIGURE 2

Figure 2. (a) The representative structures of Smo-Vismod (forest) and Smo-CLR (marine) at the free-energy minima. (b,c) Hydrophobic interaction formed among CRD groove, HD, TM6 and ECL3 of Smo-CLR and Smo-Vismod respectively; (d) the conformational comparison of CRDs between Smo-Vismod and Smo-CLR; (e) conformational comparison of TMD-sites between Smo-Vismod and Smo-CLR.

To investigate the details of conformational difference between the Smo-Vismod and Smo-CLR, we compared the CRD groove sites and TMD sites of the Smo-Vismod and Smo-CLR (Figures 2b–e). Firstly, we studied the cholesterol binding site, a hydrophobic pocket formed by the hydrophobic residues in CRD (residues V107, L108, L112, I156, V157), HD (residue V210), TM6 (residues V488, L489) and ECL3 (residues V494, I496) (Supplementary Figure S2B). And the hydrophobic pocket is still maintained by cholesterol binding in Smo-CLR (Figure 2b). However, in the stable state of Smo-Vismod, the extracellular extension of TM6 and ECL3 occupy the CRD groove (Figure 2c) and the CRD groove site collapses. The CRD groove, HD, TM6 and ECL3 form strong hydrophobic interactions. Therefore, TM6 and ECL3 occupy the CRD groove in Smo-Vismod, which hinders the cholesterol binding to groove site.

The aligned CRDs of stable states of Smo-Vismod and Smo-CLR (Figure 2d), the CRD groove of Smo-Vismod is open, while the CRD groove of Smo-CLR is closed relatively. However, in the crystal structures, the conformations of CRDs of Smo-Vismod and Smo-CLR are similar (Supplementary Figure S2C), the RMSD is about 0.30 Å, in the closed state. The open conformation of CRD groove induced by the vismodegib bound in our simulation, except the CRD-apo structure of Xenopus laevis Smo has a similar open conformation (Supplementary Table 1) (Huang et al., 2016), exists neither in other CRD crystal structures nor in the multiple-domain crystal structures that resolved previously. Obviously, the CRD groove of inactivated Smo takes the open conformation, which accommodates TM6 and ECL3 by forming strong hydrophobic interactions. In the TMD-site (Figure 2e), the notable variations are the movements of TM5, TM6, and TM7. Compared with Smo-CLR, the TM6 of Smo-Vismod moves outward from the helix bundle, while TM5 and TM7 move inward (Figure 2e).

The Conformational Dynamics of the CRD Groove Site

The Smo-Vismod undergoes significant conformational changes during the simulation (Figures 1, 2). We analyzed the conformational dynamics process of Smo-Vismod from the FES (Figure 1). Based on the hydrophobicity of the CRD groove site, and the interaction among the hydrophobic residues of CRD groove, HD, TM6, and ECL3, the solvent-accessible surface (SASA) of the hydrophobic pocket composed of them is observed over time (Figures 3C,D), along with the evolution of d1 (Figures 3A,B) to study the conformational dynamics of the ECD of Smo-Vismod.

FIGURE 3

Figure 3. The evolution of d1 and d2 along time of Smo-Vismod (A) and Smo-CLR (B). The solvent accessible surface area (SASA) of the hydrophobic pockets consisted of hydrophobic residues of CRD groove, HD, helix VI and ELC3 evolve with time for Smo-Vismod (C) and Smo-CLR (D).

At the ~30 ns of the start of the simulation, the d1 is <15 Å (Figure 3A), and the SASA of the hydrophobic pocket during the simulation is not significant change (Figure 3C), indicating that the CRD groove is still in the closed state (Supplementary Table 1). The corresponding energy basin is 4 (Figure 1A). The representative conformation of the energy basin 4 was superimposed with the crystal structure of the Smo-Vismod (Figure 4A). As seen that the CRD groove is in a closed state as in the initial state. The extracellular extension of TM6 and ECL3 do not occupy the CRD groove, nor form hydrophobic interaction. Compared to the initial structure, TM6 and ECL3 have begun to approach the groove.

FIGURE 4

Figure 4. (A,B) The conformational dynamics of ECDs of basins 4 (orange), 3 (yellow-orange), 2 (light green), and 1 (forest) of Smo-Vismod compared with the crystal structure (white). (C) The conformational changes of the CRD. (D) The conformational changes of HD.

The second obvious fluctuation is from ~30 to ~50 ns, the d1 increases rapidly and fluctuates in the range of ~15 to ~20 Å. The corresponding SASA of the hydrophobic pocket decreases rapidly, indicating that the CRD groove is beginning to open. The appropriate energy basin is 3 (Figure 1A). As shown in Figure 4A, we found that TM6 and ECL3 enter into the CRD groove. The hydrophobic pockets formed among the hydrophobic residues of CRD, HD, TM6, and ECL3 begin to be closed. Several hydrophilic residues (E158-G162) of the side of CRD groove begin to leave the groove site and hence the CRD groove gets gradually open.

At the end of 50 ns, d1 increases and fluctuates between 20 and 25 Å during the subsequent simulation (Figure 3A). The SASA of the hydrophobic pocket also fluctuates at a steady level in the last 50 ns (Figure 3C), indicating that the hydrophobic pocket maintains in compact contact. The representation of energy basins are 1 and 2 (Figure 1A). Compared the representative structures of the two energy basins (Figure 4B), the TM6 and ECL3 have completely occupied the CRD groove and form strong hydrophobic interactions with the CRD groove. CRD groove presents a significant open conformation. With the TM6 and ECL3 of Smo-Vismod gradually approach to CRD groove, CRD is away from the membrane plane (Figure 4C), the HD moves toward the direction of TM6 (Figure 4D). The synergistic movement of the extracellular domain of Smo-Vismod contributes to stabilize CRD in inactivated conformation.

In contrast, the d1 smoothly fluctuates at the beginning of ~25 ns; after ~25ns, mainly fluctuates between 10 and 15 Å (Figure 3B); the SASA shows increased at beginning ~ 20 ns (Figure 3D) in Smo-CLR. We compared the representative structures of Smo-CLR in the energy basins (Supplementary Figure S3). It is clear that CRD tilts toward the membrane plane during the simulations, the TM6 shifts toward the flank of CRD, simultaneously, the ECL3 deviates from CRD. And the HD moves away from TM6 (Supplementary Figure S3A). However, the CRD groove does not change significantly (Figure S3B). Noticeably, the activated conformational dynamics of ECD is in contrast to the inhibited conformation of Smo-Vismod.

At this point, we have proved that the movement of TM6 and the hydrophobic residues of CRD groove, HD, TM6, and ECL3 play crucial roles in the deactivation of Smo by vismodegib binding. Indeed, the TM6 and ECL3 occupy the CRD groove in Smo-Vismod, which hinders the cholesterol binding to groove site.

The Conformational Dynamics of the TMD-Site

As observed above, the movement of TM6 is crucial in the conformational dynamics of the CRD groove site. Considering the binding mode of vismodegib that directly interacts with TM5, TM6, and TM7, we firstly studied the interaction between vismodegib and TMD-site. As shown in Figure 5A, we observed that the chlorophenyl–methylsulfone moiety of vismodegib overturns nearly 90° in the stable inactivated state and closes to HD-ECL2 compared with the crystal structure. The atom O3 of the methylsulfone moiety forms H-bond with atom ND2 of side-chain of N219. The RMSD of vismodegib was calculated throughout the simulation (Figure 5B). We caught sight of the obvious remolding, that is, the RMSD suddenly increases from ~0.4 to ~1.4 Å at the beginning ~5 ns, and fluctuates stably between 1.2 and 1.6 Å, indicating that the binding mode of vismodegib remodels immediately during the simulation. We also monitored the distance between atom O3 of the methylsulfone moiety and atom ND2 of N219 (Figure 5C). The distance rapidly drops to ~3 Å from around 5.2 Å at the beginning simulation of ~5 ns and eventually runs aground. This means that the vismodegib rapidly leaves the initial binding mode shortly after the start of the simulation, afterwards, stabilizes in the new binding mode, forming steady H-bond with side-chain of N219, which also indicates the importance of HD.

FIGURE 5

Figure 5. (A) The binding mode of vismodegib of the represent structure in free energy minima (forest) varied compared with the crystal structure (white). (B) The RMSD of heavy atoms of vismodegib evolved along with time. (C) The distance between the atom O3 of vismodegib and the atom ND2 of N219 evolved along with time.

At the other side of TMD-site, the charged residue R400 of TM5, D473 of TM6, and E518 of TM7 contribute to the binding of vismodegib (Wang et al., 2014; Byrne et al., 2016). The movement of TM5, TM6, and TM7 is obvious in the stable states of Smo-Vismod and Smo-CLR (Figure 2e). Therefore, we analyzed the conformational dynamics of the three polar residues (Figure 6). As seen in Figure 6A, the distance between R400 and E518 of Smo-Vismod is much smaller than Smo-CLR throughout the simulation. This indicates that the two residues going close to each other at the beginning of the simulation and form stable electrostatic interaction Smo-Vismod, which results in TM5 and TM7 moving inward (Figures 6B,F). TM6 is extruded from the helix bundle (Figure 6F). At the same time, the two residues going close to each other results in a remolding binding mode of vismodegib. Since in the crystal structure, the negatively charged E518 is near the amide oxygen atom of vismodegib with the distance of 3.1 Å from OE2 of E518 (Supplementary Figure S4), the electrostatic repulsion causes overturning of vismodegib.

FIGURE 6

Figure 6. The conformational changes of TMD-site. (A) The center mass distance between the side chain of R400 and E518 evolved along with time in the Smo-Vismod (forest) and Smo-CLR (marine). (B) The distance of the extracellular end of TM5 (residues 397-401) and TM7 (residues 515-519) evolved along with time. (C) The scatter plot of dihedral of C, CA, CB, and CG of D473 between the Smo-Vismod and Smo-CLR. (D) The distance of the extracellular end of TM3 (residues 315-319) and TM6 (residues 488-492) evolved along with time. (E) The conformational changes of R400, D473 and E518. (F) The movements of TM5, TM6, and TM7.

In addition, due to R400 and E518 going close to each other, the electrostatic repulsion of the E518 and the electrostatic attraction of the R400 result in the D473 deflected (Figures 6C,E). The dihedral angles of the atoms C, CA, CB, and CG of the D473 of the Smo-Vismod is far less than the Smo-CLR (Figure 6C). D473 points toward R400 in the representative state of Smo-Vismod (Figure 6E) instead of pointing toward E518 in the crystal structure (Supplementary Figure S4).

From these results, the remodeled binding model of vismodegib is stabilized by forming H-bonds with N219 of HD, D384, and Y394 of ECL2, respectively. Additionally, R400, E518, and D473 constitute the electrostatic interface to form polar interaction with amide linker of vismodegib (Figures 5A, 6E). R400, D473, and E518 play a vital role in the conformational dynamics of helices. Moreover, the outward movement of TM6 makes the extracellular extension of TM6 close to the CRD groove, and thus along with ECL3 interact with hydrophobic residues of CRD groove and HD, which hinders the binding of the cholesterol. On the contrary, cholesterol occupied the groove to push the ECL3 away from the groove and keep the TM6 out of the groove, allowing the CRD to tilt toward membrane plane.

The Conformational Dynamics of the ICD

Compared with d1, the fluctuation of d2 is not significant throughout the simulation in both of systems (Figures 3A,B). Nevertheless, three major energy basins of Smo-Vismod along d2, termed to 4, 2-3 and 1, were obtained (Figure 1A). We compared the representative structures of the intracellular domain of Smo-Vismod and Smo-CLR, the ICLs give notable variations (Figure 7). Compared with Smo-CLR, the intracellular end of TM7 of Smo-Vismod shifts outward. Therefore, there is no communication between ICL1 and W535, which indicates an inhibited state of Smo. Furthermore, the significant conformational changes are the ICLs. As shown in Figure 7A, the intracellular end of TM3 shifts inward and ICL2 covers over the central interface of the helix bundle showing closed in the four representative structures of Smo-Vismod. The ICL3 is opening with the simulation due to the shifts of intracellular tips of TM5 and TM6. While, the ICL2 is open, and the ICL3 is closed in the Smo-CLR (Figure 7B). The importance of ICLs of Smo in regulating Smo signaling was confirmed by peptide mimics of ICL2 and ICL3 that suppress Smo induced tumor cell proliferation (Remsberg et al., 2007). In our simulations, the closed ICL2 and the open ICL3 are only seen in the Smo-Vismod, which proves that the intracellular loops potentially affect the interaction of downstream effectors and Smo. Compared the representative conformations of Smo-Vismod and Smo-CLR, we have got an outline of inhibited Smo induced by vismodegib. The vismodegib binding to TMD site leads to the movement of TM6, which promotes the extracellular extension of TM6 and ECL3 to occupy the groove by forming strong hydrophobic interactions (Figure 2c). Hence the vismodegib binding blocks the binding of cholesterol. At the intracellular side, the ICL2 shows closed owing to the shifts of TM3. ICL3 is open caused by the shifts of TM5 and TM6. The ICL1 and W535 are away from each other due to the TM7 shifts outward. The interaction of the multi-domain of Smo induced by vismodegib hinders the binding of the cholesterol, and destroys the intracellular interactional interface of the effector.

FIGURE 7

Figure 7. (A) The conformational dynamics of ICDs of basins 4 (orange), 3 (yellow-orange), 2 (light green), and 1 (forest) in the Smo-Vismod. (B) The conformational dynamics of ICDs of basins 1 (marine) and 2 (cyan) in Smo-CLR.

The Communication of Multi-Domain of Smo

We have made it clear so far that the antagonist vismodegib binding leads to the movements of TM5, TM6, and TM7, and then stabilizes Smo in an inactive state through coordinated movement between multiple domains. We further carried out the dynamical network analysis (Sethi et al., 2009) to investigate the communication of multi-domain of Smo in the inactive state and active state. As shown in Figure 8, the HD, ECLs spread in multiple communities in both of the systems, indicating these domains play a key role in coordinating the interaction of CRD and TMD. Which owe to the inherent flexibility and naturally ingenious arrangement of these domains. And the Smo-Vismod identifies more communities, and the communities more weakly connect each other compared with Smo-CLR, suggesting that the construction of Smo is looser induced by vismodegib. While the Smo-CLR takes a compact and ordered construction ensuring the extracellular signal transmits to the intracellular side and permitting the downstream effectors to contact.

FIGURE 8

Figure 8. The communication between the multiple domains of Smo-Vismod (A) and Smo-CLR (B).

We have demonstrated the vismodegib binding triggers the shifts of TM5, TM6, and TM7, leading the extracellular extension of TM6 and ECL3 to enter the CRD groove and block the binding of cholesterol. Simultaneously, the significant conformational changes of ICLs are attributed to the movements of TM3, TM5, TM6, and TM7. The communication of multi-domain of Smo proves that HD and ECLs play roles in the synergistic movements of CRD and TMD. The looser construction and weak communication provide the view of the inactive state of Smo.

Conclusions

In this work, the mechanism of the vismodegib allosterically inhibits Smo activation and hinders cholesterol binding is revealed by metadynamics simulations. We revealed that the vismodegib binding leads movements of TM5, TM6, and TM7, and the shift of TM6 triggers the entrance extracellular extension of TM6 and ECL3. Therefore, the TM5, TM6, and TM7 are key factors to the deactivation of Smo upon vismodegib binding. Moreover, we also found an inhibited open conformation of CRD groove, which is not shown in crystal structures. The open CRD groove accommodates the TM6 and ECL3 so that hinders cholesterol to bind and holds the CRD stacked atop the TMD. The strong hydrophobic interaction of CRD groove, HD, TM6 and ECL3 stabilizes the interaction of them. Furthermore, the HD and ECLs play a key role in the coordinated interaction of CRD, TMD and themselves. Therefore, blocking the coordinated movement of CRD, HD, and ECLs may potentially inhibit the activation of Smo. And we observed a remolding of vismodegib binding, stabilized by hydrogen bond formed with N219, D384 and Y394. These results can be taken into account for the design and discovery of novel Smo inhibitors, or providing structural information for discovering potential active sites.

Author Contributions

XA, QB, FB, and DS designed the research, performed the dynamic simulations and analyzed the data. XA and XY wrote the manuscript. HL and XY directed the project.

Funding

National Natural Science Foundation of China (Grant Nos. 2175060, 21605066).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fchem.2019.00406/full#supplementary-material

References

Arensdorf, A. M., Marada, S., and Ogden, S. K. (2016). Smoothened regulation: a tale of two signals. Trends Pharmacol. Sci. 37, 62–72. doi: 10.1016/j.tips.2015.09.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Bai, Q., Shen, Y., Jin, N., Liu, H., and Yao, X. (2014). Molecular modeling study on the dynamical structural features of human smoothened receptor and binding mechanism of antagonist LY2940680 by metadynamics simulation and free energy calculation. Biochim. Biophys. Acta 1840, 2128–2138. doi: 10.1016/j.bbagen.2014.03.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Bender, M. H., Hipskind, P. A., Capen, A. R., Cockman, M., Credille, K. M., Gao, H., et al. (2011). Identification and characterization of a novel smoothened antagonist for the treatment of cancer with deregulated hedgehog signaling. Cancer Res. 71, 2819–2819. doi: 10.1158/1538-7445.AM2011-2819

CrossRef Full Text | Google Scholar

Byrne, E. F. X., Sircar, R., Miller, P. S., Hedger, G., Luchetti, G., Nachtergaele, S., et al. (2016). Structural basis of Smoothened regulation by its extracellular domains. Nature 535, 517–522. doi: 10.1038/nature18934

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, J. K., Taipale, J., Cooper, M. K., and Beachy, P. A. (2002a). Inhibition of Hedgehog signaling by direct binding of cyclopamine to Smoothened. Genes Dev. 16, 2743–2748. doi: 10.1101/gad.1025302

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, J. K., Taipale, J., Young, K. E., Maiti, T., and Beachy, P. A. (2002b). Small molecule modulation of Smoothened activity. Proc. Natl. Acad. Sci. U.S.A. 99, 14071–14076. doi: 10.1073/pnas.182542899

PubMed Abstract | CrossRef Full Text | Google Scholar

Cooper, M. K., Wassif, C. A., Krakowiak, P. A., Taipale, J., Gong, R., Kelley, R. I., et al. (2003). A defective response to Hedgehog signaling in disorders of cholesterol biosynthesis. Nat. Genet. 33, 508–513. doi: 10.1038/ng1134

PubMed Abstract | CrossRef Full Text | Google Scholar

Darden, T., York, D., and Pedersen, L. (1993). Particle mesh Ewald: An N·log(N) method for Ewald sums in large systems. J. Chem. Phys. 98, 10089–10092. doi: 10.1063/1.464397

CrossRef Full Text | Google Scholar

Denef, N., Neubüser, D., Perez, L., and Cohen, S. M. (2000). Hedgehog Induces opposite changes in turnover and subcellular localization of patched and smoothened. Cell 102, 521–531. doi: 10.1016/S0092-8674(00)00056-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Dijkgraaf, G. J., Alicke, B., Weinmann, L., Januario, T., West, K., Modrusan, Z., et al. (2011). Small molecule inhibition of GDC-0449 refractory smoothened mutants and downstream mechanisms of drug resistance. Cancer Res. 71, 435–444. doi: 10.1158/0008-5472.CAN-10-2876

PubMed Abstract | CrossRef Full Text | Google Scholar

Dwyer, J. R., Sever, N., Carlson, M., Nelson, S. F., Beachy, P. A., and Parhami, F. (2007). Oxysterols are novel activators of the hedgehog signaling pathway in pluripotent mesenchymal cells. J. Biol. Chem. 282, 8959–8968. doi: 10.1074/jbc.M611741200

PubMed Abstract | CrossRef Full Text | Google Scholar

Fox, T., and Kollman, P. A. (1998). Application of the RESP Methodology in the Parametrization of Organic Solvents. J. Phys. Chem. B. 102, 8070–8079. doi: 10.1021/jp9717655

CrossRef Full Text | Google Scholar

Frisch, M. J., Trucks, G. W., Schlegel, H. B., Scuseria, G. E., Robb, M. A., Cheeseman, J. R., et al. (2009). Gaussian 09. Wallingford, CT: Gaussian, Inc.

Hornak, V., Abel, R., Okur, A., Strockbine, B., Roitberg, A., and Simmerling, C. (2006). Comparison of multiple Amber force fields and development of improved protein backbone parameters. Proteins 65, 712–725. doi: 10.1002/prot.21123

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, H. C., and Klein, P. S. (2004). The Frizzled family: receptors for multiple signal transduction pathways. Genome Biol. 5:234. doi: 10.1186/gb-2004-5-7-234

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, P., Nedelcu, D., Watanabe, M., Jao, C., Kim, Y., Liu, J., et al. (2016). Cellular cholesterol directly activates smoothened in hedgehog signaling. Cell 166:e1114. doi: 10.1016/j.cell.2016.08.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983). Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 79, 926–935. doi: 10.1063/1.445869

CrossRef Full Text | Google Scholar

Kalderon, D. (2000). Transducing the hedgehog signal. Cell 103, 371–374. doi: 10.1016/S0092-8674(00)00129-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, J., Aftab, B. T., Tang, J. Y., Kim, D., Lee, A. H., Rezaee, M., et al. (2013). Itraconazole and arsenic trioxide inhibit Hedgehog pathway activation and tumor growth associated with acquired resistance to smoothened antagonists. Cancer Cell 23, 23–34. doi: 10.1016/j.ccr.2012.11.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Laio, A., and Gervasio, F. L. (2008). Metadynamics: a method to simulate rare events and reconstruct the free energy in biophysics, chemistry and material science. Rep. Prog. Phys. 71:126601. doi: 10.1088/0034-4885/71/12/126601

CrossRef Full Text | Google Scholar

Laio, A., and Parrinello, M. (2002). Escaping free-energy minima. Proc. Natl. Acad. Sci. U.S.A. 99, 12562–12566. doi: 10.1073/pnas.202427399

PubMed Abstract | CrossRef Full Text | Google Scholar

Lam, C. W., Xie, J., To, K. F., Ng, H. K., Lee, K. C., Yuen, N. W., et al. (1999). A frequent activated smoothened mutation in sporadic basal cell carcinomas. Oncogene 18, 833–836. doi: 10.1038/sj.onc.1202360

PubMed Abstract | CrossRef Full Text | Google Scholar

Latorraca, N. R., Venkatakrishnan, A. J., and Dror, R. O. (2017). GPCR Dynamics: structures in motion. Chem. Rev. 117, 139–155. doi: 10.1021/acs.chemrev.6b00177

PubMed Abstract | CrossRef Full Text | Google Scholar

Lee, J., Cheng, X., Swails, J. M., Yeom, M. S., Eastman, P. K., Lemkul, J. A., et al. (2016). CHARMM-GUI Input Generator for NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM simulations using the CHARMM36 additive force field. J. Chem. Theory Comput. 12, 405–413. doi: 10.1021/acs.jctc.5b00935

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, J., Jonsson, A. L., Beuming, T., Shelley, J. C., and Voth, G. A. (2013). Ligand-dependent activation and deactivation of the human adenosine A(2A) receptor. J. Am. Chem. Soc. 135, 8749–8759. doi: 10.1021/ja404391q

PubMed Abstract | CrossRef Full Text | Google Scholar

Lomize, M. A., Lomize, A. L., Pogozheva, I. D., and Mosberg, H. I. (2006). OPM: orientations of proteins in membranes database. Bioinformatics 22, 623–625. doi: 10.1093/bioinformatics/btk023

PubMed Abstract | CrossRef Full Text | Google Scholar

Lum, L., and Beachy, P. A. (2004). The hedgehog response network: sensors, switches, and routers. Science 304, 1755–1759. doi: 10.1126/science.1098020

PubMed Abstract | CrossRef Full Text | Google Scholar

Madhavi Sastry, G., Adzhigirey, M., Day, T., Annabhimoju, R., and Sherman, W. (2013). Protein and ligand preparation: parameters, protocols, and influence on virtual screening enrichments. J. Comput. Aided Mol. Des. 27, 221–234. doi: 10.1007/s10822-013-9644-8

PubMed Abstract | CrossRef Full Text | Google Scholar

McRobb, F. M., Negri, A., Beuming, T., and Sherman, W. (2016). Molecular dynamics techniques for modeling G protein-coupled receptors. Curr. Opin. Pharmacol. 30, 69–75. doi: 10.1016/j.coph.2016.07.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Miao, Y., and McCammon, J. A. (2016). G-protein coupled receptors: advances in simulation and drug discovery. Curr. Opin. Struct. Biol. 41, 83–89. doi: 10.1016/j.sbi.2016.06.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Micheletti, C., Laio, A., and Parrinello, M. (2004). Reconstructing the density of states by history-dependent metadynamics. Phys. Rev. Lett. 92:170601. doi: 10.1103/PhysRevLett.92.170601

PubMed Abstract | CrossRef Full Text | Google Scholar

Miller-Moslin, K., Peukert, S., Jain, R. K., McEwan, M. A., Karki, R., Llamas, L., et al. (2009). 1-Amino-4-benzylphthalazines as orally bioavailable smoothened antagonists with antitumor activity. J. Med. Chem. 52, 3954–3968. doi: 10.1021/jm900309j

PubMed Abstract | CrossRef Full Text | Google Scholar

Nachtergaele, S., Mydock, L. K., Krishnan, K., Rammohan, J., Schlesinger, P. H., Covey, D. F., et al. (2012). Oxysterols are allosteric activators of the oncoprotein Smoothened. Nat. Chem. Biol. 8, 211–220. doi: 10.1038/nchembio.765

PubMed Abstract | CrossRef Full Text | Google Scholar

Nachtergaele, S., Whalen, D. M., Mydock, L. K., Zhao, Z., Malinauskas, T., Krishnan, K., et al. (2013). Structure and function of the Smoothened extracellular domain in vertebrate Hedgehog signaling. eLife 2:e01340. doi: 10.7554/eLife.01340

PubMed Abstract | CrossRef Full Text | Google Scholar

Nedelcu, D., Liu, J., Xu, Y., Jao, C., and Salic, A. (2013). Oxysterol binding to the extracellular domain of Smoothened in Hedgehog signaling. Nat. Chem. Biol. 9, 557–564. doi: 10.1038/nchembio.1290

PubMed Abstract | CrossRef Full Text | Google Scholar

Pasca di Magliano, M., and Hebrok, M. (2003). Hedgehog signalling in cancer formation and maintenance. Nat. Rev. Cancer 3, 903–911. doi: 10.1038/nrc1229

PubMed Abstract | CrossRef Full Text | Google Scholar

Phillips, J. C., Braun, R., Wang, W., Gumbart, J., Tajkhorshid, E., Villa, E., et al. (2005). Scalable molecular dynamics with NAMD. J. Comput. Chem. 26, 1781–1802. doi: 10.1002/jcc.20289

PubMed Abstract | CrossRef Full Text | Google Scholar

Rana, R., Carroll, C. E., Lee, H. J., Bao, J., Marada, S., Grace, C. R., et al. (2013). Structural insights into the role of the Smoothened cysteine-rich domain in Hedgehog signalling. Nat. Commun. 4:2965. doi: 10.1038/ncomms3965

PubMed Abstract | CrossRef Full Text | Google Scholar

Reifenberger, J., Wolter, M., Weber, R. G., Megahed, M., Ruzicka, T., Lichter, P., et al. (1998). Missense mutations in SMOH in sporadic basal cell carcinomas of the skin and primitive neuroectodermal tumors of the central nervous system. Cancer Res. 58, 1798–1803.

PubMed Abstract | Google Scholar

Remsberg, J. R., Lou, H., Tarasov, S. G., Dean, M., and Tarasova, N. I. (2007). Structural analogues of smoothened intracellular loops as potent inhibitors of hedgehog pathway and cancer cell growth. J. Med. Chem. 50, 4534–4538. doi: 10.1021/jm0705657

PubMed Abstract | CrossRef Full Text | Google Scholar

Robarge, K. D., Brunton, S. A., Castanedo, G. M., Cui, Y., Dina, M. S., Goldsmith, R., et al. (2009). GDC-0449-a potent inhibitor of the hedgehog pathway. Bioorg. Med. Chem. Lett. 19, 5576–5581. doi: 10.1016/j.bmcl.2009.08.049

PubMed Abstract | CrossRef Full Text | Google Scholar

Rohatgi, R., and Scott, M. P. (2007). Patching the gaps in Hedgehog signalling. Nat. Cell Biol. 9, 1005–1009. doi: 10.1038/ncb435

PubMed Abstract | CrossRef Full Text | Google Scholar

Rudin, C. M., Hann, C. L., Laterra, J., Yauch, R. L., Callahan, C. A., Fu, L., et al. (2009). Treatment of medulloblastoma with hedgehog pathway inhibitor GDC-0449. N. Engl. J. Med. 361, 1173–1178. doi: 10.1056/NEJMoa0902903

PubMed Abstract | CrossRef Full Text | Google Scholar

Ryckaert, J. P., Ciccotti, G., and Berendsen, H. J. C. (1977). Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys. 23, 327–341. doi: 10.1016/0021-9991(77)90098-5

CrossRef Full Text | Google Scholar

Saleh, N., Ibrahim, P., Saladino, G., Gervasio, F. L., and Clark, T. (2017). An efficient metadynamics-based protocol to model the binding affinity and the transition state ensemble of G-protein-coupled receptor ligands. J. Chem. Inf. Model. 57, 1210–1217. doi: 10.1021/acs.jcim.6b00772

PubMed Abstract | CrossRef Full Text | Google Scholar

Schneider, S., Provasi, D., and Filizola, M. (2015). The dynamic process of drug-GPCR binding at either orthosteric or allosteric sites evaluated by metadynamics. Methods Mol. Biol. 1335, 277–294. doi: 10.1007/978-1-4939-2914-6_18

PubMed Abstract | CrossRef Full Text | Google Scholar

Schulte, G. (2010). International union of basic and clinical pharmacology. LXXX. The class Frizzled receptors. Pharmacol. Rev. 62, 632–667. doi: 10.1124/pr.110.002931

PubMed Abstract | CrossRef Full Text | Google Scholar

Sethi, A., Eargle, J., Black, A. A., and Luthey-Schulten, Z. (2009). Dynamical networks in tRNA:protein complexes. Proc. Natl. Acad. Sci. U.S.A. 106, 6620–6625. doi: 10.1073/pnas.0810961106

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharpe, H. J., Pau, G., Dijkgraaf, G. J., Basset-Seguin, N., Modrusan, Z., Januario, T., et al. (2015a). Genomic analysis of smoothened inhibitor resistance in basal cell carcinoma. Cancer Cell 27, 327–341. doi: 10.1016/j.ccell.2015.02.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharpe, H. J., Wang, W., Hannoush, R. N., and de Sauvage, F. J. (2015b). Regulation of the oncoprotein Smoothened by small molecules. Nat. Chem. Biol. 11, 246–255. doi: 10.1038/nchembio.1776

PubMed Abstract | CrossRef Full Text | Google Scholar

Søndergaard, C. R., Olsson, M. H., Rostkowski, M., and Jensen, J. H. (2011). Improved treatment of ligands and coupling effects in empirical calculation and rationalization of pKa values. J. Chem. Theory Comput. 7, 2284–2295. doi: 10.1021/ct200133y

PubMed Abstract | CrossRef Full Text | Google Scholar

Sweeney, R. T., McClary, A. C., Myers, B. R., Biscocho, J., Neahring, L., Kwei, K. A., et al. (2014). Identification of recurrent SMO and BRAF mutations in ameloblastomas. Nat. Genet. 46, 722–725. doi: 10.1038/ng.2986

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, C., Wu, H., Evron, T., Vardy, E., Han, G. W., Huang, X. P., et al. (2014). Structural basis for Smoothened receptor modulation and chemoresistance to anticancer drugs. Nat. Commun. 5:4355. doi: 10.1038/ncomms5355

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, C., Wu, H., Katritch, V., Han, G. W., Huang, X. P., Liu, W., et al. (2013). Structure of the human smoothened receptor bound to an antitumour agent. Nature 497, 338–343. doi: 10.1038/nature12167

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, J., Wolf, R. M., Caldwell, J. W., Kollman, P. A., and Case, D. A. (2004). Development and testing of a general amber force field. J. Comput. Chem. 25, 1157–1174. doi: 10.1002/jcc.20035

PubMed Abstract | CrossRef Full Text | Google Scholar

Weierstall, U., James, D., Wang, C., White, T. A., Wang, D., Liu, W., et al. (2014). Lipidic cubic phase injector facilitates membrane protein serial femtosecond crystallography. Nat. Commun. 5:3309. doi: 10.1038/ncomms4309

PubMed Abstract | CrossRef Full Text | Google Scholar

Xie, J., Murone, M., Luoh, S. M., Ryan, A., Gu, Q., Zhang, C., et al. (1998). Activating Smoothened mutations in sporadic basal-cell carcinoma. Nature 391, 90–92. doi: 10.1038/34201

PubMed Abstract | CrossRef Full Text | Google Scholar

Yauch, R. L., Dijkgraaf, G. J., Alicke, B., Januario, T., Ahn, C. P., Holcomb, T., et al. (2009). Smoothened mutation confers resistance to a hedgehog pathway inhibitor in medulloblastoma. Science 326, 572–574. doi: 10.1126/science.1179386

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, X., Zhao, F., Wu, Y., Yang, J., Han, G. W., Zhao, S., et al. (2017). Crystal structure of a multi-domain human smoothened receptor in complex with a super stabilizing ligand. Nat. Commun. 8:15383. doi: 10.1038/ncomms15383

PubMed Abstract | CrossRef Full Text | Google Scholar