Molecular Modeling on Inhibitor Complexes and Active-Site Dynamics of Cytochrome P450 C17, a Target for Prostate Cancer Therapy

A molecular model for the P450 enzyme cytochrome P450 C17 (CYP17) is presented based on sequence alignments of multiple template structures and homology modeling. This enzyme plays a central role in the biosynthesis of testosterone and is emerging as a major target in prostate cancer, with the recently developed inhibitor abiraterone currently in advanced clinical trials. The model is described in detail, together with its validation, by providing structural explanations to available site-directed mutagenesis data. The CYP17 molecule in this model is in the form of a triangular prism, with an edge of ∼ 55 Å and a thickness of ∼ 37 Å. It is predominantly helical, comprising 13 α helices interspersed by six 3₁₀ helices and 11 β-sheets. Multinanosecond molecular dynamics simulations in explicit solvent have been carried out, and principal components analysis has been used to reveal the details of dynamics around the active site. Coarse-grained methods have also been used to verify low-frequency motions, which have been correlated with active-site gating. The work also describes the results of docking synthetic inhibitors, including the drug abiraterone and the natural substrate pregnenolone, in the CYP17 active site together with molecular dynamics simulations on the complexes.     

Strength of hydrogen bond network takes crucial roles in the dissociation process of inhibitors from the HIV-1 protease binding pocket

To understand the underlying mechanisms of significant differences in dissociation rate constant among different inhibitors for HIV-1 protease, we performed steered molecular dynamics (SMD) simulations to analyze the entire dissociation processes of inhibitors from the binding pocket of protease at atomistic details. We found that the strength of hydrogen bond network between inhibitor and the protease takes crucial roles in the dissociation process. We showed that the hydrogen bond network in the cyclic urea inhibitors AHA001/XK263 is less stable than that of the approved inhibitor ABT538 because of their large differences in the structures of the networks. In the cyclic urea inhibitor bound complex, the hydrogen bonds often distribute at the flap tips and the active site. In contrast, there are additional accessorial hydrogen bonds formed at the lateral sides of the flaps and the active site in the ABT538 bound complex, which take crucial roles in stabilizing the hydrogen bond network. In addition, the water molecule W301 also plays important roles in stabilizing the hydrogen bond network through its flexible movement by acting as a collision buffer and helping the rebinding of hydrogen bonds at the flap tips. Because of its high stability, the hydrogen bond network of ABT538 complex can work together with the hydrophobic clusters to resist the dissociation, resulting in much lower dissociation rate constant than those of cyclic urea inhibitor complexes. This study may provide useful guidelines for design of novel potent inhibitors with optimized interactions.     

CYP17 inhibitors for prostate cancer therapy

Prostate cancer (PC) is now the second most prevalent cause of death in men in the USA and Europe. At present, the major treatment options include surgical or medical castration. These strategies cause ablation of the production of testosterone (T), dihydrotestosterone (DHT) and related androgens by the testes. However, because these procedures do not affect adrenal, prostate and other tissues' androgen production, they are often combined with androgen receptor antagonists to block their action. Indeed, recent studies have unequivocally established that in castration-resistant prostate cancer (CRPC) many androgen-regulated genes become re-expressed and tissue androgen levels increase despite low serum levels. Clearly, inhibition of the key enzyme which catalyzes the biosynthesis of androgens from pregnane precursors, 17α-hydroxy/17,20-lyase (hereafter referred to as CYP17) could prevent androgen production from all sources. Thus, total ablation of androgen production by potent CYP17 inhibitors may provide effective treatment of prostate cancer patients. This review highlights the role of androgen biosynthesis in the progression of prostate cancer and the impact of CYP17 inhibitors, such as ketoconazole, abiraterone acetate, VN/124-1 (TOK-001) and TAK-700 in the clinic and in clinical development. Article from the special issue on Targeted Inhibitors.     

Molecular insight into pseudolysin inhibition using the MM-PBSA and LIE methods

