CO migration pathways in cytochrome P450cam studied by molecular dynamics simulations

Previous laser flash photolysis investigations between 100 and 300 K have shown that the kinetics of CO rebinding with cytochrome P450(cam)(camphor) consist of up to four different processes revealing a complex internal dynamics after ligand dissociation. In the present work, molecular dynamics simulations were undertaken on the ternary complex P450(cam)(cam)(CO) to explore the CO migration pathways, monitor the internal cavities of the protein, and localize the CO docking sites. One trajectory of 1 nsec with the protein in a water box and 36 trajectories of 1 nsec in the vacuum were calculated. In each trajectory, the protein contained only one CO ligand on which no constraints were applied. The simulations were performed at 200, 300, and 320 K. The results indicate the presence of seven CO docking sites, mainly hydrophobic, located in the same moiety of the protein. Two of them coincide with xenon binding sites identified by crystallography. The protein matrix exhibits eight persistent internal cavities, four of which corresponding to the ligand docking sites. In addition, it was observed that water molecules entering the protein were mainly attracted into the polar pockets, far away from the CO docking sites. Finally, the identified CO migration pathways provide a consistent interpretation of the experimental rebinding kinetics.     

Spring-loading the active site of cytochrome P450cam

A hydrogen bond network has been identified that adjusts protein-substrate contacts in cytochrome P450(cam) (CYP101A1). Replacing the native substrate camphor with adamantanone or norcamphor causes perturbations in NMR-detected NH correlations assigned to the network, which includes portions of a β sheet and an adjacent helix that is remote from the active site. A mutation in this helix reduces enzyme efficiency and perturbs the extent of substrate-induced spin state changes at the haem iron that accompany substrate binding. In turn, the magnitude of the spin state changes induced by alternate substrate binding parallel the NMR-detected perturbations observed near the haem in the enzyme active site.     

Progress in cytochrome P450 active site modeling

Models capable of predicting the possible involvement of cytochromes P450 in the metabolism of drugs or drug candidates are important tools in drug discovery and development. Ideally, functional information would be obtained from crystal structures of all the cytochromes P450 of interest. Initially, only crystal structures of distantly related bacterial cytochromes P450 were available-comparative modeling techniques were used to bridge the gap and produce structural models of human cytochromes P450, and thereby obtain some useful functional information. A significant step forward in the reliability of these models came four years ago with the first crystal structure of a mammalian cytochrome P450, rabbit CYP2C5, followed by the structures of two human enzymes, CYP2C8 and CYP2C9, and a second rabbit enzyme, CYP2B4. The evolution of a CYP2D6 model, leading to the validation of the model as an in silico tool for predicting binding and metabolism, is presented as a case study.     

Insight into the key active sites of afChiA1 based on molecular dynamics simulations and free energy calculations

In this study, the binding of the enzyme chitinase A1 (afChiA1) from the plant-type Aspergillus fumigatus with four potent inhibitors, allosamidin (ASM), acetazolamide (AZM), 8-chloro-theophylline (CTP) and kinetin (KIT) is investigated by molecular docking, molecular dynamics simulation and binding free energy calculation. The results reveal that the electrostatic interactions play an important role in the stabilisation of the binding of afChiA1 with inhibitors. Based on the binding energy of afChiA1-ligands, the key residues (Gln37 and Trp312) in the active binding pocket of the complex systems are confirmed by molecular mechanics/Poisson–Boltzmann surface area method, and the active inhibitors, ASM and AZM, both could form strong interaction with Gln37 and Trp312, and the non-active ligands, CTP and KIT, could not interact with these two residues, which is consistent with the result of experimental report. Then, it is identified that Gln37 and Trp312 should be one of the important active site residues of afChiA1.     

Molecular dynamics simulations of phenylimidazole inhibitor complexes of cytochrome P450 cam

