Understanding Microscopic Binding of Human ... - ACS Publications

May 14, 2008 - Microsomal prostaglandin E synthase-1 (mPGES-1) is a promising target for development of next-generation anti-inflammatory drugs...
2 downloads 0 Views 960KB Size
7320

J. Phys. Chem. B 2008, 112, 7320–7329

Understanding Microscopic Binding of Human Microsomal Prostaglandin E Synthase-1 with Substrates and Inhibitors by Molecular Modeling and Dynamics Simulation Adel Hamza, Mohamed Diwan M. AbdulHameed, and Chang-Guo Zhan* Department of Pharmaceutical Sciences, College of Pharmacy, UniVersity of Kentucky, 725 Rose Street, Lexington, Kentucky 40536 ReceiVed: January 25, 2008; ReVised Manuscript ReceiVed: March 25, 2008

Microsomal prostaglandin E synthase-1 (mPGES-1) is a promising target for development of next-generation anti-inflammatory drugs. It is crucial for rational design of the next-generation anti-inflammatory drugs to know the three-dimensional (3D) structure of mPGES-1 trimer and to understand how mPGES-1 binds with substrates and inhibitors. In the current work, a 3D structural model of human mPGES-1 trimer has been developed, for the first time, by performing combined homology modeling, molecular docking, and molecular dynamics simulation. The 3D structural model enables us to understand how mPGES-1 binds with its substrates/ inhibitors, and the key amino acid residues for the mPGES-1 binding with ligands have been identified. The detailed 3D structures and calculated binding free energies for mPGES-1’s binding with substrates and inhibitors are all consistent with available experimental data, suggesting that the 3D model of the mPGES-1 trimer and the enzyme-ligand binding modes are reasonable. The new structural insights obtained from this study should be valuable for rational design of next-generation anti-inflammatory drugs. Introduction Prostaglandin E2 (PGE2) is one of the most important prostanoids with diverse biological activity.1 The biosynthetic pathway of PGE2 has been well-characterized and involves three sequential enzymatic actions.2 The first step in this pathway involves the release of arachidonic acid (AA) from the membrane by the action of phospholipase A2 (PLA2).2 This is followed by the conversion of AA to prostaglandin H2 (PGH2, depicted in Figure 1) by the action of cyclooxygenase (COX)-1 or COX-2.2 Finally, PGH2 is converted to PGE2 by the action of terminal prostaglandin E synthase (PGES) enzymes.3 Research reported so far has identified three types of PGES; namely, microsomal prostaglandin E synthase-1 (mPGES-1), mPGES-2, and cPGES.4 Among the three, mPGES-1 has been identified as the therapeutically important enzyme. This enzyme plays a key role in various physiological processes.5 mPGES-1 is an inducible, membrane-bound enzyme.5 It is localized in the perinuclear membrane.2 mPGES-1 requires the tripeptide glutathione (GSH) as an essential cofactor for the catalytic activity.6 It is induced by proinflammatory stimuli such as interleukins and down-regulated by anti-inflammatory glucocorticoids.7 mPGES-1 was reported to be functionally coupled to COX-2 in preference to COX-1.3 A more recent report shows that mPGES-1 in vascular smooth muscles is not coupled to a particular COX isoenzyme.8 Studies on mPGES-1 have established its role in a number of disease conditions, including inflammation, arthritis, fever, pain, cancer, stroke, and bone disorders.9–15 Further, mPGES-1 knockout mice showed a reduction in the production of inflammatory PGE2 and a decrease in inflammatory response in a collagen-induced arthritis model.16 mPGES-1 has also been reported to be overexpressed in several types of cancers, including non-small-cell lung cancers.17 Upregulation of mPGES-1 has been reported in premalignant and malignant breast disease.18 The knockout studies have identified mPGES-1 as a central switch in immune-induced pyresis.19 * Phone: 859-323-3943. Fax: 859-323-3575. E-mail: [email protected]

Figure 1. Molecular structures of PGH2 (left side) and mPGES-1 inhibitors 1 to 6 (right side). In the inhibitor structures, X ) H (1), tBu (2), Ph (3), m-diphenyl (4), p-diphenyl (5), and 3-F, p-diphenyl (6).

mPGES-1-deficient mice were reported to be viable and fertile and develop normally when compared to the wild-type controls.16 PGE2 produced by mPGES-1 is a major mediator of inflammation and pain.19 The currently available nonsteroidal anti-inflammatory drugs (NSAIDS) in the market target either COX-2 (for selective COX-2 inhibitors) or COX-1 + COX-2 (for nonselective COX inhibitors).20 But these COX inhibitors have a number of side effects, such as ulcers and cardiovascular effects.21–23 More recently, two selective COX-2 inhibitors, that is, rofecoxib (Vioxx) and valdecoxib (Bextra), were withdrawn from the market due to severe side effect profiles. Grosser et al. showed that the biological basis for the cardiovascular side effects of COX-2 inhibitors is target-specific for COX-2.24 By inhibiting COX-2, the functions of all downstream PG synthases are blocked, including prostacyclin synthase (PGIS) for the conversion of PGH2 to prostaglandin I2 (PGI2). Blocking the production of PGI2 has been reported to play a role in cardiovascular side effects.25 In this situation, inhibition of the downstream enzyme mPGES-1 has emerged as a novel strategy that will target only the PGE2 pathway and without the adverse effects expected in the COX-2 inhibition.26 Thus, mPGES-1 is a very promising target for development of next-generation antiinflammatory drugs. Development of mPGES-1 inhibitors is also expected to be useful in treating other disease conditions, such as cancer. Although mPGES-1 inhibitors are expected to be potentially valuable therapeutic agents, very few inhibitors of mPGES-1 were identified in experimental screening efforts. The COX-2

10.1021/jp8007688 CCC: $40.75  2008 American Chemical Society Published on Web 05/14/2008