Pseudolysin, the extracellullar elastase of Pseudomonas aeruginosa (EC: 3.4.24.26) plays an important role in the pathogenesis of P. aeruginosa infections. In the present study, molecular dynamics simulations and theoretical affinity predictions were used to gain molecular insight into pseudolysin inhibition. Four low molecular weight inhibitors were docked at their putative binding sites and molecular dynamics (MD) simulations were performed for 5.0 ns, and the free energy of binding was calculated by the linear interaction energy method. The number and the contact surface area of stabilizing hydrophobic, aromatic, and hydrogen bonding interactions appears to reflect the affinity differences between the inhibitors. The proteinaceous inhibitor, Streptomyces metalloproteinase inhibitor (SMPI) was docked in three different binding positions and MD simulations were performed for 3.0 ns. The MD trajectories were used for molecular mechanics-Poisson-Boltzmann surface area analysis of the three binding positions. Computational alanine scanning of the average pseudolysin-SMPI complexes after MD revealed residues at the pseudolysin-SMPI interface giving the main contribution to the free energy of binding. The calculations indicated that SMPI interacts with pseudolysin via the rigid active site loop, but that also contact sites outside this loop contribute significantly to the free energy of association.     

MD simulations reveal alternate conformations of the oxyanion hole in the Zika virus NS2B/NS3 protease

Recent crystallography studies have shown that the binding site oxyanion hole plays an important role in inhibitor binding, but can exist in two conformations (active/inactive). We have undertaken molecular dynamics (MD) calculations to better understand oxyanion hole dynamics and thermodynamics. We find that the Zika virus (ZIKV) NS2B/NS3 protease maintains a stable closed conformation over multiple 100-ns conventional MD simulations in both the presence and absence of inhibitors. The S1, S2, and S3 pockets are stable as well. However, in two of eight simulations, the A132-G133 peptide bond in the binding pocket of S1' spontaneously flips to form a 3₁₀-helix that corresponds to the inactive conformation of the oxyanion hole, and then maintains this conformation until the end of the 100-ns conventional MD simulations without inversion of the flip. This conformational change affects the S1' pocket in ZIKV NS2B/NS3 protease active site, which is important for small molecule binding. The simulation results provide evidence at the atomic level that the inactive conformation of the oxyanion hole is more favored energetically when no specific interactions are formed between substrate/inhibitor and oxyanion hole residues. Interestingly, however, transition between the active and inactive conformation of the oxyanion hole can be observed by boosting the valley potential in accelerated MD simulations. This supports a proposed induced-fit mechanism of ZIKV NS2B/NS3 protease from computational methods and provides useful direction to enhance inhibitor binding predictions in structure-based drug design.     

Insights into the function of the zinc hydroxide-Thr199-Glu106 hydrogen bonding network in carbonic anhydrases

The exact functional role of the zinc hydroxide (water)-Thr199-Glu106 hydrogen bond network in the carbonic anhydrases is unknown. However, from the results of molecular dynamics simulations (MD) we are able to better define its function. From computer graphics analysis and MD simulations on the zinc hydroxide form of human carbonic anhydrase II we find that this interaction forces the hydroxide hydrogen atom to be in a "down" position relative to the deep water-binding pocket. From previous work we have found that this pocket is a high-affinity binding site for CO2. We also note that during the timescale of our simulation (126 ps) the hydrogen bonds between the hydroxide hydrogen atom and Thr199 and the one between Thr199 and Glu106 are not fluxional. We propose that the role of the zinc hydroxide (water)-Thr199-Glu106 hydrogen bond network is to lock the hydrogen atom in the down position in order to expose the CO2 molecule bound in the deep water pocket to a lone pair of the hydroxide oxygen atom. This would allow for the rapid reaction of the CO2 molecule around the zinc ion. Furthermore, if the hydroxide hydrogen atom were not locked in the down position the binding of CO2 to the deep water pocket could be interfered with by the unrestrained hydroxide hydrogen atom (e.g. the N-Zn-O-H torsion could undergo rotational transitions that would partially block the deep water pocket). In summary, the roles we ascribe to this hydrogen bonding network are (1) to allow for facile access of CO2 to the deep water pocket and (2) to allow for maximal exposure of a hydroxide oxygen lone pair to the CO2 carbon atom.     

Molecular dynamics simulations of the Apo-, Holo-, and acyl-forms of Escherichia coli acyl carrier protein