Molecular dynamics simulations have been performed on three phenylimidazole inhibitor complexes ofP450 _cam, utilizing the X-ray structures and the AMBER suite of programs. Compared to their corresponding optimized X-ray structures, very similar features were observed for the 1-phenylimidazole (1-PI) and 2-phenylimidazole (2-PI) complexes during a 100 ps MD simulation. The 1-PI inhibitor binds as a Type II complex with the imidazole nitrogen as a ligand of the heme iron. Analysis of the inhibitor-enzyme interctions during the MD simulations reveals that electrostatic interactions of the imidazole with the heme and van der Waals interactions of the phenyl ring with nearby hydrophobic residues are dominant. By contrast, 2-PI binds as a Type I inhibitor in the substrate binding pocket, but not as a ligand of the iron. The interactions of this inhibitor are qualitatively different from that of the Type II 1-PI, being mainly electrostatic/H-bonding interactions with a bound water and polar residues. Although the third compound, 4-PI, in common with 1-PI, also binds as a Type II inhibitor, with one nitrogen of the imidazole as a ligand to the iron, the MD average binding orientation deviates significantly from the X-ray structure. The most important changes observed include: (1) the rotation of the imidazole ring of this inhibitor by about 90° to enhance electrostatic interactions of the imidazole NH group with the carbonyl group of LEU244, and (2) the rotation of the carbonyl group of ASP251 to form a H-bond with VAL254. An analysis of the H-bonding network surrounding this substrate in the optimized crystal structure revealed that there is no H-bonding partner either for the free polar NH group in the imidazole ring of 4-phenylimidazole or for the polar carbonyl group of the nearby ASP251 residue. The deviation of the dynamically averaged inhibitor-enzyme structure of the 4-PI complex from the optimized crystal structure can therefore be rationalized as a consequence of the optimization of the electrostatic interactions among the polar groups.     

Molecular dynamics characterization of the active cavity of carboxypeptidase A and some of its inhibitor adducts.

Molecular dynamics (MD) calculations have been performed on carboxypeptidase A and on its adducts with inhibitors, such as d-phenylalanine (dPhe) and acetate. The catalytically essential zinc ion present in the protein was explicitly included in all the simulations. The simulation was carried out over a sphere of 15 A centered on the zinc ion. The crystallographic water molecules were explicitly taken into account; then the protein was solvated with a 18 A sphere of water molecules. MD calculations were carried out for 45-60 ps. There is no large deviation from the available X-ray structures of native and the dPhe adduct for the MD structures. Average MD structures were calculated starting from the X-ray structure of the dPhe adduct, and, from a structure obtained by docking the inhibitor in the native structure. Comparison between these two structures and with that of the native protein shows that some of the key variations produced by inhibitor binding are reproduced by MD calculations. Addition of acetate induces structural changes relevant for the understanding of the interaction network in the active cavity. The structural variations induced by different inhibitors are examined. The effects of these interactions on the catalytic mechanism and on the binding of substrate are discussed.     

Position-resolved free energy of solvation for amino acids in lipid membranes from molecular dynamics simulations

Studies of insertion and interactions of amino acids in lipid membranes are pivotal to our understanding of membrane protein structure and function. Calculating the insertion cost as a function of transmembrane helix sequence is thus an important step towards improved membrane protein prediction and eventually drug design. Here, we present position-dependent free energies of solvation for all amino acid analogs along the membrane normal. The profiles cover the entire region from bulk water to hydrophobic core, and were produced from all-atom molecular dynamics simulations. Experimental differences corresponding to mutations and costs for entire segments match experimental data well, and in addition the profiles provide the spatial resolution currently not available from experiments. Polar side-chains largely maintain their hydration and assume quite ordered conformations, which indicates the solvation cost is mainly entropic. The cost of solvating charged side-chains is not only significantly lower than for implicit solvation models, but also close to experiments, meaning these could well maintain their protonation states inside the membrane. The single notable exception to the experimental agreement is proline, which is quite expensive to introduce in vivo despite its hydrophobicity--a difference possibly explained by kinks making it harder to insert helices in the translocon.     

Controlling the regiospecificity and coupling of cytochrome P450cam: T185F mutant increases coupling and abolishes 3-hydroxynorcamphor product.

Cytochrome P450cam (P450CIA1) catalyzes the hydroxylation of camphor and several substrate analogues such as norcamphor and 1-methyl-norcamphor. Hydroxylation was found experimentally at the 3, 5, and 6 positions of norcamphor, but only at the 5 and 6 positions of 1-methyl-norcamphor. In the catalytic cycle, the hydroxylation of substrate is coupled to the consumption of NADH. For camphor, the degree of coupling is 100%, but for both norcamphor and 1-methyl-norcamphor, the efficiency is dramatically lowered to 12% and 50%, respectively. Based on an examination of the active site of P450cam, it appeared that mutating position 185 might dramatically alter the product specificity and coupling of hydroxylation of norcamphor by P450cam. Analysis of molecular dynamics trajectories of norcamphor bound to the T185F mutant of cytochrome P450cam predicted that hydroxylation at the 3 position should be abolished and that the coupling should be dramatically increased. This mutant was constructed and the product profile and coupling experimentally determined. The coupling was doubled, and hydroxylation at the 3 position was essentially abolished. Both of these results are in agreement with the prediction.     