Binding of mPGES-1 with Ligands inhibitor NS-398, 5-lipoxygenase activating protein (FLAP) inhibitor MK-886, and the active metabolite of another NSAID (sulindac) were found to inhibit mPGES-1 with an IC50 of 20, 1.6, and 80 µM, respectively.27,28 Leukotriene C4 was reported to inhibit mPGES-1 with micromolar IC50, probably by competing with glutathione (GSH).27 In addition to small molecules, several polyunsaturated fatty acids and stable analogs of PGE2 were reported to inhibit mPGES-1.29 Riendeau et al.30 recently reported a series of mPGES-1 inhibitors. These compounds were synthesized on the basis of the scaffold of MK-886 (FLAP inhibitor). Some of these newly synthesized mPGES-1 inhibitors are potent, with an IC50 value of a few nanomolar in vitro. However, none of these inhibitors is specific for mPGES-1, such that all of these compounds have a very low activity against mPGES-1 in vivo. It is highly desirable to design and discover novel, mPGES-1-specific inhibitors with different scaffolds for development of a new generation of therapeutics. One of the important impeding factors for rational design and discovery of novel mPGES-1 inhibitors is the absence of a detailed three-dimensional (3D) crystal structure of mPGES-1. mPGES-1 is recognized as a member of a group of membraneassociated proteins involved in eicosanoid and glutathione metabolism (MAPEG) family.31 This protein family consists of membrane-bound proteins with diverse functions, such as leukotriene C4 synthase, FLAP, and microsomal glutathione transferase-1 (MGST1), MGST2, and MGST3.32 mPGES-1 shows significant homology with other MAPEG proteins.32 There is nearly 39% sequence identity between mPGES-1 and MGST1.33 These two proteins form a separate group (group I) under the MAPEG superfamily.32 Hydropathy analysis suggests that all the MAPEG proteins have similar 3D structures and membrane-spanning topographic properties.33 Because mPGES-1 is a membrane protein, it is generally difficult to produce a qualified single crystal for X-ray diffraction. But initial electron crystallography studies have revealed that mPGES-1 is organized as a homotrimer and should be similar to the homotrimer structure of MGST1.33 Various hydrodynamic studies further support the trimeric structure of mPGES-1.33 Site-directed mutagenesis studies have revealed that Arg-110 has an essential role in catalytic activity.7 mPGES-1 is also reported to catalyze the conjugation of 1-chloro-2,4-dinitrobenzene to GSH.33 This reaction shows the evolutionary relationship between mPGES-1 and MGST1.33 We recently reported a 3D structural model of the substrate-binding domain of mPGES-1.34 This first-ever 3D model of mPGES-1 was constructed by developing and using an ab initio approach of the protein structure prediction. The 3D model was also used to explore how mPGES-1 binds with substratePGH2 andglutathione(GSH).Thepredictedenzyme-substrate binding mode was examined by performing site-directed mutagenesis and the enzyme activity tests on various mPGES-1 mutants. The data obtained from the wet experimental tests are qualitatively consistent with the 3D model of mPGES-1. The studies also revealed that Tyr130 and Thr114 are essential for the activity of mPGES-1.34 However, without a 3D structural model of the mPGES-1 trimer, we were unable to refine the microscopic structure of the protein-substrate interaction through molecular dynamics (MD) simulation. Very recently, a crystal structure of MGST1 complexed with GSH became available.35 The 3D crystal structure reveals that the glutathione binding site in MGST1 is different from the cytosolic GSTs. It also shows the trimeric structure of MGST1. Given the ∼39% sequence identity between mPGES-1 and MGST1, this MGST1 structure serves as a good starting point for developing homology-based modeling of mPGES-1. Here,

J. Phys. Chem. B, Vol. 112, No. 24, 2008 7321 for the first time, we report a 3D structural model of the mPGES-1 trimer and its binding with substrates and inhibitors that is based on the homology modeling and further computational studies. The obtained structural and energetic data are in good agreement with available experimental biochemical data and provide valuable new insights into the detailed microscopic binding of mPGES-1 with its substrates and inhibitors. METHODS Construction of Initial 3D Models of mPGES-1. The amino acid sequence of mPGES-1 [AC: O14684] was generated from the GenBank database.36 The search for sequence similarities with several members of the MAPEG family within the Protein Data Bank (PDB)37 database was performed with the BLAST program.38 The crystal structure of the rat microsomal glutathione transferase 1 (MGST1) [ID: 2H8A]35 was used as the template to build the initial mPGES-1 model. The multiplesequence alignment was performed using the Homology module of InsightII software (Accelrys, Inc.). The 3D models of the transmembrane (TM) regions of the mPGES-1 protein were generated using the automated homology modeling tool Modeler (version 9v1)39–41 with default parameters, whereas the N- and C-terminal domains and the remaining loops were built de novo and refined to access proper folding with minimum steric clash by using the loop subroutine of Modeler. Atomic coordinates of the side chains were added from the standard residues library of Modeler. The Modeler is a well-known comparative modeling tool that generates a refined 3D homology model of a protein sequence automatically and rapidly, on the basis of a given sequence alignment to a known 3D protein structure template. It employs probability density functions (PDFs) as the spatial restraints rather than energy.39–41 The Modeler generates a large number of template-derived restraints to force the model toward the structure of the template and converges to the best possible structure by simultaneously satisfying a network of spatial restraints and molecular geometry using the CHARMM forcefield.42 The optimization process to generate the 3D model consists of applying the variable target function as well as conjugate gradients and molecular dynamics with simulated annealing. Thus, to construct a preliminary 3D structural model as an input of further refinement, a set of 150 conformations of the 3D structure were generated using the fast simulatedannealing procedure implemented in the Modeler. It is wellknown that MD simulations are valuable in the improvement of medium resolution (∼4 Å) structures and homology models. As well as restricted MD having been used in some solved structure refinement programs,43 it has also been used in the completion/refinement of homology models.44 Validation of the Models. In general, every homology model contains errors. The number of errors (for a given method) is mainly dependent on the percentage sequence identity between the template and target and on the number of errors in the template itself. An essential step in the homology modeling process is to verify/validate the model. We used several steps to estimate errors in our 3D models. First, the top 10 models obtained with the corresponding lowest PDFs were selected. The stereochemical quality of these structures was evaluated to determine whether the bond lengths and bond angles are within their normal ranges and whether there are lots of bumps in the models (corresponding to a high van der Waals energy). This was done by calculating the G-factor using Procheck.45,46 Second, the atomic contact analysis was also performed by using the program to identify bad packing of side chain atoms or unusual residue contacts. Then the normality indices (z-score)

7322 J. Phys. Chem. B, Vol. 112, No. 24, 2008 and the atomic contact analysis that describe how well a given characteristic of the model resembles the same characteristic in the template structure were determined by using Whatif.47 Finally, the 3D model with the best score (in terms of G-factor and z-score) was considered the best initial structure for further refinement as described below. Model Refinement. Previous structural comparison studies with DALI server48 revealed that MGST1 has a striking correspondence to subunit I of ba3-cytochrome c oxidase.35 It was reported that the TM regions, the connectivity between and the order of R-helices of subunit I of ba3-cytochrome c oxidase is in good agreement with that in the MGST1 structure. The report further demonstrated the similarity in the TM regions as well as in the connectivity by superimposition of the two structures.35 The root-mean-square deviation (rmsd) between the main-chain atoms in the corresponding TM regions of the molecule was reported as 2.04 Å. Thus, ba3-cytochrome c oxidase has a protein structure with a high structural similarity to MGST1. To build the mPGES-1 trimer, the above-validated 3D model obtained from the homology modeling was superimposed on subunit I of ba3-cytochrome c oxidase (PDB code: 1EHK), then the two other monomers were added in the same manner. Finally, the constructed mPGES-1 trimer model was solvated (with phospholipid bilayer and water) and subjected to a short run of energy minimization (3000 cycles) to relieve possibly unfavorable steric interactions and to optimize the stereochemistry (see below for the MD simulation procedure). This refined 3D model of mPGES-1 was again evaluated periodically for its stereochemical quality as well as its residue packing and atomic contacts, as described above. Molecular Docking with GSH, PGH2, and Inhibitors. The geometries of the ligands (i.e., cofactor GSH, substrate PGH2, and inhibitors depicted in Figure 1) were optimized by performing ab initio electronic structure calculations using the Gaussian03 program49 at the HF/6-31G* level. The optimized geometries were used to calculate the electrostatic potentials (at the same HF/6-31G* level) used to determine the partial atomic charges according to the standard restrained electrostatic potential (RESP) fitting procedure.50 The determined RESP charges were used in the calculations of electrostatic energy terms in the docking and MD simulation processes. The missing force field parameters for the ligands were taken from the General Amber Force Field (GAFF) implemented in the Amber9 program.51 There is a dilemma in accounting for solvent effects during molecular docking simulation. On one hand, it would be theoretically better to describe the solvent effects by using explicit water molecules. On the other hand, technically, the existence of water molecules within the active site can hinder the positioning of the ligand and, thus, decrease the efficiency of the conformational search during the molecular docking. To extensively search for the possible binding modes and conformations of ligands in the protein pocket, we decided to first use an implicit solvent model with a distance-dependent dielectric constant function (ε ) 4r)52–54 during the molecular docking. The obtained initial protein-ligand binding structures were subjected to energy minimizations and MD simulations on the explicitly solvated protein-ligand complexes (see below for the MD simulations) that can better account for the solvent effects on the protein-ligand binding. For GSH binding with the protein, the GSH molecule was added to the validated mPGES-1 trimer model by superimposing with the MGST1-GSH template structure. The refinement of the initial mPGES-1-GSH complex was obtained by performing