Acyl carrier protein (ACP) is an essential co-factor protein in fatty acid biosynthesis that shuttles covalently bound fatty acyl intermediates in its hydrophobic pocket to various enzyme partners. To characterize acyl chain-ACP interactions and their influence on enzyme interactions, we performed 19 molecular dynamics (MD) simulations of Escherichia coli apo-, holo-, and acyl-ACPs. The simulations were started with the acyl chain in either a solvent-exposed or a buried conformation. All four short-chain (< or = C10) and one long-chain (C16) unbiased acyl-ACP MD simulation show the transition of the solvent-exposed acyl chain into the hydrophobic pocket of ACP, revealing its pathway of acyl chain binding. Although the acyl chain resides inside the pocket, Thr-39 and Glu-60 at the entrance stabilize the phosphopantetheine linker through hydrogen bonding. Comparisons of the different ACP forms indicate that the loop region between helices II and III and the prosthetic linker may aid in substrate recognition by enzymes of fatty acid synthase systems. The MD simulations consistently show that the hydrophobic binding pocket of ACP is best suited to accommodate an octanoyl group and is capable of adjusting in size to accommodate chain lengths as long as decanoic acid. The simulations also reveal a second, novel binding mode of the acyl chains inside the hydrophobic binding pocket directed toward helix I. This study provides a detailed dynamic picture of acyl-ACPs that is in excellent agreement with available experimental data and, thereby, provides a new understanding of enzyme-ACP interactions.     

Theoretical studies on the interaction of biphenyl inhibitors with Mycobacterium tuberculosis protein tyrosine phosphatase MptpB

MptpB is an essential secreted virulence factor for M. tuberculosis. Inhibition of MptpB impairs mycobacterial survival in host macrophages and thus helps reduce tuberculosis infections. However, the binding mode of the biphenyl inhibitors, which are known as some of the most potent MptpB inhibitors, remains unclear. In this study, to understand the interactions between biphenyl inhibitors and MptpB, docking and molecular dynamics simulations were carried out using AutoDock and GROMACS softwares. Calculation results show that all the biphenyl inhibitors can be docked to the binding site of MptpB, with the acid warheads forming a hydrogen bond network at the active site. But the binding modes of other terminals of these inhibitors are different. The cyclohexyl and trifluoromethyl substituents at R1 and R2 sites are necessary for the inhibitors to adopt their double-site binding mechanism. The estimated binding affinities are basically consistent with the experimental results. MD simulations show that these binding complexes display different stability.     

Binding free energies for nicotine analogs inhibiting cytochrome P450 2A6 by a combined use of molecular dynamics simulations and QM/MM-PBSA calculations

Molecular dynamics (MD) simulations and hybrid quantum mechanical/molecular mechanical (QM/MM) calculations have been performed to explore the dynamic behaviors of cytochrome P450 2A6 (CYP2A6) binding with nicotine analogs (that are typical inhibitors) and to calculate their binding free energies in combination with Poisson–Boltzmann surface area (PBSA) calculations. The combined MD simulations and QM/MM-PBSA calculations reveal that the most important structural parameters affecting the CYP2A6-inhibitor binding affinity are two crucial internuclear distances, that is, the distance between the heme iron atom of CYP2A6 and the coordinating atom of the inhibitor, and the hydrogen-bonding distance between the N297 side chain of CYP2A6 and the pyridine nitrogen of the inhibitor. The combined MD simulations and QM/MM-PBSA calculations have led to dynamic CYP2A6-inhibitor binding structures that are consistent with the observed dynamic behaviors and structural features of CYP2A6-inhibitor binding, and led to the binding free energies that are in good agreement with the experimentally-derived binding free energies. The agreement between the calculated binding free energies and the experimentally-derived binding free energies suggests that the combined MD and QM/MM-PBSA approach may be used as a valuable tool to accurately predict the CYP2A6-inhibitor binding affinities in future computational design of new, potent and selective CYP2A6 inhibitors.     

Regioselectivity of CYP2B6: homology modeling,molecular dynamics simulation,docking