Computational insights into the binding pattern of mitochondrial calcium uniporter inhibitor through homology modeling,molecular dynamics simulation,binding free energy prediction and density functional theory calculation

Abstract

The mitochondrial calcium uniporter (MCU) is the critical protein of the inner mitochondrial membrane that is the primary mediator for calcium uptake into the mitochondrial matrix. Herein we built the optimal homology model of human MCU which was refined through all-atom molecular dynamics simulation. Then, the binding mode of known inhibitor was predicted through molecular docking method, along with molecular dynamics simulation and binding free energy calculation to verify the docking result and stability of the protein-inhibitor complex. Finally, density functional theory (DFT) calculation enhanced our understanding of the molecular interaction of MCU inhibitor. Our research would provide a deeper insight into the interactions between human MCU and its inhibitor, which boosts to develop novel therapy against MCU related disease.

Communicated by Ramaswamy H. Sarma     

Structural evidence for a functionally relevant second camphor binding site in P450cam: model for substrate entry into a P450 active site

P450cam has long served as a prototype for the cytochrome P450 (CYP) gene family. But, little is known about how substrate enters its active site pocket, and how access is achieved in a way that minimizes exposure of the reactive heme. We hypothesize that P450cam may first bind substrate transiently near the mobile F-G helix that covers the active site pocket. Such a two-step binding process is kinetically required if P450cam rarely populates an open conformation-as suggested by previous literature and the inability to obtain a crystal structure of P450cam in an open conformation. Such a mechanism would minimize exposure of the heme by allowing P450cam to stay in a closed conformation as long as possible, since only brief flexing into an open conformation would be required to allow substrate entry. To test this model, we have attempted to dock a second camphor molecule into the crystal structure of camphor-bound P450cam. The docking identified only one potential entry site pocket, a well-defined cavity on the F-helix side of the F-G flap, 16 A from the heme iron. Location of this entry site pocket is consistent with our NMR T1 relaxation-based measurements of distances for a camphor that binds in fast exchange (active site camphor is known to bind in slow exchange). Presence of a second camphor binding site is also confirmed with [(1)H-(13)C] HSQC titrations of (13)CH3-threonine labeled P450cam. To confirm that camphor can bind outside of the active site pocket, (13)CH3-S-pyridine was bound to the heme iron to physically block the active site, and to serve as an NMR chemical shift probe. Titration of this P450cam-pyridine complex confirms that camphor can bind to a site outside the active site pocket, with an estimated Kd of 43 microM. The two-site binding model that is proposed based on these data is analogous to that recently proposed for CYP3A4, and is consistent with recent crystal structures of P450cam bound to tethered-substrates, which force a partially opened conformation.     

Computational analysis of binding of P1 variants to trypsin

The binding of P1 variants of bovine pancreatic trypsin inhibitor (BPTI) to trypsin has been investigated by means of molecular dynamics simulations. The specific interaction formed between the amino acid at the primary binding (P1) position of the binding loop of BPTI and the specificity pocket of trypsin was estimated by use of the linear interaction energy (LIE) method. Calculations for 13 of the naturally occurring amino acids at the P1 position were carried out, and the results obtained were found to correlate well with the experimental binding free energies. The LIE calculations rank the majority of the 13 variants correctly according to the experimental association energies and the mean error between calculated and experimental binding free energies is only 0.38 kcal/mole, excluding the Glu and Asp variants, which are associated with some uncertainties regarding protonation and the possible presence of counter-ions. The three-dimensional structures of the complex with three of the P1 variants (Asn, Tyr, and Ser) included in this study have not at present been solved by any experimental techniques and, therefore, were modeled on the basis of experimental data from P1 variants of similar size. Average structures were calculated from the MD simulations, from which specific interactions explaining the broad variation in association energies were identified. The present study also shows that explicit treatment of the complex water-mediated hydrogen bonding network at the protein-protein interface is of crucial importance for obtaining reliable binding free energies. The successful reproduction of relative binding energies shows that this type of methodology can be very useful as an aid in rational design and redesign of biologically active macromolecules.     

Investigation on the binding mechanism of loratinib with the c-ros oncogene 1 (ROS1) receptor tyrosine kinase via molecular dynamics simulation and binding free energy calculations