Hamza et al. energy minimization in vacuum (converged to 0.001 kcal mol-1 Å-1) followed by MD simulation for 100 ps. During this process, only the side chains of residues in the binding pocket were kept free to move. A 10.0 kcal mol-1 Å-2 positional constraint was applied to the other residues of the protein. The cutoff used for nonbonded interactions was 12 Å. The temperature of the system was maintained at 298 K with a 0.2 ps coupling constant. To explore the possible mode of mPGES-1 binding with a ligand (i.e., PGH2 in the presence of GSH or inhibitor), the first step was to identify by virtue the key residues in the interface between two mPGES-1 monomers in the trimer and to dock the ligand. We aimed to find where the ligand could be inserted most comfortably. The molecular docking for each ligand binding was carried out in the same way as we recently did for studying other protein-ligand binding systems.55,56 Briefly, a ligand-binding site was defined as that consisting of the residues at the interface of each mPGES-1 monomer and facing to the crucial Arg110 residue of a monomer in the mPGES-1 trimer. The initial docking calculations were performed on the ligand with the mPGES-1 trimer binding site using the “automatic docking” affinity module of the InsightII package (Accelrys, Inc.). The Affinity methodology uses a combination of Monte Carlo type and simulated annealing procedures to dock the guest molecule (ligand) to the host (the receptor; i.e., the mPGES-1 trimer).57 During the simulation, the side chains of residues in the specified binding site and PGH2 molecule were allowed to move while the others of the mPGES-1 trimer and the GSH cofactor were kept rigid. In the docking calculations for mPGES-1 binding with PGH2, 100 docking runs were carried out, with each run producing a single binding structure. For each binding structure, we estimated the receptor-ligand interaction energy as the sum of the electrostatic interaction energy and van der Waals interaction energy between the ligand and the mPGES-1 trimer, then the 100 possible binding structures obtained from all of the docking runs were first sorted geometrically in terms of their structural similarities. The possible binding structures having a rmsd value (in atomic positions) of -5.0). The overall root-mean-square deviation of the positions of the monomer backbone atoms in the energyminimized mPGES-1 trimer model from those in the crystal structure of MGST1 was found to be ∼0.7 Å, highlighting the overall similarity of the folding pattern, especially for the highly conserved TM domain (Figure S4 in the Supporting Information). Finally, the resulting mPGES-1 homotrimer model does not exhibit steric clashes involving main chains (see Figures S4 and S6 in the Supporting Information). To validate the stability of the structural model and to further test whether the modeling of this homotrimeric membrane protein based on the rat MGST1 structure is reasonable, we carried out the fully relaxed MD simulation on the mPGES-1 trimer structure in phospholipid bilayer surrounded by water. The MD simulation confirmed that the modeled mPGES-1 trimer structure was very stable (see Figure S5 in the Supporting Information for the time-dependence of the rmsd values), with all TM helices for each monomer structure remaining almost intact after the MD simulation for ∼7 ns. As expected, the TM helices of the trimer displayed a low rmsd fluctuation during the MD simulation, with a flat curve after ∼3.2 ns, indicating a minor structural change from the initial structure. Not surprisingly, the extramembraneous loops connecting with the

TM helices displayed relatively larger structural changes from the initial structure (data not shown). We also examined the rmsd of the TM helices for each monomer in the mPGES-1 trimer and obtained similar rmsd curves for the three monomers during the simulation on the trimer; after ∼3.2 ns, the rmsd values were converged to ∼1.5-1.9 Å (see the Supporting Information). In the MD-simulated mPGES-1 trimer structure, each monomer has a TM helix (residues 62-92) facing to the center of the trimer to contribute to the formation of a strong hydrophobic core (Figure S6 in the Supporting Information). The interactions within the core are characterized mainly by a network of dipole-quadrupole interactions between the phenyl ring of Phe84, Phe87, Phe90, and Tyr80 of one monomer and the equivalent residues of the other two monomers. In addition, the hydrophobic core of the trimer is also stabilized by the hydrogen bonds of the hydroxyl hydrogen of Tyr80, Tyr80′, and Tyr80′′ with the backbone oxygen atoms of the TM helices facing the core. Examination of the surface of the trimer revealed the presence of a crevice located at the interface between two neighboring monomers. Interestingly, superposition of the crystal structure of MGST1 with the modeled structure of mPGES-1 trimer brings the bound GSH into the crevice of the modeled trimer. A GSH molecule has been fitted in the crevice of the modeled trimer, in stacking interaction with Arg73. In addition to this interaction with Arg73, hydrogen bonds could be established with Asn74. Moreover, the Arg110 side chain forms the lower wall of the crevice between monomers and is the only charged residue to be largely accessible to the solvent in the mPGES-1 trimer. Following the exhaustive conformational sampling of our refined homology-built mPGES-1 trimer model, we assessed its quality and reliability by analyzing the generated 3D structure with respect to the available site-directed mutagenesis data34 and by examining simulated side chain interactions of the ligands’ (GSH and PGH2) binding sites. Binding of mPGES-1 with GSH and PGH2. Molecular docking revealed only five clusters of the mPGES1-GSH-PGH2 binding modes, but only the primary cluster is associated with the expected strong binding between PGH2 and the Arg110 side chain. This binding mode corresponds to the lowest interaction energy, so the binding of PGH2 with mPGES1-GSH can be regarded as strong and specific because it forms two hydrogen bonds with the Arg110 side chain. The mPGES1-GSH-PGH2 binding structure obtained from the docking was refined by performing MD simulation. The computational results clearly reveal how GSH and PGH2 bind with the receptor and predict the corresponding binding free energies. As shown in Figure 3, substrate PGH2 was docked into a cavity formed at the interface between neighboring monomers and near the GSH cofactor. The average PGH2 structure obtained from the

Binding of mPGES-1 with Ligands

J. Phys. Chem. B, Vol. 112, No. 24, 2008 7325

Figure 3. Ribbon representation of the structure of the mPGES-1 trimer existing in the phospholipid bilayer model (A). Solvent accessible surface area of the active site of the mPGES-1-GSH complex (B). Views from two different orientations of the binding mode from the MDsimulated structures of GSH and PGH2 in the catalytic active site of the mPGES-1 trimer (C and D). The ribbons of three monomers in the trimer were highlighted in different colors, and the H-bonds between the ligands and the key residues were highlighted in a dashed line.