Human cytochrome P450 (CYP) 2B6 activates the anticancer prodrug cyclophosphamide (CPA) by 4-hydroxylation. In contrast, the same enzyme catalyzes N-deethylation of a structural isomer, the prodrug ifosfamide (IFA), thus causing severe adverse drug effects. To model the molecular interactions leading to a switch in regioselectivity, the structure of CYP2B6 was modeled based on the structure of rabbit CYP2C5. We modeled the missing 22-residue loop in CYP2C5 between helices F and G (the F-G loop), which is not resolved in the X-ray structure, by molecular dynamics (MD) simulations using a simulated annealing protocol. The modeled conformation of the loop was validated by unconstrained MD simulations of the complete enzymes (CYP2C5 and CYP2B6) in water for 70 and 120 ps, respectively. The simulations were stable and led to a backbone r.m.s. deviation of 1.7 A between the two CYPs.The shape of the substrate binding site of CYP2B6 was further analyzed. It consists of three well-defined hydrophobic binding pockets adjacent to the catalytic heme. Size, shape and hydrophobicity of these pockets were compared to the shapes of the two structurally isomeric substrates. In their preferred orientation in the binding site, both substrates fill all three binding pockets without repulsive interactions. The distance to the heme iron is short enough for 4-hydroxylation and N-deethylation to occur for CPA and IFA, respectively. However, if the substrates are docked in the non-preferred orientation (such that 4-hydroxylation and N-deethylation would occur for IFA and CPA, respectively), one pocket is left empty, and clashes were observed between the substrates and the enzyme.     

The flexibility of P-glycoprotein for its poly-specific drug binding from molecular dynamics simulations

The multidrug efflux pump P-glycoprotein (P-gp) contributes to multidrug resistance in about half of human cancers. Recently, high resolution X-ray crystal structures of mouse P-gp (inward-facing) were reported, which significantly facilitates the understanding of the function of P-gp and the structure-based design of inhibitors for P-gp. Here we perform 20?ns molecular dynamics simulations of inward-facing P-gp with/without ligand in explicit lipid and water to investigate the flexibility of P-gp for its poly-specific drug binding. By analyzing the interactions between P-gp and QZ59-RRR or QZ59-SSS, we summarize the important residues and the flexibility of different parts of P-gp. Particularly, the flexibility of the side chains of aromatic residues (Phe and Tyr) allows them to form rotamers with different orientations in the binding pocket, which plays a critical role for the poly-specificity of the drug-binding cavity of P-gp. MD simulations reveal that trans-membrane (TM) TM12 and TM6 are flexible and contribute to the poly-specific drug binding, while TM4 and TM5 are rigid and stabilize the whole structure. We also construct outward-facing P-gp based on the MsbA structure and perform 20?ns MD simulations. The comparison between the MD results for outward-facing P-gp and those for inward-facing P-gp shows that the TM parts in outward-facing P-gp undergo significant conformational change to facilitate the export of small molecules.     

Why human cytochrome P450c21 is a progesterone 21-hydroxylase

Human cytochrome P450c21 (steroid 21-hydroxylase, CYP21A2) catalyzes the 21-hydroxylation of progesterone (P4) and its preferred substrate 17α-hydroxyprogestrone (17OHP4). CYP21A2 activities, which are required for cortisol and aldosterone biosynthesis, involve the formation of energetically disfavored primary carbon radicals. Therefore, we hypothesized that the binding of P4 and 17OHP4 to CYP21A2 restricts access of the reactive heme-oxygen complex to the C-21 hydrogen atoms, suppressing oxygenation at kinetically more favorable sites such as C-17 and C-16, which are both hydroxylated by cytochrome P450c17 (CYP17A1). We reasoned that expansion of the CYP21A2 substrate-binding pocket would increase substrate mobility and might yield additional hydroxylation activities. We built a computer model of CYP21A2 based principally on the crystal structure of CYP2C5, which also 21-hydroxylates P4. Molecular dynamics simulations indicate that binding of the steroid nucleus perpendicular to the plane of the CYP21A2 heme ring limits access of the heme oxygen to the C-21 hydrogen atoms. Residues L107, L109, V470, I471, and V359 were found to contribute to the CYP21A2 substate-binding pocket. Mutation of V470 and I471 to alanine or glycine preserved P4 21-hydroxylase activity, and mutations of L107 or L109 were inactive. Mutations V359A and V359G, in contrast, acquired 16α-hydroxylase activity, accounting for 40% and 90% of the P4 metabolites, respectively. We conclude that P4 binds to CYP21A2 in a fundamentally different orientation than to CYP17A1 and that expansion of the CYP21A2 substrate-binding pocket allows additional substrate trajectories and metabolic switching.     