The c-ros oncogene 1 (ROS1) has proven to be an important cancer target for the treatment of various human cancers. The anaplastic lymphoma kinase inhibitor crizotinib has been granted approval for the treatment of patients with ROS1 positive metastatic non-small-cell lung cancer by the Food and Drug Administration on 2016. However, serious resistance due to the secondary mutation of glycine 2032 to arginine (G2032R) was developed in clinical studies. Loratinib (PF-06463922), a macrocyclic analog of crizotinib, showed significantly improved inhibitory activity against wild–type (WT) ROS1 and ROS1^G2032R mutant. To provide insights into the inhibition mechanism, molecular dynamics simulations and free energy calculations were carried out for the complexes of loratinib with WT and G2032R mutated ROS1. The apo-ROS1^WT and apo-ROS1^G2032R systems showed similar RMSF distributions, while ROS1^G2032R-loratinib showed significantly higher than that of WT ROS1-loratinib, which revealed that the binding of loratinib to ROS1^G2032R significantly interfered the ?uctuation of protein. Calculations of binding free energies indicate that G2032R mutation significantly reduces the binding affinity of loratinib for ROS1, which arose mostly from the increase of conformation entropy and the decrease of solvation energy. Furthermore, detailed per-residue binding free energies highlighted the increased and decreased contributions of some residues in the G2032R mutated systems. The present study revealed the detailed inhibitory mechanism of loratinib as potent WT and G2032R mutated ROS1 inhibitor, which was expected to provide a basis for rational drug design.     

Mechanism of the plant cytochrome P450 for herbicide resistance: a modelling study

Plant cytochrome P450 is a key enzyme responsible for the herbicide resistance but the molecular basis of the mechanism is unclear. To understand this, four typical plant P450s and a widely resistant herbicide chlortoluron were analysed by carrying out homology modelling, molecular docking, molecular dynamics simulations and binding free energy analysis. Our results demonstrate that: (i) the putative hydrophobic residues located in the F-helix and polar residues in I-helix are critical in the herbicide resistance; (ii) the binding mode analysis and binding free energy calculation indicate that the distance between catalytic site of chlortoluron and heme of P450, as well as the binding affinity are key elements affecting the resistance for plants. In conclusion, this work provides a new insight into the interactions of plant P450s with herbicide from a molecular level, offering valuable information for the future design of novel effective herbicides which also escape from the P450 metabolism.     

Segment swapping aided the evolution of enzyme function: The case of uroporphyrinogen III synthase

In an earlier study, we showed that two‐domain segment‐swapped proteins can evolve by domain swapping and fusion, resulting in a protein with two linkers connecting its domains. We proposed that a potential evolutionary advantage of this topology may be the restriction of interdomain motions, which may facilitate domain closure by a hinge‐like movement, crucial for the function of many enzymes. Here, we test this hypothesis computationally on uroporphyrinogen III synthase, a two‐domain segment‐swapped enzyme essential in porphyrin metabolism. To compare the interdomain flexibility between the wild‐type, segment‐swapped enzyme (having two interdomain linkers) and circular permutants of the same enzyme having only one interdomain linker, we performed geometric and molecular dynamics simulations for these species in their ligand‐free and ligand‐bound forms. We find that in the ligand‐free form, interdomain motions in the wild‐type enzyme are significantly more restricted than they would be with only one interdomain linker, while the flexibility difference is negligible in the ligand‐bound form. We also estimated the entropy costs of ligand binding associated with the interdomain motions, and find that the change in domain connectivity due to segment swapping results in a reduction of this entropy cost, corresponding to ～20% of the total ligand binding free energy. In addition, the restriction of interdomain motions may also help the functional domain‐closure motion required for catalysis. This suggests that the evolution of the segment‐swapped topology facilitated the evolution of enzyme function for this protein by influencing its dynamic properties. Proteins 2016; 85:46–53. © 2016 Wiley Periodicals, Inc.     

Comparative study of the binding mode between cytochrome P450 17A1 and prostate cancer drugs in the absence of haem iron

Abstract

According to the X-ray crystal structures of CYP17A1 (including its complexes with inhibitors), it is shown that a hydrogen bond exists between CYP17A1 and its inhibitors (such as abiraterone and TOK-001). Previous short MD simulations (50?ns) suggested that the binding of abiraterone to CYP17A1 is stronger than that of TOK-001. In this work, by carrying out long atomistic MD simulations (200?ns) of CYP17A1 and its complexes with abiraterone and TOK-001, we observed a binding mode between CYP17A1 and abiraterone, which is different from the binding mode between CYP17A1 and TOK-001. In the case of abiraterone binding, the unfilled volume in the active site cavity increases the freedom of movement of abiraterone within CYP17A1, leading to the collective motions of the helices G and B′ as well as the breaking of hydrogen bond existing between the 3β-OH group of abiraterone and N202 of CYP17A1. However, the unfilled volume in the active site cavity can be occupied by the benzimidazole ring of TOK-001, restraining the motion of TOK-001. By pulling the two inhibitors (abiraterone and TOK-001) out of the binding pocket in CYP17A1, we discovered that abiraterone and TOK-001 were moved from their binding sites to the surface of protein similarly through the channels formed by the helices G and B′. In addition, based on the free energy calculations, one can see that it is energetically favorable for the two inhibitors (abiraterone and TOK-001) to enter into the binding pocket in CYP17A1.     