MD simulation shows that the carboxylate moiety of the substrate faces the side chains of Arg110 and Thr114 and is nearly parallel to the scaffold of cofactor GSH. GSH is bound in an extended conformation to the side of the TM2 helix by forming several hydrogen bonds to the atoms of the protein side chain, which include contacts with Arg70 (Arg73), Arg73 (Leu76), and Asn74 (Asn77) of mPGES-1 (MGST1). Such a conformation was also observed in other GSH binding protein.76,77 The Arg73 residue of mPGES-1 is unique, whereas the other two residues are highly conserved among the members of the MAPEG family. We conclude that the architecture of the active site of mPGES-1 may be distinctly different from that of others GSTs. Indeed, in the crystal structure of MGST1, GSH is associated with a monomer through only one hydrogen bond between the γ-glutamyl carboxyl group of the ligand and the Arg73 side chain, whereas the glycine carboxyl moiety of the GSH is in a relatively unfavorable environment including Glu80 residue. In contrast, the extended conformation of GSH in mPGES-1 is maintained by the hydrogen bonds between the two carboxyl groups of the cofactor and the side chains of the Arg73 and Arg70 residues. The modeled binding mode is consistent with available experimental observations reported in the literature. For ex-

ample, our previous site-directed mutagenesis studies34 have demonstrated that the enzymatic activity of mPGES-1 is mainly influenced by mutations on the Arg110, Thr114, and Tyr130 residues.34 Arg110 and Thr114 are essential to confer the mPGES-1 activity. An analysis of the disposition of these residues in the mPGES-1 structure reveals that the Arg110, Thr114, and Tyr130 residues are all in close interaction with PGH2. The average PGH2 structure from the MD simulation reveals an interaction of its carboxyl moiety with the guanidinium side chain of Arg110, whereas the ω-chain of PGH2 faces the core of the mPGES-1 trimer (Figure 3, also Figure S7 in the Supporting Information), the stability of the model was evaluated by monitoring some key distances between the key residues of mPGES-1 and its cofactor GSH with substrate PGH2. Thus, the key residues in the PGH2 binding site include Arg110 and Thr129 of the first monomer, and Leu69′ and His72′ of the second monomer (i.e., residues Leu69 and His72 of another monomer in the trimer). In the bound PGH2, the C-terminal carboxyl group of the R chain forms the salt bridges to Arg110, and the side chain of the latter is stabilized in an extended conformation interacting with the His113 side chain. In our 3D model, both Arg110 and His113 side chains are stabilized by a π-π stacking interaction and are involved in a network of hydrogen bonds with the PGH2 carboxyl group, Thr129, and Thr114 residues. Indeed, the hydroxyl group of Thr114 interacts strongly with the Nδ1 atom of His113 through a hydrogen bond, while the Hε2 atom of this latter contributes to a hydrogen bond with the Thr129 hydroxyl group. As depicted in Figure 3, because Arg110 is involved in an ionic interaction with the PGH2 carboxyl group, it is clear that the Arg110Thr mutation drastically decreases the affinity of substrate PGH2 with mPGES-1 due to the shift of the bidentate salt bridge interaction (Arg110NH1/NH2-PGH2carboxyl) to a monodental bifurcated salt bridge (Thr129OH-PGH2carboxyl). The results from the MM/PBSA binding free energy calculation for the Arg110Thr mutant are given in Table 1. The electrostatic component of ∆EMM is weaker than in the wild type (WT), while the van der Waals components are similar to the WT. The net values of ∆G bind are thus -3.1 for the mutant and -6.4 kcal/ mol for the WT. The mutation-caused binding free energy change is 3.3 kcal/mol. Thus, the change in the binding free energy is associated with the loss of the electrostatic interaction. Alignment of the known MAPEG family members78 demonstrates that Arg110 in mPGES-1 is strictly conserved. As demonstrated in previous experimental studies, replacement of Arg110 by Thr or Ser indeed significantly decreases the catalytic function of mPGES-1, implying a crucial role of this residue.7,34 In addition, the Thr114Val mutation significantly decreases the interaction between PGH2 and mPGES-1, which is consistent

TABLE 1: Energetic Results (kcal/mol) Obtained from MM/PBSA Binding Free Energy Calculations on PGH2 Binding with Wild-Type (WT) mPGES-1 and Its Mutants in the Presence of GSH

a

energy

mPGES-1-WT

mPGES-1-R110T

mPGES-1-T114V

mPGES-1-Y130I

ele ∆Eint vdW ∆Eint ∆EMM ∆Esol ele ∆Etot ∆Ebind -T∆S ∆Gbind exptla

-381.5 -49.4 -430.9 400.1 23.1 -30.7 24.3 -6.4 -6.6 to -5.2

-326.3 -53.6 -379.9 352.5 30.6 -27.4 24.3 -3.1

-374.7 -50.8 -425.5 399.3 29.09 -26.2 24.3 -1.9

-374.0 -47.9 -421.9 394.6 24.8 -27.3 24.3 -3.0

Experimental values from ref 34.

7326 J. Phys. Chem. B, Vol. 112, No. 24, 2008

Hamza et al.

Figure 4. Ribbon diagram for the binding mode of the MD-simulated structure of compound 3 in the active site of the mPGES-1 trimer (A). Superposition of the mPGES-1 binding site complexed with the MDsimulated structures of compounds 1-3 (B) and 4-6 (C). The H-bonds between the ligands and key residues of mPGES-1 are highlighted in a dashed line. Plots of MD-simulated internuclear distances and rmsd for atomic positions of the ligands 6 versus the simulation time for mPGES-1 trimer (D). Traces D1 and D2 represent the H-bond distances between the ligand carboxyl group and the Arg110 guanidinium and T129 hydroxyl side chains, respectively. Trace D3 represents the distance between the indole ring of the ligand and the imidazole ring of the H72′ side chain. Trace D4 describes the distance of the planes between the Cl-phenyl ring of the ligand and the guanidinium Arg73 side chain. Trace D5 represents the H-bond distance between the hydrogen HH12 atom of the R126 side chain and the fluorine atom of compound 6.

with the experimental observation.34 This is most likely due to the steric hindrance induced by a bulky Val side chain, possibly accompanied by a charge-charge interaction. As shown in Table 1, the interaction energy, ∆EMM decreases by ∼5 kcal/mol for the Thr114Val mutant due to the effect of the charge repulsion. The binding free energy calculations show that PGH2 has a lower binding affinity to the Thr114Val mutant. The Tyr130Ile mutation also decreases the binding affinity of PGH2 to mPGES-1. The decrease of the binding affinity is mainly due to a change in the electrostatic component of the ∆EMM term, as seen in Table 1. Indeed, a structural analysis of the simulated mutant showed a steric interaction of the Ile130 side chain with the carboxyl group of the PGH2, thus destroying the salt bridge between the PGH2 and the Arg110 residue. In the wild type, the side chain of Tyr130 stays toward the solvent,

