- 1The Key Laboratory of Natural Medicine and Immuno-Engineering, Henan University, Kaifeng, China
- 2School of Biological and Chemical Engineering, Zhejiang University of Science and Technology, Hangzhou, China
Butyrylcholinesterase (BChE) is a non-specific enzyme with clinical pharmacological and toxicological significance, which was a renewed interest as therapeutic target in Alzheimer's disease (AD) nowadays. Here, all-atom molecular dynamics simulations of butyrylcholinesterase with tacrine complex were designed to characterize inhibitor binding modes, strengths, and the hydrogen-bond dependent non-covalent release mechanism. Four possible release channels were identified, and the most favorable channel was determined by random acceleration molecular dynamics molecular dynamics (RAMD MD) simulations. The thermodynamic and dynamic properties as well as the corresponding Detour-forward delivery mechanism were determined according to the classical molecular dynamics (MD) simulations accompanied with umbrella sampling. The free energy barrier of the tacrine release process for the most beneficial pathway is about 10.95 kcal/mol, which is related to the non-covalent interactions from the surrounding residues, revealing the specific binding characteristics in the active site. The residues including Asp70, Ser79, Trp82, Gly116, Thr120, Tyr332, and His438 were identified to play major roles in the stabilization of tacrine in the pocket of BChE, where hydrogen bonding and π-π interactions are significant factors. Tyr332 and Asp70, which act as gate keepers, play crucial roles in the substrate delivery. The present results provide a basic understanding for the ligand transport mechanism depending on the BChE enzymatic environment, which is useful for the design of BChE inhibitors in the future.
Introduction
Alzheimer's disease (AD), a progressive neurodegenerative disorder, is characterized by generalized dementia such as memory impairment, aphasia, apraxia, visual spatial skill impairment, executive dysfunction, etc. (Tschanz and Andersen, 2011). With the aging of the population, age-dependent neurodegenerative disorders have become a serious medical problem in modern society. Thus, tremendous material and financial resources are gradually devoted into treatment of AD. Whereas, the etiology of AD is not completely understood due to its being complex and multifactorial, low level of acetylcholine (ACh), β-amyloid (Aβ) aggregation, tau-protein hyperphosphorylation, oxidative stress, etc. (Grundke-Iqbal et al., 1986; Talesa, 2001; Hardy, 2009; Rosini et al., 2014) all have been correlated at some extent to AD.
Cholinesterases (ChEs), including acetylcholinesterase (AChE) and butyrylcholinesterase (BChE), play important roles in cholinergic transmission by hydrolysis of the neurotransmitter acetylcholine (ACh) (Masson et al., 1999). Many studies have found that AD is related to the decrease in acetylcholine (ACh), a neurotransmitter in the brain (Anand and Singh, 2013a; Sun et al., 2015; Zhu et al., 2017). Actually, the shortfall in cholinergic activity affects the activity of AChE and BChE, both of which regulate acetylcholine in the human brain (Giacobini, 2003). It is worth noting that the level of AChE in the brain declines to 55–67%, while BChE increases to 120% of normal levels during the progressed AD, indicating that BChE plays a critical role for ACh hydrolysis in the late stage of AD (Greig et al., 2001; Mushtaq et al., 2014). Hence, the BChE is considered to be a potentially useful target (Giacobini, 2003; Hartmann et al., 2007; Nordberg et al., 2013). The importance of BChE has also been demonstrated by the AChE knockout mice model, in which BChE compensates for the lack of AChE, maintaining normal cholinergic pathways in AChE nullizygous animals (Xie et al., 2000; Mesulam M. M. et al., 2002). It has been discovered that BChE can replace AChE to hydrolyze ACh in the presence of a specific AChE inhibitor in the human brain (Mesulam M. et al., 2002). In recent years, increasing the level of ACh in the brain by inhibiting the biological activity of AChE is a significant therapeutic approach for the treatment of AD. Accordingly, inhibition of BChE can be viewed as an alternative for the restoration of cholinergic activity and improvement of cognitive performance for patients with AD (Mushtaq et al., 2014).
Nowadays, there are five anti-Alzheimer's disease drugs approved for clinical application (Figure 1), namely (a) tacrine, (b) galantamine, (c) memantine, (d) donepezil, and (e) rivastigmine, respectively. Among them, tacrine is a reversible inhibitor and a classical pharmacophore as the first generation of anti-AD agent, which was approved by the Food and Drug Administration (FDA) in 1993 (Summers, 2006). However, it was withdrawn (Qizilbash et al., 1999) due to its clinical shortcomings instigating hepatotoxicity (Lagadic-Gossmann et al., 1998). In spite of severe hepatotoxicity, tacrine is still considered as a good scaffold to design multitarget drugs for Alzheimer's disease (Bartolini and Marco-Contelles, 2019). Thus, the development of tacrine derivatives has drawn immense attention since tacrine is a classical pharmacophore and inhibits both AChE and BChE at micromolar scale (Sameem et al., 2017). Plenty of tacrine derivatives have been developed with the aim of improving inhibitory activity, eliminating toxicity, and enlarging properties as cholinesterase inhibitors (ChEIs) (Anand and Singh, 2013b). More recently, a new strategy has been proposed with the aim of designing compounds that can simultaneously target different targets that are closely related to AD. The compounds are referred to as multitarget-directed ligands (MTDLs), and tacrine is preferred in the designing of MTDLs (Bajda et al., 2011; Zhang et al., 2015; Lin et al., 2017; Chalupova et al., 2019). Moreover, the selective or non-selective inhibition of BChE may raise neuroprotective and disease-modifying effects (Greig et al., 2005). Therefore, the structural data of tacrine complex with BChE is extremely importance for developing the hopeful ChEIs, and it is recently resolved at 2.1 Å resolution (PDB ID: 4BDS) (Nachon et al., 2013). The X-ray cocrystal structure is shown in Figure 2. As we can see, key interactions mainly come from (1) the aromatic π-π stacking with Trp82; (2) hydrogen bonding between the N7 (as shown in Figure 1A) and the main chain carbonyl of His438; and (3) two hydrogen bonds between the amino group of tacrine and two water molecules in the surrounding water molecular network (Nachon et al., 2013).
Figure 1. Major drugs approved by the US Food and Drug Administration (FDA) for treatment of Alzheimer's disease (AD). (A) Tacrine, (B) galantamine, (C) memantine, (D) donepezil, and (E) rivastigmine.
Figure 2. X-ray crystal structure of BChE complexed with tacrine (PDB ID: 4BDS). Key residues and ligand are represented as sticks, and waters are shown as turquoise spheres.
The statistic crystal structure provides the intuitive structural model, while it is also significant to understand the interactions between BChE and tacrine, and several crucial questions are still open: (i) What is the specific contribution of key residues for tacrine binding in the active site? (ii) What are possible channels for the substrate tacrine delivery and which one is the most advantageous? (iii) How about the thermodynamic and dynamics properties and mechanisms during the tacrine release process? (iv) Which residues play crucial roles in the tacrine delivery process? Here, in order to clarify these issues, classical molecular dynamics simulations and extensive combined random acceleration molecular dynamics molecular dynamics (RAMD MD) simulations combined with the umbrella sampling technique have been performed. Such efforts can be served as reasonable initial stages to explore novel potential specific BChE inhibitors, which provide a molecular basis for designing new BChE inhibitors with better kinetic properties selectivity. The modeling of tacrine entrance into the active site is difficult to build for its uncertain positions outside the protein. Moreover, according to our previous studies, the process of ligand entrance into the active site is similarly reversible to the ligand release in a certain extent, especially for the possible pathways, delivery mechanisms, and conformation changes of protein (Zhao et al., 2016). Hence, here, we mainly focus on the tacrine release process from the active site to outside of BChE.
Computational Methods and Details
Setup of the Enzyme–Substrate Complex Model
In the present study, the initial model for the simulations was acquired from the high-resolution crystal structures of BChE in complex with tacrine (PDB ID: 4BDS; resolution: 2.1 Å) (Nachon et al., 2013). The missing fragment of 378–379 residues was complemented by the biopolymer module in the (Sybyl-X 2.1, 2013) software. The proton state of ionizable residues was determined at pH 6.5 by PROPKA 3.0 (Rostkowski et al., 2011) as well as their surrounding environment. The partial atomic charges of tacrine were derived by the restrained electrostatic potential (RESP) (Cornell et al., 1993) charge at the HF/6-31G* level with Gaussian03 software (Frisch et al., 2004). The protein and water molecules were described by the ff14SB force field (Maier et al., 2015; Shao and Zhu, 2018) and the TIP3P model, respectively (Jorgensen et al., 1983). The complex was dissolved in a cubic water box of 86 × 87 × 102 Å with a buffer distance of 10 Å on each side between the box wall and the nearest solute atoms, and then, the whole model box was neutralized with chlorine ions. The tleap modules in the AMBER 16 (Case et al., 2016) software was used to generate the topology parameters and initial coordinates of the whole system. The process of the model building is depicted in Figure 3.
Figure 3. The structure of the butyrylcholinesterase (BChE) protein and the water box of the building model.
Molecular Mechanical Molecular Dynamics Simulations
In order to remove poor interactions between atoms, the whole system was minimized for 10,000 steps by using the steepest descent method, followed by a conjugate gradient method for another 10,000 steps minimization. Afterwards, the system was heated from 0 to 300 K gradually for 100 ps, and then, another 100 ps NPT MD simulations were carried out to relax the density to ~1.0 g/cm3. Finally, 100 ns classical MD simulations with NPT ensemble by using periodic boundary condition were performed. Langevin method was employed to maintain the temperature at 300 K (Pastor et al., 1988). Finally, 100 ns classical MD simulations with NPT ensemble by using periodic boundary condition were performed, which adopt the nearest mirror image principle. The 12 Å cutoff value was set to calculate the van der Waals force and the electrostatic interaction given the desired accuracy and computing costs, which has been successfully used in several other approximations even in bigger systems (Zhao et al., 2018; Fan et al., 2019). All of the bonds with hydrogen atoms in the model were constrained with the SHAKE algorithm (Ryckaert et al., 1977). The stability of the backbone of BChE was evaluated using root-mean-square deviation (RMSD) analysis. Then, the hydrogen bonds were analyzed using the last 10 ns trajectories of the MD simulations. The whole MD simulations were performed using the AMBER 16 package (Case et al., 2016).
Molecular Mechanics/Generalized Born Surface Area
Molecular mechanics/generalized born surface area (MM/GBSA) (Kollman et al., 2000; Wang et al., 2001, 2006) has been widely applied to calculate the free energy for various protein–ligand (Kuhn and Kollman, 2000; Huo et al., 2002; Brown and Muchmore, 2006; Hou and Yu, 2007), protein–protein, and protein–peptide (Wang and Kollman, 2000; Gohlke and Case, 2004; Hou et al., 2006; She et al., 2019; Wei et al., 2019) complexes successfully. The binding free energy between the ligand and receptor was calculated with the following equation. Moreover, MM/GBSA allows for rigorous free energy decomposition into contributions originating from different groups of atoms or types of interaction (Gohlke et al., 2003; Hou et al., 2008, 2011).
where ΔGcomplex, ΔGprotein, and ΔGligand stand for the free energies of complex, receptor, and ligand, respectively (Srivastava and Sastry, 2012; Munnaluri et al., 2015). The free energy also can be estimated as the sum of three terms:
where ΔEmm is the total gas phase energy expressed as the sum of the internal (int), electrostatic (ele), and van der Waals (vdW) energies:
where ΔGsol accounts for the polar (ΔGGB) and non-polar (ΔGSA) solvation energies. ΔGGB is the polar contribution to the solvation free energy, described by the Generalized Born (GB) calculation, while ΔGSA, the non-polar solvation energy, is calculated from the solvent accessible surface area (SASA). TΔS is the conformational entropy of binding and usually calculated by the normal-mode analysis. The aim of this method is to further understand the effect of each residue on the release process of tacrine, and the main results show that van der Waals, electrostatic interaction, and non-polar solvent effects provide the main driving force for substrate binding, while polar solvent energy plays a negative role.
RAMD MD Simulations for Detecting Tacrine Release Channels
In order to identify the possible channels for ligand transportation from active site to the outside of the protein, the combined random acceleration molecular dynamics and MD simulations (RAMD MD) (Lüdemann et al., 2000; Vashisth and Abrams, 2008) have been carried out by using NAMD 2.9 software (Phillips et al., 2005). Proteins and ligands were described with ff14SB (Maier et al., 2015; Shao and Zhu, 2018) and GAFF force fields (Wang et al., 2004), respectively. The initial structure of the BChE–tacrine system was derived from the equilibrium state of MM MD simulations. In the RAMD MD approach, a small randomly oriented force is added to the center of the mass of tacrine to identify the possible channels of the ligand release. The direction of the random force is kept for a certain period of time. During this time, if tacrine displacement reaches the threshold parameter, the direction of the force is maintained, otherwise, a new direction will be selected randomly. When the ligand escapes from its initial position, the classic MD simulation will be performed to balance the sampling, which avoided the error arising from the higher random force. In the present work, random accelerations of 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40, 0.45, and 0.50 kcal/Å/g and thresholds of 0.5, 0.4, 0.3, and 0.2 Å were applied to the initial model, respectively. Totally, 144 RAMD MD trajectories were acquired to identify possible delivery channels of tacrine by statistics analysis.
Umbrella Sampling
Based on the statistics of the RAMD MD simulations, four possible channels were generated. The umbrella sampling (US) technique (Kästner, 2011) was applied to map out the free energy profiles of the most favorable pathway to further identify the thermodynamic and dynamics properties. In this process, the conformational changes of the protein can also be observed. The initial model was constructed based on the stable snapshot structure in the conventional MD simulations. Based on the most probable channel determined from RAMD MD simulations, the distance between the C10 atom of tacrine and the Cα of Ile442 was chosen as the reaction coordinate (RC) for the substrate transportation (see Figure 4), which varies from 12.5 to 32.5 Å with a 0.5-Å interval for two adjacent windows. Afterwards, for each window, 40 ns MM MD simulations with an appropriate biasing harmonic potential were performed along the reaction coordinate to ensure a reasonable overlap for all windows. The last 7 ns trajectory data from the all windows were analyzed by the weighted histogram analysis method (WHAM) to generate the potential of mean force (PMF) (Kumar et al., 1992). Interactions between the ligand and the protein as well as the role of the key residues were investigated.
Figure 4. Definitions of reaction coordinate for the tacrine delivery. RC, the distance between the C10 of tacrine and the Cα of Ile442.
Results and Discussions
Equilibrium Structure Analysis
In order to investigate the stability of the enzyme–substrate complex structure for subsequent detection of possible substrate release pathways, long time MM MD simulations have been carried out. On the basis of the root-mean-square deviation (RMSD) results of the backbone atoms, the system approaches equilibrium after 40 ns (see Figure S1), which suggests that the interaction between the ligand and BChE reaches an equilibrium state at that time. Thus, the last 20 ns sampling snapshots have been used for the following analysis. The representative snapshots from the final MD simulation are displayed in Figure 5, which shows that tacrine has been immobilized in the active site of BChE. At that time, strong face-to-face π-π stacking interactions from the aromatic planes of Trp82 have been discovered to stable the substrate in the active site, which plays an essential role in the whole binding model. As shown in Figure 5, it is found that not only a face-to-face π-π stacking interactions from the aromatic planes of Trp82 but also an edge-to-face π-π stacking was formed between the benzene ring from tacrine and Tyr332 immediately. Furthermore, there existed a series of sound van der Waals interaction between the tacrine and the residues of Asp70, Trp82, Gly115, Gly116, Gly117, Thr120, Ala328, Phe329, Tyr332, and His438, which play crucial roles in stabilizing the substrate in the active site. In addition, two localized water molecules (WatA and WatB) in the binding gorge region have been observed to form hydrogen bonds with the amino group of the aromatic ring. Meanwhile, we found that another water molecule (WatC) also plays a crucial role in protein–tacrine association by mediating interactions between the aromatic nitrogen and the main chain carbonyl of His438. Beyond that, a significant number of water molecules have been observed to stay close to the tacrine binding site and form a network of hydrogen bonds to assist the tacrine stabilization. Therefore, the tacrine was stabilized by those abundant hydrogen bonds and strong aromatic stacking aromatic stacking. After 100 ns simulations, the interatomic distances between key residues and tacrine are collected, which is listed in Figure S2, showing good consistence with the X-ray structure. Moreover, the protein and location of key residues in substrate binding, such as Trp82, Gly116, Ser198, Pro285, Ala328, Phe329, Tyr332, and His438, are also extracted to compare with that in the crystal structure by overlap technology. As shown in Figure 6, both protein and active site observed in theoretical calculations and experiment are overlapping well, which means that our calculated model is reasonable.
Figure 5. The equilibrium structure of butyrylcholinesterase (BChE)–tacrine complex in the active site after molecular mechanics–molecular dynamics (MM MD) simulations. Key residues are represented as sticks; all oxygen and nitrogen atoms are colored red and blue, respectively. Carbons are colored pink in the ligand and green in residues of the complex. Waters are represented by red spheres, hydrogen bonds are drawn as black dashes, and aromatic stacking are drawn as yellow dashes.
Figure 6. Overlap of initial crystal structure and stable state after molecular mechanics–molecular dynamics (MM MD) simulations. Tacrine in initial crystal structure is represented as green sticks; in MD simulations, structure is represented as orange sticks.
Free Energy Decomposition Identifies “Efficient Amino Acids”
The MM/GBSA method has been used to analyze the contribution of all amino acids to the binding free energy of BChE–tacrine complex systematically by employing 500 snapshots from the last 5 ns in MM MD simulations. The total binding free energy of the complex was −16.90 kcal/mol that is constituted by the integrated effects of van der Waals (ΔEvdW), electrostatic (ΔEele), non-polar solvation (ΔGSA), and polar solvation interaction (ΔGGB) (see Table S2), and the ΔEvdW is the most contributor that is benefit for tacrine binding. Moreover, the binding energy contributions of individual amino acids were also calculated, as we can see in Figure 7A, the residues supported interactions with the ligand directly and outstanding contributions to the binding energy have been labeled specifically. We divide van der Waals (ΔEvdW), electrostatic (ΔEele), non-polar solvation (ΔGSA), and polar solvation interaction (ΔGGB) into ΔGnon−pol and ΔGpol. ΔGnon−pol equals to the sum of ΔEvdw and ΔGSA, and ΔGpol is comprised of the electrostatic contribution (ΔEele) and polar solvation interaction (ΔGGB). Meanwhile, the energy decomposition of these main residues including Trp82, Gly116, Ser198, Pro285, Ala328, Phe329, Tyr332, and His438 has been depicted in Figure 7B. ΔGnon−pol of these main residues are −2.93, −0.82, −0.43, −0.78, −0.74, −1.33, −0.97, and −1.22 kcal/mol; ΔGpol of these main residues are 0.74, 0.17, −0.08, 0.01, 0.25, 0.28, 0.49, and 0.42 kcal/mol. As demonstrated in Table 1, for these main residues, van der Waals, electrostatic, and non-polar solvated interactions of BChE–tacrine are favorable in the binding of the complex. Polar interaction is 15.16 kcal/mol, which consists of electrostatic interaction and polar solvation effect, playing an unfavorable role for binding. However, it can be counteracted by non-polar interaction (−32.06 kcal/mol) that is much favorable for the binding, which comprises of van der Waals and non-polar solvation energy. Therefore, the BChE and tacrine can be bound with each other at last considering these above factors. Among them, Trp82, Tyr332, and His438 have a significant contribution to the binding. Trp82 exhibited a particularly high correlation with tacrine especially, which can be justified by the fact that this residue forms parallel strong π-π stacking with tacrine, and the van der Waals interactions and non-polar solvent effects are favorable for the substrate binding. Tacrine forms an edge-to-face π-π stacking with Tyr332. Meanwhile, His438 provides a hydrogen bond with tacrine molecule to stabilize it in the active site. The equilibrium structure shows that Gly115, Gly116, and Gly117 also form a hydrogen bond network around the ligand, which is the main source of their binding energy contribution. At the same time, Asp70 forms a water-mediated hydrogen bond with the ligand. These interactions can be also acquired clearly by snapshots sampling in different windows during the release process of tacrine. In addition, the possibility of hydrogen bonds forming between tacrine and BChE as well as surrounding water molecules is also analyzed in detail by the last 10 ns MD simulation trajectory, as shown in Table 2 and Figure 8. The serial number of atoms in tacrine is listed in Figure S3. The probability of hydrogen bonding between Tyr332 and Asp70 is 68.82%. There is a water-mediated hydrogen bond between tacrine and His438 with the probability of 28.78 and 15.35% for each of them, and Gly116 and Thr120 form a hydrogen bond network around the ligand with the probability of 58.22 and 59.87%, which means that they also make a contribution to stabilize the ligand.
Figure 7. (A) Contribution of each residue to butyrylcholinesterase (BChE)–tacrine binding. (B) The binding free energy decomposition of the main residues in the BChE–tacrine systems.
Table 1. The mean and standard deviation (M ± SD) of binding free energy decomposition of key residues in butyrylcholinesterase (BChE)–tacrine complexes (kcal/mol).
Table 2. Hydrogen bonds formed between tacrine and butyrylcholinesterase (BChE) and the corresponding percent analysis with the last 10 ns trajectories in the MD simulations.
Figure 8. The hydrogen bond network in the active site and the probability of bonding are labeled in red.
Furthermore, to identify the role of the specific residues, four key residues of Trp82, Pro285, Phe329, and His438 are selected to take site-directed mutations of alanine, since they contributed the most to the entire substrate and protein binding process according to the MM/GBSA calculation results. A single virtual alanine scanning calculations was made for each of them by 100 snapshots from last 5 in 100 ns MM MD simulations for wild system to verify their contribution to protein–ligand complex binding. These mutant systems are named as Trp82Ala, Pro285Ala, Phe329Ala, and His438Ala, respectively. Detailed results are summarized in Table 3 and Table S3. More importantly, if we carefully analyze all data presented there, we could find that the van der Waals interactions between ligand and protein dominate most to the binding free energy. To illuminate this, the comparison between calculated GvdW and Gtotal is plotted in Figure 9. Table 3 described the difference value of binding free energy between the wide and mutation type. Here, Trp82 is the biggest contributor among the four residues involved in the binding of ligand through the strong π-π stacking. Compared to the wild type, the energies of Trp82Ala decrease to −13.48 kcal/mol, 4.43 kcal/mol less than the wild type. For Phe329, coming in the second, the energies decrease to −17.32 kcal/mol. His438 decreases the binding energy of 2.30 kcal/mol after mutating to Ala residue. It can be found that there is an aromatic stacking with Trp82 and a weak hydrogen bonding between the N7 atom of tacrine and His438 (3.2 Å) experimentally (Nachon et al., 2013). Therefore, our calculated results are in good agreement with the experimental data, suggesting that residues Trp82 and His438 contribute significantly to binding of tacrine. The mutation of these residue leads to a significant decrease in binding affinity. In addition, the mutations of Pro285Ala also cause a tiny reduction in binding energy. It is probably due to the reason that this residue does not interact directly with ligands and would be neglected by the calculation error. Moreover, there still exist a weak interaction between the mutant alanine and the ligand. If one residue was changed, the interactions around tacrine may be changed to some extent.
Table 3. The mean and standard deviation (M ± SD) of binding free energy of wild and mutant systems (kcal/mol).
Figure 9. Comparison of ΔGvdW and ΔGtotal for Ser79Ala, Trp82Ala, Thr120Ala, Tyr332Ala, and His438Ala.
Transport Mechanism of Tacrine Delivery
According to RAMD MD simulations, 144 trajectories were obtained, and four possible pathways were identified, named as P1 (among helix 1 and loop 1), P2 (between helix 1 and loop 3), P3 (through loops 1 and 2), and P4 (between loops 2 and 3) in Figure 10 and Table 4. The possibilities of these channels are 61.11, 29.86, 6.25, and 2.78% for P1, P2, P3, and P4, respectively. As we can see from Figure 10, there is a cavity in the P1 direction, which will facilitate the release of the ligand. Therefore, as we can see in Table 4, path P1, located between helix 1 and loop 1, is recognized as the predominant pathway for the release of substrate, which will be discussed in detail.
Figure 10. Definition of the possible pathways P1, P2, P3, and P4 for the substrate release based on the random acceleration molecular dynamics molecular dynamics (RAMD MD) simulations.
Table 4. The location, number, and possibility of trajectories for tacrine delivery pathways based on random acceleration molecular dynamics molecular dynamics (RAMD MD) simulations.
According to the escaping direction of tacrine for the most favorable channel P1, the distance between the C10 atom of tacrine and the Cα atom of Ile442 was selected as the reaction coordinate (RC) (see Figure 4). In order to obtain a more reliable substrate release mechanism and visualized potential of the mean force (PMF), a total of 1,640 ns MD simulations combined with the umbrella sampling technique along the reaction coordinate from 12.5 to 32.5 Å were carried, and 41 windows were collected. For each window, a series of biasing harmonic potential of 30 kcal/mol was added to ensure full sampling and enough overlapping between neighboring windows (Figure S4). Figure S5 shows the fluctuations of the backbone of BChE, indicating that the trajectories reach a steady state after 20 ns. Moreover, in order to verify the convergence of PMF curves (Magistrato et al., 2017), different time periods and districts are considered (see Figures S6, S7), different time periods and districts are considered (see Figures S6, S7), which implies that the free energy profile probably reach convergence after 25 ns. Therefore, the energy profiles of different sampling time durations of 26–40, 28–40, 27–39, and 26–38 ns have been compared in Figure 11. Clearly, the relative free energy profiles for different sampling time durations are essentially the same with the maximal SD of 0.38 kcal/mol, showing reliable convergence of MD simulations for the PMF calculations. Based on the free energy and the conformational changes of residues and proteins shown in Figures 11–13, the escape process of substrates from the proteins is divided into four consecutive stages.
Figure 12. Structure of the protein changes along P1 in the four phases of the tacrine release determined by umbrella sampling simulations.
Figure 13. Key residues and their conformations in the active site from different windows along the reaction coordinate. Key residues and waters are represented as sticks. Hydrogen bonds are drawn as black dashes; aromatic stacking are drawn as blue dashes.
In the first stage (12.5 Å ≤ RC < 16.0 Å), the delivery channel is totally closed (shown in Figure 14). At 13.0 Å, Tyr332 forms a hydrogen bond with Asp70 to assist the “gate” totally close, and a stable face-to-face π-π stacking formed between tacrine and Trp82 indole ring in the active site. At the same time, the amino group of tacrine forms two water-mediated hydrogen bonds with the water. In addition, van der Waals interactions can be detected from the main chains of Thr120 and Gly115. Both hydrogen bond and van der Waals interactions provide the necessary force to stabilize tacrine at this stage. Then, the hydrogen bond between Tyr332 and Asp70 is broken, causing the “gate” beginning to open; the side chain of Asp70, one of the door keepers, provides a driving force through a new water-mediated hydrogen bond with the substrate at 14.5 Å of reaction coordinate, resulting in the tacrine beginning to move to the “door” as seen in Figure 13. In the second stage (16.0 Å ≤ RC < 20.0 Å), with the reaction coordinate changing, the substrate tries to break away from the active site, leading the water-mediated hydrogen bond between tacrine and Asp70 broken as shown in RC = 16.5 Å and forming a new aromatic hydrogen bond between the amino and Tyr332. We can observe that Tyr332, as another door keeper, turns on one side gradually as a “swinging gate.” At this stage, these two door keepers (Asp70 and Tyr332) interact with the substrate successively and get the hydrogen bond between Tyr332 and Asp70 broken, causing the “gate” beginning to open. At the same time, a new water-mediated hydrogen bond is formed between the carbonyl of Phe329 and substrate, which also provides a driving force for tacrine delivery. Then, water provides a driving force for the substrate escape through a new hydrogen bond. Meanwhile, the “swinging gate” gets back to the original “close” state; the substrate moves toward the protein “cavity” where the number of residues around tacrine is obviously reduced. In the third stage (20.0 Å ≤ RC < 23.5 Å), the newborn hydrogen bond with Tyr332 is broken, following the formation of a water-mediated hydrogen bond with Tyr332 for a short time, and the aromatic nitrogen atom of the tacrine forms a water-mediated hydrogen bond with Ser287. Meanwhile, the main chain of Pro285 forms two hydrogen bonds with tacrine. Then, the conjugate ring flips and tries to move toward outside. The amino group of tacrine forms a strong hydrogen bond with Asn289; the aromatic nitrogen atom of the tacrine forms a hydrogen bond with the water. This flipping resulted to the fluctuation in the third stage in the free energy profiles along the reaction coordinate. At this stage, tacrine tries to escape from the active site; Asn289, closer to the entrance of the pocket, provides a sound driving force by a hydrogen bond. Especially, Asn289 forms a hydrogen bond with tacrine directly. The specific protein environment makes it difficult for the substrate opening the “swinging gate” and escaping from the pocket facilely. Thus, the detour-forward strategy was adopted for the substrate to release. In the last stage (RC ≥ 23.5 Å), tacrine gets rid of the constraint of the hydrogen bonding with Asn289 and gradually flees away from the constraint of BChE. The interactions around tacrine are obviously reduced and weaken, and thus, the tacrine fluctuates stronger, which caused its conformation changes apparently. The inhibitor looks like traveling a long distance in small difference between two adjacent windows. For example, at 28.0 Å, a water bridge has been identified between Pro281 and N7 of tacrine, and at 28.5 Å, a direct hydrogen bond has been found between Asn68 and –NH2 of tacrine. After 32.5 Å, the tacrine escapes from BChE completely, which floats freely that totally exposed to the water environment.
To give a clearer indication of the release process of tacrine and the opening and closing process of the “swinging gate,” we made the illustration of the “gate” during the substrate release, as shown in Figure 14. In addition, to make the release more intuitive, only the key residues, Asp70 and Tyr332, acted as gatekeepers as shown in Figure 14. The conceptual graphs revealing the release process of tacrine are shown Figure 12. As we can see from Figure 14, in the first stage (12.5 Å ≤ RC < 16.0 Å), Tyr332 forms a stable hydrogen bond with Asp70, and the “swinging gate” keeps the tight closing state. The substrate begins to move outward, trying to break away from the active site. Tacrine approaches the “swinging gate” guided by the doorkeeper Asp70, causing the hydrogen bond between two doorkeepers (Asp70 and Tyr332) broken. In the second stage (16.0 Å ≤ RC < 20.0 Å), the “swinging gate” totally open, although that is does not keep a long time. The substrate moves out of the protein under the guidance of doorkeeper Asp70. Then, the “door” gets back to the “close” state. In the last two stages (RC ≥ 20.0 Å), tacrine has reached the cavity outside the “door.” The “door” no longer constrains the release of tacrine, and the “swinging gate” rotates back to the original position.
As shown in Figure 11, this release process is exothermic, while there is a need to overcome an energy barrier of 10.95 kcal/mol. First, the barrier keeps arising from the 13.0 Å of the reaction coordination until almost close to the top at 23.0 Å. The previous three stages are a series of energy increasing steps over constantly. In the first stage, tacrine began to move toward the exit and given rise to the π-π stacking with Trp82 and the hydrogen bond with His438 broken simultaneously. Thus, this is probably the main cause of the barrier when the reaction coordination changes from the 13.0 to 15.0 Å. However, the broken hydrogen bond interaction with Asp70 and Tyr332 may be the primary reason of the risen energy at the second stage. When the reaction coordination moves to the third stage, the free energy fluctuates from 6.2 to 10.1 kcal/mol and reaches the peak of the profiles. It probably comes from tacrine flipping; the flipping leads to the breaking of all non-bond interactions. As to the last stage, tacrine gets rid of the constraint of the hydrogen bonding with Asn289 and Asn68. Thus, there generated an unimpeded path for the substrate delivery, which is probably the main reason of the fluctuation in the free energy profile in Figure 11. Generally, according to the different stages in the delivery mechanism and their thermodynamic properties, the tacrine escaping from BChE will need to overcome a barrier of ~1,0.95 kcal/mol, and π-π stacking interactions with Trp82 are the major factors that affect the overall ligand transportation process. Besides, the hydrogen bond interactions with the residues around the ligand surely affect the transportation process as well.
Conclusion
In this work, the binding characteristics of BChE and antagonist tacrine complex as well as the release mechanism of ligand have been studied by using extensive MM MD simulations combined with several methods and technologies such as MM/GBSA, RAMD, umbrella sampling, etc., and some significant findings were obtained. First, Trp82, Tyr332, His438, and water bridges have been observed to play a crucial role in substrate binding by hydrogen bonds and π-π stacking interactions, which can be also verified by mutate calculations of key residues. Moreover, the non-polar interaction and side chain effects are discovered as the main contributing factors for tacrine and BChE binding. Second, four possible channels have been identified, and the most favorable channel P1 that located around helix 1 and loop 1 has been chosen to describe the mechanism and thermodynamic properties of substrate delivery process. The barrier for the tacrine release is measured to be around 10.57–10.95 kcal/mol, which mainly comes from the flip of the substrate as well as the broken π-π stacking interactions between Trp82 and tacrine.
Moreover, the hydrogen bond formed between Tyr332 and tacrine also contributes to the barrier. Besides that, the substrate release is a hydrogen-bond-dependent process, especially for the hydrogen bond interactions coming from Tyr332 and Asp70, which are the main factors that influence the delivery of tacrine. It is an intricate and flexural process for the substrate to escape from the deep pocket; thus, a detour-forward strategy was adopted to cross the “swinging gate” and release. Lastly, the key residues Tyr332 and Asp70 forms a “swinging gate,” which shows a dynamical switch mechanism, probably dominating the transportation of tacrine. In summary, the results of this work elucidate the antagonist-binding mechanism in the pocket of BChE on the basis of tacrine and supply detailed information about its dissociation pathway, mechanism, and thermodynamic properties at the atomic level. Nowadays, many derivatives are designed based on tacrine, which have good performance in inhibitor activity (Santos et al., 2016; Jerabek et al., 2017; Okten et al., 2019). Therefore, we hope that this work could provide more dynamic information to design novel BChE inhibitors with better kinetic properties by improving the drug residence time, which will be useful for carrying out rational clinical drug treatment for AD.
Data Availability Statement
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.
Author Contributions
ZZ and FF performed calculations. ZZ, FF, and YZ analyzed the data. WL, YZ, and CW designed the research. ZZ, FF, YZ, and CW wrote the paper. All authors have approved the final version of the manuscript.
Funding
This work was supported by the National Science Foundation of China (Nos. 21603057 and U1704176), the Project funded by the China Post-doctoral Science Foundation (Nos. 2017M622324 and 2018T110721), and the Project funded by the Henan Post-doctoral Science Foundation (No. 001702017).
Conflict of Interest
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.
Acknowledgments
We thank Prof. Zexing Cao at Xiamen University for his kind help. We also thank the Supercomputer Center in Wuhan University and National Supercomputer Center in Changsha and Guangzhou for providing the computational resources.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fchem.2020.00730/full#supplementary-material
References
Anand, P., and Singh, B. (2013a). Flavonoids as lead compounds modulating the enzyme targets in Alzheimer's disease. Med. Chem. Res. 22, 3061–3075. doi: 10.1007/s00044-012-0353-y
Anand, P., and Singh, B. (2013b). A review on cholinesterase inhibitors for Alzheimer's disease. Arch. Pharm. Res. 36, 375–399. doi: 10.1007/s12272-013-0036-3
Bajda, M., Guzior, N., Ignasik, M., and Malawska, B. (2011). Multi-target-directed ligands in Alzheimer's disease treatment. Curr. Med. Chem. 18, 4949–4975. doi: 10.2174/092986711797535245
Bartolini, M., and Marco-Contelles, J. (2019). Tacrines as therapeutic agents for Alzheimer's disease. IV. The tacripyrines and related annulated tacrines. Chem. Rec. 19, 927–937. doi: 10.1002/tcr.201800155
Brown, S. P., and Muchmore, S. W. (2006). High-throughput calculation of protein–ligand binding affinities: modification and adaptation of the MM-PBSA protocol to enterprise grid computing. J. Chem. Inf. Model. 46, 999–1005. doi: 10.1021/ci050488t
Case, D. A., Betz, R. M., Cerutti, D. S., Cheatham, T. E., Darden, T. A., Duke, R. E., et al. (2016). AMBER 16. San Francisco, CA: University of California.
Chalupova, K., Korabecny, J., Bartolini, M., Monti, B., Lamba, D., Caliandro, R., et al. (2019). Novel tacrine-tryptophan hybrids: multi-target directed ligands as potential treatment for Alzheimer's disease. Eur. J. Med. Chem. 168, 491–514. doi: 10.1016/j.ejmech.2019.02.021
Cornell, W. D., Cieplak, P., Bayly, C. I., and Kollman, P. A. (1993). Application of RESP charges to calculate conformational energies, hydrogen bond energies, and free energies of solvation. J. Am. Chem. Soc. 115, 9620–9631. doi: 10.1021/ja00074a030
Fan, F., Zhao, Y., and Cao, Z. (2019). Insight into the delivery channel and selectivity of multiple binding sites in bovine serum albumin towards naphthalimide-polyamine derivatives. Phys. Chem. Chem. Phys. 21, 7429–7439. doi: 10.1039/C9CP00527G
Frisch, M. J., Trucks, G. W., Schlegel, H. B., Scuseria, G. E., Robb, M. A., Cheeseman, J. R., et al. (2004). Gaussian 03 Rev. CO2. Wallingford, CT: Gaussian Inc.
Giacobini, E. (2003). Cholinesterases: new roles in brain function and in Alzheimer's disease. Neurochem. Res. 28, 515–522. doi: 10.1023/A:1022869222652
Gohlke, H., and Case, D. A. (2004). Converging free energy estimates: MM-PB(GB)SA studies on the protein–protein complex ras–raf. J. Comput. Chem. 25, 238–250. doi: 10.1002/jcc.10379
Gohlke, H., Kiel, C., and Case, D. A. (2003). Insights into protein-protein binding by binding free energy calculation and free energy decomposition for the ras-raf and Ras-RalGDS complexes. J. Mol. Biol. 330, 891–913. doi: 10.1016/S0022-2836(03)00610-7
Greig, N. H., Utsuki, T., Ingram, D. K., Wang, Y., Pepeu, G., Scali, C., et al. (2005). Selective butyrylcholinesterase inhibition elevates brain acetylcholine, augments learning and lowers Alzheimer β-amyloid peptide in rodent. Proc. Natl. Acad. Sci. U.S.A. 102, 17213–17218. doi: 10.1073/pnas.0508575102
Greig, N. H., Utsuki, T., Yu, Q.-S., Zhu, X., Holloway, H. W., Perry, T., et al. (2001). A new therapeutic target in Alzheimer's disease treatment: attention to butyrylcholinesterase. Curr. Med. Res. Opin. 17, 159–165. doi: 10.1185/03007990152673800
Grundke-Iqbal, I., Iqbal, K., Tung, Y. C., Quinlan, M., Wisniewski, H. M., and Binder, L. I. (1986). Abnormal phosphorylation of the microtubule-associated protein tau (tau) in Alzheimer cytoskeletal pathology. Proc. Natl. Acad. Sci. U.S.A. 83, 4913–4917. doi: 10.1073/pnas.83.13.4913
Hardy, J. (2009). The amyloid hypothesis for Alzheimer's disease: a critical reappraisal. J. Neurochem. 110, 1129–1134. doi: 10.1111/j.1471-4159.2009.06181.x
Hartmann, J., Kiewert, C., Duysen, E. G., Lockridge, O., Greig, N. H., and Klein, J. (2007). Excessive hippocampal acetylcholine levels in acetylcholinesterase-deficient mice are moderated by butyrylcholinesterase activity. J. Neurochem. 100, 1421–1429. doi: 10.1111/j.1471-4159.2006.04347.x
Hou, T., Chen, K., McLaughlin, W. A., Lu, B., and Wang, W. (2006). Computational analysis and prediction of the binding motif and protein interacting partners of the Abl SH3 domain. PLoS Comput. Biol. 2:e1. doi: 10.1371/journal.pcbi.0020001
Hou, T., Wang, J., Li, Y., and Wang, W. (2011). 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. 51, 69–82. doi: 10.1021/ci100275a
Hou, T., and Yu, R. (2007). Molecular dynamics and free energy studies on the wild-type and double mutant HIV-1 protease complexed with amprenavir and two amprenavir-related inhibitors: mechanism for binding and drug resistance. J. Med. Chem. 50, 1177–1188. doi: 10.1021/jm0609162
Hou, T., Zhang, W., Case, D. A., and Wang, W. (2008). Characterization of domain–peptide interaction interface: a case study on the amphiphysin-1 SH3 domain. J. Mol. Biol. 376, 1201–1214. doi: 10.1016/j.jmb.2007.12.054
Huo, S., Wang, J., Cieplak, P., Kollman, P. A., and Kuntz, I. D. (2002). Molecular dynamics and free energy analyses of cathepsin D-inhibitor interactions: insight into structure-based ligand design. J. Med. Chem. 45, 1412–1419. doi: 10.1021/jm010338j
Jerabek, J., Uliassi, E., Guidotti, L., Korabecny, J., Soukup, O., Sepsova, V., et al. (2017). Tacrine-resveratrol fused hybrids as multi-target-directed ligands against Alzheimer's disease. Eur. J. Med. Chem. 127, 250–262. doi: 10.1016/j.ejmech.2016.12.048
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
Kästner, J. (2011). Umbrella sampling. Wiley Interdiscipl. Rev. Comput. Mol. Sci. 1, 932–942. doi: 10.1002/wcms.66
Kollman, P. A., Massova, I., Reyes, C., Kuhn, B., Huo, S., Chong, L., et al. (2000). Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models. Acc. Chem. Res. 33, 889–897. doi: 10.1021/ar000033j
Kuhn, B., and Kollman, P. A. (2000). Binding of a diverse set of ligands to avidin and streptavidin: an accurate quantitative prediction of their relative affinities by a combination of molecular mechanics and continuum solvent models. J. Med. Chem. 43, 3786–3791. doi: 10.1021/jm000241h
Kumar, S., Rosenberg, J. M., Bouzida, D., Swendsen, R. H., and Kollman, P. A. (1992). THE weighted histogram analysis method for free-energy calculations on biomolecules. I. The method. J. Comput. Chem. 13, 1011–1021. doi: 10.1002/jcc.540130812
Lagadic-Gossmann, D., Rissel, M., Le Bot, M. A., and Guillouzo, A. (1998). Toxic effects of tacrine on primary hepatocytes and liver epithelial cells in culture. Cell Biol. Toxicol. 14, 361–373. doi: 10.1023/A:1007589808761
Lin, H., Li, Q., Gu, K., Zhu, J., Jiang, X., Chen, Y., et al. (2017). Therapeutic agents in Alzheimer's disease through a multi-targetdirected ligands strategy: recent progress based on tacrine core. Curr. Top. Med. Chem. 17, 3000–3016. doi: 10.2174/1568026617666170717114944
Lüdemann, S. K., Lounnas, V., and Wade, R. C. (2000). How do substrates enter and products exit the buried active site of cytochrome P450cam? 1. Random expulsion molecular dynamics investigation of ligand access channels and mechanisms. J. Mol. Biol. 303, 797–811. doi: 10.1006/jmbi.2000.4154
Magistrato, A., Sgrignani, J., Krause, R., and Cavalli, A. (2017). Single or multiple access channels to the CYP450s active site? An answer from free energy simulations of the human aromatase enzyme. J. Phys. Chem. Lett. 8, 2036–2042. doi: 10.1021/acs.jpclett.7b00697
Maier, J. A., Martinez, C., Kasavajhala, K., Wickstrom, L., Hauser, K. E., and Simmerling, C. (2015). ff14SB: improving the accuracy of protein side chain and backbone parameters from ff99SB. J. Chem. Theory Comput. 11, 3696–3713. doi: 10.1021/acs.jctc.5b00255
Masson, P., Xie, W., Froment, M.-T., Levitsky, V., Fortier, P.-L., Albaret, C., et al. (1999). Interaction between the peripheral site residues of human butyrylcholinesterase, D70 and Y332, in binding and hydrolysis of substrates. Biochim. Biophys. Acta Protein Struct. Mol. Enzymol. 1433, 281–293. doi: 10.1016/S0167-4838(99)00115-6
Mesulam, M., Guillozet, A., Shaw, P., and Quinn, B. (2002). Widely spread butyrylcholinesterase can hydrolyze acetylcholine in the normal and Alzheimer brain. Neurobiol. Dis. 9, 88–93. doi: 10.1006/nbdi.2001.0462
Mesulam, M. M., Guillozet, A., Shaw, P., Levey, A., Duysen, E. G., and Lockridge, O. (2002). Acetylcholinesterase knockouts establish central cholinergic pathways and can use butyrylcholinesterase to hydrolyze acetylcholine. Neuroscience 110, 627–639. doi: 10.1016/S0306-4522(01)00613-3
Munnaluri, R., Sivan, S. K., and Manga, V. (2015). Molecular docking and MM/GBSA integrated protocol for designing small molecule inhibitors against HIV-1 gp41. Med. Chem. Res. 24, 829–841. doi: 10.1007/s00044-014-1185-8
Mushtaq, G., Greig, N. H., Khan, J. A., and Kamal, M. A. (2014). Status of acetylcholinesterase and butyrylcholinesterase in Alzheimer's disease and type 2 diabetes mellitus. CNS Neuro. Disord. Drug Targets 13, 1432–1439. doi: 10.2174/1871527313666141023141545
Nachon, F., Carletti, E., Ronco, C., Trovaslet, M., Nicolet, Y., Jean, L., et al. (2013). Crystal structures of human cholinesterases in complex with huprine W and tacrine: elements of specificity for anti-Alzheimer's drugs targeting acetyl- and butyryl-cholinesterase. Biochem. J. 453, 393–399. doi: 10.1042/BJ20130013
Nordberg, A., Ballard, C., Bullock, R., Darreh-Shori, T., and Somogyi, M. (2013). A review of butyrylcholinesterase as a therapeutic target in the treatment of Alzheimer's disease. Prim. Care Companion CNS Disord. 15:PCC.12r01412. doi: 10.4088/PCC.12r01412
Okten, S., Ekiz, M., Kocyigit, U. M., Tutar, A., Celik, I., Akkurt, M., et al. (2019). Synthesis, characterization, crystal structures, theoretical calculations and biological evaluations of novel substituted tacrine derivatives as cholinesterase and carbonic anhydrase enzymes inhibitors. J. Mol. Struct. 1175, 906–915. doi: 10.1016/j.molstruc.2018.08.063
Pastor, R. W., Brooks, B. R., and Szabo, A. (1988). An analysis of the accuracy of langevin and molecular dynamics algorithms. Mol. Phys. 65, 1409–1419. doi: 10.1080/00268978800101881
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
Qizilbash, N., Birks, J., López Arrieta, J., Lewington, S., and Szeto, S. (1999). Tacrine for Alzheimer's disease. Cochrane Database Syst. Rev. 2000:CD000202. doi: 10.1002/14651858.CD000202
Rosini, M., Simoni, E., Milelli, A., Minarini, A., and Melchiorre, C. (2014). Oxidative stress in Alzheimer's disease: are we connecting the dots? J. Med. Chem. 57, 2821–2831. doi: 10.1021/jm400970m
Rostkowski, M., Olsson, M. H., Søndergaard, C. R., and Jensen, J. H. (2011). Graphical analysis of pH-dependent properties of proteins predicted using PROPKA. BMC Struct. Biol. 11:6. doi: 10.1186/1472-6807-11-6
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
Sameem, B., Saeedi, M., Mahdavi, M., and Shafiee, A. (2017). A review on tacrine-based scaffolds as multi-target drugs (MTDLs) for Alzheimer's disease. Eur. J. Med. Chem. 128, 332–345. doi: 10.1016/j.ejmech.2016.10.060
Santos, M. A., Chand, K., and Chaves, S. (2016). Recent progress in repositioning Alzheimer's disease drugs based on a multitarget strategy. Future Med. Chem. 8, 2113–2142. doi: 10.4155/fmc-2016-0103
Shao, Q., and Zhu, W. (2018). Assessing AMBER force fields for protein folding in an implicit solvent. Phys. Chem. Chem. Phys. 20, 7206–7216. doi: 10.1039/C7CP08010G
She, N., Zhao, Y., Hao, J., Xie, S., and Wang, C. (2019). Uridine diphosphate release mechanism in O-N-acetylglucosamine (O-GlcNAc) transferase catalysis. Biochim. Biophys. Acta Gen. Subj. 1863, 609–622. doi: 10.1016/j.bbagen.2018.12.005
Srivastava, H. K., and Sastry, G. N. (2012). Molecular dynamics investigation on a series of HIV protease inhibitors: assessing the performance of MM-PBSA and MM-GBSA approaches. J. Chem. Inf. Model. 52, 3088–3098. doi: 10.1021/ci300385h
Summers, W. K. (2006). Tacrine, and Alzheimer's treatments. J. Alzheimers Dis. 9, 439–445. doi: 10.3233/JAD-2006-9S350
Sun, X., Chen, W.-D., and Wang, Y.-D. (2015). β-Amyloid: the key peptide in the pathogenesis of Alzheimer's disease. Front. Pharmacol. 6:221. doi: 10.3389/fphar.2015.00221
Talesa, V. N. (2001). Acetylcholinesterase in Alzheimer's disease. Mech. Ageing Dev. 122, 1961–1969. doi: 10.1016/S0047-6374(01)00309-8
Tschanz, J. T., and Andersen, A. (2011). “Alzheimer's Dementia,” in Encyclopedia of Clinical Neuropsychology, eds J. S. Kreutzer, J. Deluca, B. Caplan (New York, NY: Springer New York), 99–105. doi: 10.1007/978-0-387-79948-3_489
Vashisth, H., and Abrams, C. F. (2008). Ligand escape pathways and (Un)binding free energy calculations for the hexameric insulin-phenol complex. Biophys. J. 95, 4193–4204. doi: 10.1529/biophysj.108.139675
Wang, J., Hou, T., and Xu, X. (2006). Recent advances in free energy calculations with a combination of molecular mechanics and continuum models. Curr. Comput. Aided Drug Des. 2, 287–306. doi: 10.2174/157340906778226454
Wang, J. M., 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
Wang, W., Donini, O., Reyes, C. M., and Kollman, P. A. (2001). Biomolecular simulations: recent developments in force fields, simulations of enzyme catalysis, protein-ligand, protein-protein, and protein-nucleic acid noncovalent interactions. Annu. Rev. Biophys. Biomol. Struct. 30, 211–243. doi: 10.1146/annurev.biophys.30.1.211
Wang, W., and Kollman, P. A. (2000). Free energy calculations on dimer stability of the HIV protease using molecular dynamics and a continuum solvent model. J. Mol. Biol. 303, 567–582. doi: 10.1006/jmbi.2000.4057
Wei, Y., Zhang, Z., She, N., Chen, X., Zhao, Y., and Zhang, J. (2019). Atomistic insight into the inhibition mechanisms of suppressors of cytokine signaling on janus kinase. Phys. Chem. Chem. Phys. 21, 12905–12915. doi: 10.1039/C9CP02257K
Xie, W., Stribley, J. A., Chatonnet, A., Wilder, P. J., Rizzino, A., McComb, R. D., et al. (2000). Postnatal developmental delay and supersensitivity to organophosphate in gene-targeted mice lacking acetylcholinesterase. J. Pharmacol. Exp. Ther. 293, 896–902.
Zhang, X., Wang, J., Hong, C., Luo, W., and Wang, C. (2015). Design, synthesis and evaluation of genistein-polyamine conjugates as multi-functional anti-Alzheimer agents. Acta Pharm. Sin. B 5, 67–73. doi: 10.1016/j.apsb.2014.12.008
Zhao, Y., Chen, N., Wang, C., and Cao, Z. (2016). A comprehensive understanding of enzymatic catalysis by hydroxynitrile lyases with S stereoselectivity from the α/β-hydrolase superfamily: revised role of the active-site lysine and kinetic behavior of substrate delivery and sequential product release. ACS Catal. 6, 2145–2157. doi: 10.1021/acscatal.5b02855
Zhao, Y., She, N., Ma, Y., Wang, C., and Cao, Z. (2018). A description of enzymatic catalysis in N-acetylhexosamine 1-kinase: concerted mechanism of two-magnesium-ion-assisted GlcNAc phosphorylation, flexibility behavior of lid motif upon substrate recognition, and water-assisted GlcNAc-1-P release. ACS Catal. 8, 4143–4159. doi: 10.1021/acscatal.8b00006
Keywords: Alzheimer's disease, tacrine, butyrylcholinesterase, release mechanism molecular dynamics, MM/GBSA, umbrella sampling
Citation: Zhang Z, Fan F, Luo W, Zhao Y and Wang C (2020) Molecular Dynamics Revealing a Detour-Forward Release Mechanism of Tacrine: Implication for the Specific Binding Characteristics in Butyrylcholinesterase. Front. Chem. 8:730. doi: 10.3389/fchem.2020.00730
Received: 05 June 2019; Accepted: 14 July 2020;
Published: 25 August 2020.
Edited by:
Nino Russo, University of Calabria, ItalyReviewed by:
Hujun Xie, Zhejiang Gongshang University, ChinaVicent Moliner, University of Jaume I, Spain
Copyright © 2020 Zhang, Fan, Luo, Zhao and Wang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Yuan Zhao, zhaoyuan@henu.edu.cn; Chaojie Wang, wcjsxq@henu.edu.cn
†These authors have contributed equally to this work