Exploration of the conformational landscape in pregnane X receptor reveals a new binding pocket

Ligand‐regulated pregnane X receptor (PXR), a member of the nuclear receptor superfamily, plays a central role in xenobiotic metabolism. Despite its critical role in drug metabolism, PXR activation can lead to adverse drug‐drug interactions and early stage metabolism of drugs. Activated PXR can induce cancer drug resistance and enhance the onset of malignancy. Since promiscuity in ligand binding makes it difficult to develop competitive inhibitors targeting PXR ligand binding pocket (LBP), it is essential to identify allosteric sites for effective PXR antagonism. Here, molecular dynamics (MD) simulation studies unravelled the existence of two different conformational states, namely “expanded” and “contracted”, in apo PXR ligand binding domain (LBD). Ligand binding events shifted this conformational equilibrium and locked the LBD in a single “ligand‐adaptable” conformational state. Ensemble‐based computational solvent mapping identified a transiently open potential small molecule binding pocket between α5 and α8 helices, named “α8 pocket”, whose opening‐closing mechanism directly correlated with the conformational shift in LBD. A virtual hit identified through structure‐based virtual screening against α8 pocket locks the pocket in its open conformation. MD simulations further revealed that the presence of small molecule at allosteric site disrupts the LBD dynamics and locks the LBD in a “tightly‐contracted” conformation. The molecular details provided here could guide new structural studies to understand PXR activation and antagonism.     

An In Silico Analysis of the Binding Modes and Binding Affinities of Small Molecule Modulators of PDZ-Peptide Interactions

Inhibitors of PDZ-peptide interactions have important implications in a variety of biological processes including treatment of cancer and Parkinson’s disease. Even though experimental studies have reported characterization of peptidomimetic inhibitors of PDZ-peptide interactions, the binding modes for most of them have not been characterized by structural studies. In this study we have attempted to understand the structural basis of the small molecule-PDZ interactions by in silico analysis of the binding modes and binding affinities of a set of 38 small molecules with known K_i or K_d values for PDZ2 and PDZ3 domains of PSD-95 protein. These two PDZ domains show differential selectivity for these compounds despite having a high degree of sequence similarity and almost identical peptide binding pockets. Optimum binding modes for these ligands for PDZ2 and PDZ3 domains were identified by using a novel combination of semi-flexible docking and explicit solvent molecular dynamics (MD) simulations. Analysis of the binding modes revealed most of the peptidomimectic ligands which had high K_i or K_d moved away from the peptide binding pocket, while ligands with high binding affinities remained in the peptide binding pocket. The differential specificities of the PDZ2 and PDZ3 domains primarily arise from differences in the conformation of the loop connecting βB and βC strands, because this loop interacts with the N-terminal chemical moieties of the ligands. We have also computed the MM/PBSA binding free energy values for these 38 compounds with both the PDZ domains from multiple 5 ns MD trajectories on each complex i.e. a total of 228 MD trajectories of 5 ns length each. Interestingly, computational binding free energies show good agreement with experimental binding free energies with a correlation coefficient of approximately 0.6. Thus our study demonstrates that combined use of docking and MD simulations can help in identification of potent inhibitors of PDZ-peptide complexes.     

Multiple binding modes for palmitate to barley lipid transfer protein facilitated by the presence of proline 12