allowing the formation of a tightly packet hydrophobic and dipole-quadrupole interactions between the aromatic side chains of Tyr130 and Ala31′ (i.e., the Ala31 residue of the second monomer). In addition, the ω-chain of the PGH2 moved slightly away from the side chain of His72′ (i.e., the His72 residue of the second monomer), thus preventing the formation of the π-π interaction found in the wild-type mPGES-1-GSH-PGH2 binding complex. Binding of mPGES-1 with Inhibitors. The modeled 3D structure of the mPGES-1 trimer was also used to study how mPGES-1 binds with two sets of known inhibitors depicted in Figure 1. The first set (denoted by group I) consists of compounds 1-3, whereas the second set (denoted by group II) consists of compounds 4-6.30 All of these compounds are derivatives of compound MK-886. The MD simulations on the six solvated mPGES-1-ligand complexes led to six average complexes, one average structure for each complex. The obtained average structures revealed that the carboxyl moiety of the compounds in both groups I and II make salt bridges with the Arg110 side chain and hydrogen bond with Thr129, whereas the indole ring (scaffold) binds at the interface between the first two monomers in the homotrimeric mPGES-1 enzyme. Inspection of the MD trajectory revealed that the rmsd of the ligand (compound 1, 2, or 3) heavy atoms from the initial structure was found to be stable after ∼250 ps of the MD simulation, and the rmsd value was rather small, ∼0.4 Å. The stability of the structures obtained from the MD simulations is illustrated as figures in the Supporting Information for compounds 1 and 3. The ligand fits into the pocket by establishing strong electrostatic interactions with the guanidinium of Arg110 and hydroxyl group of Thr129 side chains, then the indole ring is stabilized in the pocket by two π-π stacking interactions with His72′ (i.e., His72 of the second monomer in the trimer) and Tyr117 residues, and the Cl-phenyl group plugs into a positively charged pocket surrounded by the side chains of Arg73, Arg73′ (i.e., Arg73 of the second monomer in the trimer), Arg73′′ (i.e., Arg73 of the third monomer in the trimer), and Met76′ (i.e., Met76 of the second monomer in the trimer) (Figure 4). In addition, two significant dipole-quadrupole interactions are also seen between the indole ring and residues Arg126 and Arg73. Being consistent with these observations, the substitution of the Cl-phenyl group by a hydrogen or methyl group indeed decreases the binding affinity of compound 1.30 The MM/PBSA-calculated binding free energies are summarized in Table 2. When compared to mPGES-1 complexed with compound 1, the more favorable binding of compounds 2 and 3 with mPGES-1 is attributed to the molecular mechanics interaction energy component (∆EMM) and the nonpolar contribution to the solvent effects (∆Esol). The MM/PBSA results suggest that a crucial factor affecting the binding affinity is to

TABLE 2: Energetic Results (kcal/mol) Obtained from MM/PBSA Binding Free Energy Calculations on mPGES-1 Binding Inhibitors with a 1 to 6.

a

energy

mPGES-1-1

mPGES-1-2

mPGES-1-3

mPGES-1-4

mPGES-1-5

mPGES-1-6

∆Eeleint ∆EvdWint ∆EMM ∆Esol ∆Eeletot ∆Ebind -T∆S ∆Gbind exptla

-401.9 -45.1 -447.0 421.3 23.4 -25.7 17.6 -8.1 -7.5

-403.5 -48.5 -452.0 423.7 24.9 -28.3 19.7 -8.6 -8.84

-409.7 -48.0 -457.7 430.3 25.2 -27.4 19.2 -8.2 -8.49

-408.9 -52.6 -461.5 429.5 25.8 -32.0 23.0 -9.0 -9.27

-411.0 -53.1 -464.1 430.7 24.7 -33.4 23.0 -10.4 -10.6

-413.8 -53.2 -467.0 432.3 23.8 -34.7 23.1 -11.6 -11.11

Values converted from the experimental IC50 reported in ref 30.

Binding of mPGES-1 with Ligands achieve an optimal electrostatic interaction between the ligand and the protein active site and also to suffer less (or keep the same) desolvation penalty. The similar molecular docking and MD simulations were also performed to examine mPGES-1 binding with compounds 4-6. Each obtained enzyme-ligand structure can achieve a good equilibrium after ∼250 ps of MD simulation, giving a rmsd value of ∼0.5 Å for the positions of all ligand atoms as compared to those in the starting structure, as shown in Figure 4D. The ligand binding mode can be characterized by some internuclear distances between some ligand atoms and key amino acid residues. The results obtained from the MD simulations show that the binding modes of compounds 4-6 are similar to that of the compounds in group I. Thus, the side chain of the highly conserved Arg110 is the primary anchor point for the carboxyl group of the ligand as described by the distance D1 in Figure 4D, whose interaction is further stabilized by a hydrogen bond (distance D2) with the hydroxyl group of Thr129. These two hydrogen bonds are persistent during the MD simulation, with an average length of ∼1.7 and ∼1.9 Å for D1 and D2, respectively. Beside these important hydrogen bonds, several favorable interactions between the ligand and the enzyme are clearly demonstrated in Figure 4C. Similar to the binding of compounds 1-3, the Cl-phenyl ring and the indole group in compounds 4-6 stay in the same positively charged pocket, while the biphenyl ring is stabilized in a “sandwich” conformation by establishing a network of π-π interactions with Arg70 and Arg126 residues. Moreover, the ligands are also stabilized by forming another π-π interaction between the Tyr117 side chain and the indole ring of the ligand. Regarding the interaction with the phenyl ring in compounds 4 and 5, the para substitution increases the interaction energy term ∆EMM by ∼2 kcal/mol, giving rise to a larger binding affinity of compound 5 as compared to compound 4. Further introducing a fluorine atom at the meta position increases the ∆EMM value by establishing a favorable electrostatic interaction between the NH2 group of Arg126 and the fluorine atom in compound 6. The binding free energies calculated for compounds 4, 5, and 6 are -9.0, -10.4, and -11.6 kcal/mol, respectively. According to the energy components of the binding free energies (Table 2), the major favorable contributors to ligand binding are van der Waals and electrostatic interaction terms, whereas the solvation and entropy terms oppose the binding. Furthermore, it is encouraging to note that the order of the experimental affinities of all the six inhibitors is qualitatively consistent with that of our calculated binding free energies, suggesting that the determined enzyme-inhibitor binding structures are reasonable. CONCLUSIONS A comprehensive 3D structural model of human microsomal prostaglandin E synthase-1 homotrimer has been developed by carrying out combined homology modeling, molecular docking, and molecular dynamics simulation. The 3D structural model enables us to clearly identify the major differences and similarities between the mPGES-1 and MGST1 structures and to study how mPGES-1 binds with its substrates and inhibitors, and the key amino acid residues for the mPGES-1 binding with ligands have been identified. The simulated enzyme-ligand binding includes hydrogen bonding and aromatic interactions between mPGES-1 and ligands. The detailed 3D structures and calculated binding free energies for mPGES-1 binding with substrates and inhibitors are all consistent with available experimental data, suggesting that the 3D model of the mPGES-1 trimer and the