Evaluating the suitability of RNA intervention mechanism exerted by some flavonoid molecules against dengue virus MTase RNA capping site: a molecular docking,molecular dynamics simulation,and binding free energy study

Abstract

Dengue virus (DENV) is one of the most dangerous mosquito-borne human pathogens known to the mankind. Currently, no vaccines or standard therapy is avaliable to treate DENV infection. This makes the drug development against DENV more significant and challenging. The MTase domain of DENV RNA RdRp NS5 is a promising drug target, because this domain hosts the RNA capping process of DENV RNA to escape from human immune system. In the present study, we have analysed the RNA intervention mechanism exerted by flavoniod molecules against NS5 MTase RNA capping site by using molecular docking, molecular dynamics simulation and the binding free energy calculations. The results from the docking analysis confirmed that the RNA intervention mecanism is exerted by the quercetagetin (QGN) molecule with all necessary intermolecular interactions and high binding affinity. Notably, QGN forms strong hydrogen bonding interactions with Asn18, Leu20 and Ser150 residues and π???π stacking interaction with Phe25 residue. The apo and QGN bound NS5 MTase and QGN-NS5 MTase complex were used for MD simulation. The results of MD simulation reveal that the RMSD and RMSF values of QGN-MTase complex have increased on comparing the apo protein due to the effect of ligand binding. The binding free energy calulation includes prediction of total binding free energy of ligand-protein complex and per-residue free energy decomposition. The QGN binding to NS5 MTase affects it’s native motion, this result is found from Principal component analysis.

Communicated by Ramaswamy H. Sarma     

Insight into the binding modes and inhibition mechanisms of adamantyl-based 1,3-disubstituted urea inhibitors in the active site of the human soluble epoxide hydrolase

Soluble epoxide hydrolase (sEH) is a promising new target for treating hypertension and inflammation. Considerable efforts have been devoted to develop novel inhibitors. In this study, the binding modes and interaction mechanisms of a series of adamantyl-based 1,3-disubstituted urea inhibitors were investigated by molecular docking, molecular dynamics simulations, binding free energy calculations, and binding energy decomposition analysis. Based on binding affinity, the most favorable binding mode was determined for each inhibitor. The calculation results indicate that the total binding free energy (ΔG_TOT, the sum of enthalpy ΔG_MM-GB/SA, and entropy ?TΔS) presents a good correlation with the experimental inhibitory activity (IC₅₀, r²?=?.99). The van der Waals energy contributes most to the total binding free energy (ΔG_TOT). A detailed discussion on the interactions between inhibitors and those residues located in the active pocket is made based on hydrogen bond and binding modes analysis. According to binding energy decomposition, the residues Asp333 and Trp334 contribute the most to binding free energy in all systems. Furthermore, Hip523 plays a major role in determining this class of inhibitor-binding orientations. Combined with the results of hydrogen bond analysis and binding free energy, we believe that the conserved hydrogen bonds play a role only in anchoring the inhibitors to the exact site for binding and the number of hydrogen bonds may not directly relate to the binding free energy. The results we obtained will provide valuable information for the design of high potency sEH inhibitors.     

Molecular modeling and identification of substrate binding site of orphan human cytochrome P450 4F22

Cytochrome P450s are superfamily of heme proteins which generally monooxygenate hydrophobic compounds. The human cytochrome P450 4F22 (CYP4F22) was categorized into "orphan" CYPs because of its unknown function. CYP4F22 is a potential drug target for cancer therapy. However, three-dimensional structure, the active site topology and substrate specificity of CYP4F22 remain unclear. In this study, a three-dimensional model of human P450 4F22 was constructed by comparative modeling using Modeller 9v5. The resulting model was refined by energy minimization subjected to the quality assessment from both geometric and energetic aspects and was found to be of reasonable quality. Docking approach was employed to dock arachidonic acid into the active site of CYP4F22 in order to probe the ligand-binding modes. As a result, several key residues were identified to be responsible for the binding of arachidonic acid with CYP4F22. These findings provide useful information for understanding the biological roles of CYP4F22 and structure-based drug design.