Molecular dynamics simulations have been used to characterise the binding of the fatty acid ligand palmitate in the barley lipid transfer protein 1 (LTP) internal cavity. Two different palmitate binding modes (1 and 2), with similar protein–ligand interaction energies, have been identified using a variety of simulation strategies. These strategies include applying experimental protein–ligand atom–atom distance restraints during the simulation, or protonating the palmitate ligand, or using the vacuum GROMOS 54B7 force‐field parameter set for the ligand during the initial stages of the simulations. In both the binding modes identified the palmitate carboxylate head group hydrogen bonds with main chain amide groups in helix A, residues 4 to 19, of the protein. In binding mode 1 the hydrogen bonds are to Lys 11, Cys 13, and Leu 14 and in binding mode 2 to Thr 15, Tyr 16, Val 17, Ser 24 and also to the OH of Thr 15. In both cases palmitate binding exploits irregularity of the intrahelical hydrogen‐bonding pattern in helix A of barley LTP due to the presence of Pro 12. Simulations of two variants of barley LTP, namely the single mutant Pro12Val and the double mutant Pro12Val Pro70Val, show that Pro 12 is required for persistent palmitate binding in the LTP cavity. Overall, the work identifies key MD simulation approaches for characterizing the details of protein–ligand interactions in complexes where NMR data provide insufficient restraints.     

Molecular mechanics PBSA ligand binding energy and interaction of Efavirenz derivatives with HIV-1 reverse transcriptase

In order to evaluate the properties of several HIV-1 reverse transcripase(RT) inhibitors, Efavirenz (SUSTIVA) and a set of its derivatives (benzoxazinones) have been placed into the nonnucleoside analogue binding site of the enzyme by molecular docking. The resulting geometries were used for a molecular dynamics simulation and binding energy calculations. The enzyme-inhibitor binding energies were estimated from experimental inhibitory activities (IC90). The correlation of the predicted and experimental binding energies were satisfactory acceptable as indicated by r2 = 0.865. Based on MD simulations, the obtained results indicate that the tight association of the ligand to the HIV-1 RT binding pocket was based on hydrogen bonding between Efavirenz's N1 and the oxygen of the backbone of Lys 101, with an estimated average distance of 1.88 A. Moreover, electrostatic interaction was mainly contributed by two amino acid residues in the binding site; Lys 101 and His 235. MD simulations open the possibility to study the reaction of the flexible enzyme to those substances as well as the overall affinity.     

Synthesis, biological evaluation and molecular modelling studies of methyleneimidazole substituted biaryls as inhibitors of human 17alpha-hydroxylase-17,20-lyase (CYP17). Part I: Heterocyclic modifications of the core structure

Novel chemical entities were prepared via Suzuki and S(N) reaction as AC-ring substrate mimetics of CYP17. The synthesised compounds 1-31 were tested for activity using human CYP17 expressed in Escherichia coli. Promising compounds were tested for selectivity against hepatic CYP enzymes (3A4, 2D6, 1A2, 2C9, 2C19, 2B6). Two potent inhibitors (27, IC50 = 373 nM/28, IC50 = 953 nM) were further examined in rats regarding their effects on plasma testosterone levels and their pharmacokinetic properties. Compound 28 was similarly active as abiraterone and showed better pharmacokinetic properties (higher bioavailability, t(1/2) 9.5 h vs 1.6 h). Docking studies revealed two new binding modes different from the one of the substrates and steroidal inhibitors.     

All-atom molecular dynamics simulations reveal significant differences in interaction between antimycin and conserved amino acid residues in bovine and bacterial bc1 complexes

Antimycin A is the most frequently used specific and powerful inhibitor of the mitochondrial respiratory chain. We used all-atom molecular dynamics (MD) simulations to study the dynamic aspects of the interaction of antimycin A with the Q_i site of the bacterial and bovine bc₁ complexes embedded in a membrane. The MD simulations revealed considerable conformational flexibility of antimycin and significant mobility of antimycin, as a whole, inside the Q_i pocket. We conclude that many of the differences in antimycin binding observed in high-resolution x-ray structures may have a dynamic origin and result from fluctuations of protein and antimycin between multiple conformational states of similar energy separated by low activation barriers, as well as from the mobility of antimycin within the Q_i pocket. The MD simulations also revealed a significant difference in interaction between antimycin and conserved amino acid residues in bovine and bacterial bc₁ complexes. The strong hydrogen bond between antimycin and conserved Asp-228 (bovine numeration) was observed to be frequently broken in the bacterial bc₁ complex and only rarely in the bovine bc₁ complex. In addition, the distances between antimycin and conserved His-201 and Lys-227 were consistently larger in the bacterial bc₁ complex. The observed differences could be responsible for a weaker interaction of antimycin with the bacterial bc₁ complex.