J. Phys. Chem. B, Vol. 112, No. 24, 2008 7327 enzyme-ligand binding modes are reasonable. The new structural insights obtained from this study should be valuable for rational design of next-generation anti-inflammatory drugs. Acknowledgment. The research was supported by the Kentucky Science and Engineering Foundation (Grant 925RDE-008 to C.G.Z.). The authors also acknowledge the Center for Computational Sciences (CCS) at the University of Kentucky for supercomputing time on an IBM X-series Cluster with 34 nodes and 1360 processors. Supporting Information Available: Twelve figures concerning more detailed information about the homology modeling, MD-simulated structure of mPGES-1 trimer binding with PGH2 and GSH, and MD trajectories for mPGES-1 trimer binding with inhibitors. This information is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Serhan, C. N.; Levy, B. Success of prostaglandin E2 in structurefunction is a challenge for structure-based therapeutics. Proc. Natl. Acad. Sci. U.S.A 2003, 100, 8609–8611. (2) Kudo, I.; Murakami, M. Prostaglandin E synthase, a terminal enzyme for prostaglandin E2 biosynthesis. J. Biochem. Mol. Biol. 2005, 38, 633–638. (3) Fahmi, H. mPGES-1 as a novel target for arthritis. Curr. Opin. Rheumatol. 2004, 16, 623–627. (4) Park, J. Y.; Pillinger, M. H.; Abramson, S. B. Prostaglandin E2 synthesis and secretion: The role of PGE2 synthases. Clin. Immunol 2006, 119, 229–240. (5) Sampey, A. V.; Monrad, S.; Crofford, L. J. Microsomal prostaglandin E synthase-1: the inducible synthase for prostaglandin E2. Arthritis Res. Ther. 2005, 7, 114–117. (6) Ouellet, M.; Falgueyret, J-P.; Ear, P. H.; Pen, A.; Mancini, J. A.; Riendeau, D.; Percival, M. D. Purification and characterization of recombinant microsomal prostaglandin E synthase -1. Protein Expr. Purif. 2002, 26, 489–495. (7) Murakami, M.; Naraba, H.; Tanioka, T.; Semmyo, N.; Nakatani, Y.; Kojima, F.; Ikeda, T.; Fueki, M.; Ueno, A.; Oh-ishi, S.; Kudo, I. Regulation of prostaglandin E2 biosynthesis by inducible membrane associated prostaglandin E2 synthase that acts in concert with cyclooxygenase-2. J. Biol. Chem. 2000, 275, 32783–32792. (8) Camacho, M.; Gerboiles, E.; Escudero, J.-R.; Anton, R.; GarciaMoll, X.; Vila, L. Microsomal prostaglandin E synthase-1, which is not coupled to a particular cyclooxygenase isoenzyme, is essential for prostaglandin E2 biosynthesis in vascular smooth muscle cells. J. Thromb. Haemost. 2007, 5, 1411–1419. (9) Uematsu, S.; Matsumoto, M.; Takeda, K.; Akira, S. Lipopolysaccharide-dependent prostaglandin E2 production is regulated by the glutathione-dependent prostaglandin E2 synthase gene induced by Toll-like receptor 4/MyD88/NF-IL6 pathway. J. Immunol. 2002, 168, 5811–5816. (10) Kaemi, D.; Murakami, M.; Nakatani, Y.; Ishikawa, Y.; Ishi, T.; Kudo, I. Potential role of microsomal prostaglandin E synthase-1 in tumorigenesis. J. Biol. Chem. 2003, 278, 19396–19405. (11) Kaemi, D.; Yamakawa, K.; Takegoshi, Y.; Mikami-Nakanishi, M.; Nakatani, Y.; Oh-ishi, S.; Yasui, H.; Azuma, Y.; Hirasawa, N.; Ohuchi, K.; Kawaguchi, H.; Ishikawa, Y.; Ishii, T.; Uematsu, S.; Akira, S.; Murakami, M.; Kudo, I. Reduced pain hypersensitivity and inflammation in mice lacking microsomal prostaglandin E synthase-1. J. Biol. Chem. 2004, 279, 33684–33695. (12) Ikeda-Matsuo, Y.; Ota, A.; Fukada, T.; Uematsu, S.; Akira, S.; Sasaki, Y. Microsomal prostaglandin E2 synthase-1 is a critical factor of stroke-reperfusion injury. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 11790– 11795. (13) Murakami, M.; Kudo, I. Recent advances in the molecular biology and physiology of the prostaglandin E2 biosynthetic pathway. Prog. Lipid Res. 2004, 43, 3–35. (14) Claveau, D.; Sirinyan, M.; Guay, J.; Gordon, R.; Chan, C.-C.; Bureau, Y.; Riendeau, D.; Mancini, J. A. Microsomal prostaglandin E2 synthase is a major terminal synthase that is selectively up-regulated during cyclooxygenase-2 dependent prostaglandin E2 production in the rat adjuvantinduced arthritis model. J. Immunol. 2003, 170, 4738–4744. (15) Oshima, H.; Oshima, M.; Inaba, K.; Taketo, M. M. Hyperplastic gastric tumors induced by activated macrophages in COX-2/mPGES-1 transgenic mice. EMBO J. 2004, 23, 1669–1678. (16) Trebino, C. E.; Stock, J. L.; Gibbons, C. P.; Naiman, B. M.; Wachtmann, T. S.; Umland, J. P.; Pandher, K.; Lapointe, J.-M.; Saha, S.;

7328 J. Phys. Chem. B, Vol. 112, No. 24, 2008 Roach, M. L.; Carter, D.; Thomas, N. A.; Durtschi, B. A.; McNeish, J. D.; Hambor, J. E; Jakobsson, P.-J.; Carty, T. J; Perez, J. R; Audoly, L. P. Impaired inflammatory and pain response in mice lacking an inducible prostaglandin E synthase. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 9044– 9049. (17) Yoshimatsu, K.; Altorki, N. K.; Golijanin, D.; Zhang, F.; Jakobsson, P-J.; Dannenberg, A. J.; Subbaramaiah, K. Inducible prostaglandin E synthase is overexpressed in non-small cell lung cancer. Clin. Cancer Res. 2001, 7, 2669–2674. (18) Mehrotra, S.; Morimiya, A.; Agarwal, B.; Konger, R.; Badve, S. Microsomal prostaglandin E2 synthase-1 in breast cancer: a potential target for therapy. J. Pathol. 2006, 208, 356–363. (19) Engblom, D.; Saha, S.; Engstrom, L.; Westman, M.; Audoly, L. P.; Jakobsson, P. J.; Blomqvist, A. Microsomal prostaglandin E synthase-1 is the central switch during immune-induced pyresis. Nat. Neurosci. 2003, 6, 1137–1138. (20) Simmons, D. L.; Botting, R. M.; Hla, T. Cyclooxygenase isozymes: The biology of prostaglandin synthesis and inhibition. Pharmacol. ReV. 2004, 56, 387–437. (21) Silverstein, F. E.; Faich, G.; Goldstein, J. L.; Simon, L. S.; Pincus, T.; Whelton, A.; Makuch, R.; Eisen, G.; Agrawal, N. M.; Stenson, W. F.; Burr, A. M.; Zhao, W. W.; Kent, J. D.; Lefkowith, J. B.; Verburg, K. M.; Geis, G. S. Gastrointestinal toxicity with celecoxib vs nonsteroidal antiinflammatory drugs for osteoarthritis and rheumatoid arthritis. The CLASS study: A randomized controlled trial. J. Am. Med. Assoc. 2000, 284, 1247– 1255. (22) Caldwell, B.; Aldington, S.; Weatherall, M.; Shirtcliffe, P.; Beasley, R. Risk of cardiovascular events and celecoxib: a systematic review and meta-analysis. J. R. Soc. Med. 2006, 99, 132–140. (23) Dickman, A.; Ellershaw, J. NSAIDS: gastroprotection or selective COX-2 inhibitor. PalliatiVe Med. 2004, 18, 275–286. (24) Grosser, T.; Fries, S.; FitzGerald, G. A. Biological basis of the cardiovascular consequences of COX-2 inhibition: therapeutic challenges and opportunities. J. Clin. InVest. 2006, 116, 4–15. (25) Cheng, Y.; Wang, M.; Yu, Y.; Lawson, J.; Funk, C. D.; FitzGerald, G. A. Cyclooxygenases, microsomal prostaglandin E synthase-1, and cardiovascular function. J. Clin. InVest. 2006, 116, 1391–1399. (26) Mehrotra, S.; Morimiya, A.; Agarwal, B.; Konger, R.; Badve, S. Microsomal prostaglandin E2 synthase-1 in breast cancer: a potential target for therapy. J. Pathol. 2006, 208, 356–363. (27) Thoren, S.; Jakobsson, P. J. Coordinate up- and down regulation of glutathione-dependent prostaglandin E synthase and cyclooxygenase-2 in A549 cells. Inhibition by NS-398 and leukotriene C4. Eur. J. Biochem. 2000, 267, 6428–6434. (28) Mancini, J. A.; Blood, K.; Guay, J.; Gordon, R.; Claveau, D.; Chan, C. C.; Riendeau, D. J. Cloning, expression and up-regulation of inducible rat prostaglandin E synthase during lipopolysaccharide-induced pyresis and adjuvant-induced arthritis. J. Biol. Chem. 2001, 276, 4469–4475. (29) Quraishi, O.; Mancini, J. A.; Riendeau, D. Inhibition of inducible prostaglandin E2 synthase by 15-deoxy-D12,14-prostaglandin J2 and polyunsaturated fatty acids. Biochem. Pharmacol. 2002, 63, 1183–1189. (30) Riendeau, D.; Aspiotis, R.; Ethier, D.; Gareau, Y.; Grimm, E. L.; Guay, J.; Gurial, S.; Juteau, H.; Mancini, J. A.; Methot, N.; Rubin, J.; Friesen, R. W. Inhibitors of the inducible microsomal prostaglandin E2 synthase (mPGES-1) derived from MK-886. Bioorg. Med. Chem. Lett. 2005, 15, 3352–3355. (31) Jakobsson, P. J.; Thoren, S.; Morgenstern, R.; Samuelsson, B. Identification of human prostaglandin E synthase: A microsomal, glutathione-dependent, inducible enzyme, constituting a potential novel drug target. Proc. Natl. Acad. Sci. U.S.A. 1999, 96, 7220–7225. (32) Murakami, M.; Kudo, I. Prostaglandin E Synthase: A novel drug target for inflammation and cancer. Curr. Pharm. Des. 2006, 12, 943–954. (33) Thoren, S.; Weinander, R.; Saha, S.; Jegerschold, C.; Pettersson, P. L.; Samuelsson, B.; Hebert, H.; Hamberg, M.; Morgenstern, R.; Jakobsson, P.-J. Human microsomal prostaglandin E synthase-1. Purification, functional characterization, and projection structure determination. J. Biol. Chem. 2003, 278, 22199–22209. (34) Huang, X.; Yan, W.; Gao, D.; Tong, M.; Tai, H. H.; Zhan, C.-G. Structural and functional characterization of human microsomal prostaglandin E synthase-1 by computational modeling and site directed mutagenesis. Bioorg. Med. Chem. 2006, 14, 3553–3562. (35) Holm, P. J.; Bhakat, P.; Jegerschold, C.; Gyobu, N.; Mitsuoka, K.; Fujiyoshi, Y.; Morgenstern, R.; Hebert, H. Structural basis for detoxification and oxidative stress protection in membranes. J. Mol. Biol. 2006, 360, 934– 945. (36) Benson, D. A.; Karsch-Mizrachi, I.; Lipman, D. J.; Ostell, J.; Wheeler, D. L. Genbank. Nucleic Acids Res. 2000, 28, 15–8. (37) Bernstein, F. C.; Koetzle, T. F.; Williams, G. J.; Meyer, E. E., Jr.; Brice, M. D.; Rodgers, J. R. The Protein Data Bank. J. Mol. Biol. 1977, 112, 535–42. (38) Altschul, S. F.; Gish, W.; Miller, W.; Myers, E. W.; Lipman, D. J. Basic local alignment search tool. J. Mol. Biol. 1990, 215, 403–10.

Hamza et al. (39) Fiser, A.; Do, R. K.; Sali, A. Modeling of loops in protein structures. Protein Sci. 2000, 9, 1753–1773. (40) Marti-Renom, M. A.; Stuart, A. C.; Fiser, A.; Sanchez, R.; Melo, F.; Sali, A. Comparative protein structure modeling of genes and genomes. Annu. ReV. Biophys. Biomol. Struct. 2000, 29, 291–325. (41) Sali, A.; Blundell, T. L. Comparative protein modelling by satisfaction of spatial restraints. J. Mol. Biol. 1993, 234, 779–815. (42) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D. ; States, D. J.; Swaminathan, S.; Karplus, M. CHARMM: A program for macromolecular energy minimization and dynamics calculations. J. Comput. Chem. 1983, 4, 187–217. (43) Sutton, R. B.; Rasshauer, D.; Jahn, R.; Brunger, A. T. Crystal structure of a SNARE complex involved in synaptic exocytosis at 2.4 A resolution. Nature 1998, 395, 347–353. (44) Flohil, J. A.; Vriend, G; Berendsen, H. J. C. Completion and refinement of 3-D homology models with restricted molecular dynamics: Application to targets 47, 58, and 111 in the CASP modeling competition and posterior analysis. Proteins Struct. Funct. Genet. 2002, 48, 593–604. (45) Laskowski, R. A.; MacArthur, M. W.; Moss, D. S.; Thornton, J. M. PROCHECK: a program to check the stereochemical quality of protein structures. J. Appl. Crystallogr. 1993, 26, 283–291. (46) Morris, A. L.; MacArthur, M. W.; Hutchinson, E. G.; Thornton, J. M. Stereochemical quality of protein structure coordinates. Proteins 1992, 12, 345–364. (47) Vriend, G. WHAT IF: a molecular modeling and drug design program. J. Mol. Graph. 1990, 8, 52. (48) Holm, L.; Sander, C. Protein structure comparison by alignment of distance matrices. J. Mol. Biol. 1993, 233, 123–138. (49) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Zakrzewski, V. G.; Montgomery, J. A.; Stratmann, R. E.; Burant, J. C.; Dapprich, S.; Millan, J. M.; Daniels, A. D.; Kudin, K. N.; Strain, M. C.; Farkas, O.; Tomasi, J.; Barone, V.; Cossi, M.; Cammi, R.; Mennucci, B.; Pomelli, C.; Adamo, C.; Clifford, S.; Ochterski, J.; Petersson, G. A.; Ayala, P. Y.; Cui, Q.; Morokuma, K.; Malich, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Cioslowski, J.; Ortiz, J. V.; Baboul, A. G.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Gonzales, C.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Andreas, J. L.; HeadGordon, M.; Reploge, E. S.; Pople, J. A. Gaussian 03, reVision B-03; Gaussian, Inc.: Pittsburgh, PA, 2003. (50) (a) Cieplak, P.; Cornell, W. D.; Bayly, C. I.; Kollman, P. A. Application of the multimolecule and multiconformational RESP methodology to biopolymers: charge derivation for DNA, RNA, and proteins. J. Comput. Chem. 1995, 16, 1357. (b) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. A well-behaved electrostatic potential based method using charge restraints for deriving atomic charges: the RESP model. J. Phys. Chem. 1993, 97, 10269. (51) Case, D. A.; Darden, T. A.; Cheatham, T. E., III; Simmerling, C. L.; Wang, J.; et al. AMBER 9, University of California, San Francisco, CA, 2006. (52) Pearlman, D. A. et al. AMBER 4.1. University of California, San Francisco, CA. 1995. (53) Harvey, S. C. Treatment of electrostatic effects in macromolecular modeling. Proteins 1989, 5, 78–92. (54) Guenot, J.; Kollman, P. A. Molecular dynamics studies of a DNAbinding protein: 2. An evaluation of implicit and explicit solvent models for the molecular dynamics simulation of the Escherichia coli trp repressor. Protein Sci. 1992, 1, 1185–1205. (55) Hamza, A.; Zhan, C-G. How Can (-)-epigallocatechin gallate from green tea prevent HIV-1 infection? mechanistic insights from computational modeling and the implication for rational design of anti-HIV-1 entry inhibitors. J. Phys. Chem. B 2006, 110, 2910–2917. (56) Bargagna, M. P.; Hamza, A.; Kim, Y. E.; Khuan, Y.; Ho, A.; MorVaknin, N.; Wendschlag, N.; Liu, J.; Evans, R. M.; Markovitz, D. M.; Zhan, C.-G.; Kim, K. B.; Mohan, R. The tumor inhibitor and antiangiogenic agent withaferin A targets the intermediate filament protein vimentin. Chem. Biol. 2007, 14, 623–34. (57) Kuntz, I. D.; Meng, E. C.; Shoichet, B. K. Structure-based molecular design. Acc. Chem. Res. 1994, 27, 117. (58) Jorgensen, W. L. Quantum and statistical mechanical studies of liquids. 10. Transferable intermolecular potential functions for water, alcohols, and ethers. Application to liquid water. J. Am. Chem. Soc. 1981, 103, 335. (59) Pan, Y.; Gao, D.; Yang, W.; Cho, H.; Yang, G.-F.; Tai, H.-H.; Zhan, C.-G. Computational redesign of human butyrylcholinesterase for anticocaine medication. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 16656. (60) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: an N · log(N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089. (61) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684.

Binding of mPGES-1 with Ligands (62) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. Numerical integration of the Cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys. 1977, 23, 327. (63) Kollman, P.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Duan, Y.; Wang, W.; Donini, O.; Cieplak, P.; Srnivasan, J.; Case, D.; Cheatham, T. E. Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models. Acc. Chem. Res. 2000, 33, 889–897. (64) Kollman, P. A.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Lee, T.; Duan, Y.; Wang, W.; Donini, O.; Cieplak, P.; Srinivasan, J.; Case, D. A.; Cheatham, T. E. Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models. Acc. Chem. Res. 2000, 33, 889–897. (65) Hou, T.; Guo, S.; Xu, X. Predictions of binding of a diverse set of ligands to Gelatinase-A by a combination of molecular dynamics and continuum solvent models. J. Phys. Chem. B 2002, 106, 5527–5535. (66) Xiong, Y.; Li, Y.; He, H.; Zhan, C. G. Theoretical calculation of the binding free energies for pyruvate dehydrogenase E1 binding with ligands. Bioorg. Med. Chem. Lett. 2007, 17, 5186–5190. (67) Masukawa, K. M.; Kollman, P. A.; Kuntz, I. D. Investigation of neuraminidase-substrate recognition using molecular dynamics and free energy calculations. J. Med. Chem. 2003, 46, 5628–5637. (68) Huo, S.; Wang, J.; Cieplak, P.; Kollman, P. A.; Kuntz, I. D. Molecular dynamics and free energy analyses of cathepsin D-inhibitor interactions: insights into structure-based ligand design. J. Med. Chem. 2002, 45, 1412–1419. (69) Donini, O. A. T.; Kollman, P. A. Calculation and prediction of binding free energies for the matrix metalloproteinases. J. Med. Chem. 2000, 43, 4180–4188. (70) Rafi, S. B.; Cui, G.; Song, K.; Cheng, X.; Tonge, P. J.; Simmerling, C. Insight through molecular mechanics Poisson-Boltzmann surface area

J. Phys. Chem. B, Vol. 112, No. 24, 2008 7329 calculations into the binding affinity of triclosan and three analogues for Fab1, the E. coli enoyl reductase. J. Med. Chem. 2006, 49, 4574–4580. (71) Wang, J.; Kang, X.; Kuntz, I. D.; Kollman, P. A. Hierarchical database screenings for HIV-1 reverse transcriptase using a pharmacophore model, rigid docking, solvation docking and MM-PB/SA. J. Med. Chem. 2005, 48, 2432–2444. (72) Wang, J.; Morin, P.; Wang, W.; Kollman, P. A. Use of MM-PBSA in reproducing the binding free energies to HIV-1 RT of TIBO derivatives and predicting the binding mode to HIV-1 RT of efavirenz by docking and MM-PBSA. J. Am. Chem. Soc. 2001, 123, 5221–5230. (73) Abdulhameed, M. D. M.; Hamza, A.; Zhan, C.-G. Microscopic modes and free energies of 3-phosphoinositide dependent kinase-1 (PDK1) binding with celecoxib and other inhibitors. J. Phys. Chem. B 2006, 110, 26365–26374. (74) Williamson, M. P.; McCormick, T. G.; Nance, C. L.; Shearer, W. T. Epigallocatechin gallate, the main polyphenol in green tea, binds to the T-cell receptor, CD4: Potential for HIV-1 therapy. J. Allergy Clin. Immunol. 2006, 118, 1369–1374. (75) Fujiwara, S.; Amisaki, T. Identification of high affinity fatty acid binding sites on human serum albumin by MM-PBSA method. Biophys. J. 2008, 94, 95–103. (76) Oakley, A. J. Glutathione transferases: new functions. Curr. Opin. Struct. Biol. 2005, 15, 716–723. (77) Ji, X.; Zhang, P.; Armstrong, R. N.; Gilliland, G. L. The threedimensional structure of a glutathione S-transferase from the mu gene class. Structural analysis of the binary complex of isoenzyme 3-3 and glutathione at 2.2-Å resolution. Biochemistry 1992, 31, 10169–10184. (78) Jakobsson, P. J.; Morgenstern, R.; Mancini, J.; Ford-Hunchinton, A.; Persson, B. Common structural features of MAPEG: a widespread superfamily of membrane associated proteins with highly divergent functions in eicosanoid and glutathione metabolism. Protein Sci. 1999, 8, 689–692.

JP8007688