880
J. Phys. Chem. B 2001, 105, 880-888
Electrostatic Complementarity at Ligand Binding Sites: Application to Chorismate Mutase Erik Kangas†,‡ and Bruce Tidor†,* Departments of Chemistry and Physics, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139-4307 ReceiVed: September 23, 2000
Charge optimization methods facilitate examination and, potentially, improvement of electrostatic interactions between binding partners. Here charge optimization was applied to the chorismate mutase from Bacillus subtilis binding an endo-oxabicyclic transition-state analogue. Electrostatically optimized templates based on calculations using the X-ray crystal structure were used to define regions of the transition-state analogue whose electrostatic properties are sub-optimal for binding. Variants of the analogue that could exhibit improved electrostatic affinity for the enzyme were considered that more closely mimicked the optimal charge distributions. Results indicate that the transition-state analogue is remarkably complementary to the enzyme active site in terms of electrostatics throughout much of the binding site. Particularly good electrostatic complementarity is exhibited for most of the groups on the analogue that make hydrogen bonds with the enzyme. While some small potential opportunities for improvement exist throughout the ligand, one significant under-compensated interaction stands out. The C10 carboxylate pays substantial desolvation penalty upon binding but does not recover strong hydrogen-bonded interactions with the enzyme in the available crystal structure. Calculations suggest that replacement with a virtually isosteric nitro group may improve the electrostatic contribution to the binding free energy by 2-3 kcal/mol. The principal mechanism is a decrease in the desolvation penalty of the ligand with a smaller loss in ligand-enzyme interactions. This work shows that charge optimization techniques are capable of identifying chemically reasonable substitutions to known ligands that produce substantial improvements in computed binding affinity through electrostatic enhancement. Comparison of the results across twelve independently refined active sites confirms that the approach is not overly sensitive to subtle structural variation.
I. Introduction Molecular association in aqueous solution is essentially an exchange reaction in which interactions with solvent in the unbound state are exchanged for intermolecular interactions between partners in the bound complex. The nature of this exchange is such that progressively increasing ligand charge or polarity to improve bound-state interactions will eventually diminish binding free energy due to an overwhelming increase in the desolvation penalty. Likewise, under-charging ligands to reduce their desolvation penalty will at some point lead to insufficient interactions with a polar or charged binding site. Using continuum models it has been shown that there is an intermediate ligand charge distribution for which the electrostatic contribution to the binding free energy is as favorable as possible for a given geometry, and methodology for computing chargeoptimized ligand templates has been developed.1,2 A number of properties of optimized charge distributions have been investigated. For example, in several calculations the resulting optimized electrostatic component of the binding free energy has been favorable,2-4 in contrast to the role of electrostatics in many natural complexes, where the net effect of polar and charged moieties is frequently unfavorable due to under-compensation of desolvation effects.5-9 In fact, it has been * Author to whom correspondence should be addressed at Department of Chemistry, Room 6-135, Massachusetts Institute of Technology, Cambridge, MA 02139-4307. Phone: (617) 253-7258. Fax: (617) 2521816. E-mail:
[email protected]. † Department of Chemistry. ‡ Department of Physics.
proven that the electrostatic binding free energy for optimized complexes is favorable for many situations of biophysical interest.10 One implication of these studies is that there may be a reasonably large margin for electrostatic improvement in many binding reactions. An issue that has yet to be addressed, however, is whether charge optimization can suggest chemically appropriate modifications to known ligands to improve binding affinity and whether the expected enhancements are significant. In this work we examine the electrostatic binding properties of an endo-oxabicyclic transition-state analogue 411,12 to Bacillus subtilis chorismate mutase (BsCM; Figure 1). To the extent that the transition-state analogue is sub-optimal, the calculations may suggest variants with improved binding electrostatics. Moreover, the results provide a plausible estimate of what charge distribution the enzyme has been designed to bind best. In this sense, it is of interest to consider the results in light of available mechanistic studies. Chorismate mutase represents a rare example of an enzymecatalyzed pericyclic rearrangement. The conversion of chorismate 1 to prephenate 3 (Figure 2), formally a Claisen rearrangement, is a key step in the synthesis of the aromatic amino acids tyrosine and phenylalanine via the shikimate pathway in microorganisms and plants, but not in humans or other mammals.13 The enzyme is of interest both for an understanding of its catalytic rate enhancement and as a target for the discovery of antibiotics, antifungals, and herbicides. An endo-oxabicyclic analog 411,12 based on the presumed transition state 2 binds with micromolar affinity to the natural enzyme from a number of organisms and has been used as a hapten to generate catalytic
10.1021/jp003449n CCC: $20.00 © 2001 American Chemical Society Published on Web 01/06/2001
Electrostatic Complementarity at Ligand Binding Sites
J. Phys. Chem. B, Vol. 105, No. 4, 2001 881
Figure 2. The chorismate mutase enzyme catalyzes the conversion of substrate chorismate 1 to product prephenate 3. The endo-oxabicyclic transition-state analogue 4 resembles the presumed transition state 2.
Figure 1. Structure of B. subtilis chorismate mutase. The enzyme is a trimer with pseudo-3-fold symmetry. Three active sites, one at each dimeric interface, are present in each trimer. (a) Each of the trimer chains is shown in a different color, and the endo-oxabicyclic transitionstate analogue 4 is shown in one of the three pseudo-symmetric sites. (b) Enlargement of 4 in an active site with key side chains indicated.
antibodies, which are currently less proficient than the enzymes.14-16 Much is known about the enzyme-catalyzed mechanism.17-29 Because the uncatalyzed reaction proceeds readily in solution, it has been possible to compare enzymic to uncatalyzed processes. Both have been shown to proceed via a chairlike transition state.18,19 In solution the substrate exists as an equilibrium between pseudodiaxial and pseudodiequitorial conformations, with the less stable pseudodiaxial being in the proper conformation for reaction to proceed.21 Kinetic isotope evidence suggests that the transition state for nonenzymic and enzymic reactions is highly polarized, in which bond breaking precedes bond formation.17,27,29 Wiest and Houk used ab initio quantum mechanical calculations to examine the conversion of chorismate to prephenate.23,24 The calculated transition-state structure was similar to the endo-oxabicyclic transition-state analogue in many ways and different in others. The two carboxylate groups of
both structures occupied similar positions, the cyclohexene ring of the analogue was not as flat as the corresponding region of the transition state, and the overall charge pattern appeared similar with a few exceptions.23 Examination of the roles of the carboxylate groups showed that, while each lowered the activation barrier somewhat alone, the pair together increased the barrier due to enhanced repulsion between them in the transition state. Using model systems based on the BsCM active site, hydrogen-bonding interactions that lowered the activation barrier were identified, including Arg7 and Arg90.24 Hybrid quantum mechanical and molecular mechanical (QM/MM) simulations of the B. subtilis reaction carried out by Lyne, Mulholland, and Richards suggested the enzyme binds a somewhat strained form of the pseudodiaxial conformation of the substrate and is most complementary to the transition state, with particularly large contributions arising from Glu78 and Arg90.26 Substantial effort has been made to synthesize enzyme inhibitors, many based on substrate, presumed transition state, and product. Early studies by Andrews and co-workers produced inhibitors designed as transition-state analogues, including a number of substituted adamantanes.30 Bartlett and co-workers synthesized a somewhat higher-affinity endo-oxabicyclic transition-state analogue 4, which binds to many chorismate mutases and has been used in a number of studies, including crystal structures.11,12 Further synthetic studies have been undertaken to probe the relative roles of individual functional groups in the endo-oxabicyclic analogue, to attempt to isolate more potent inhibitors, and to search for inhibitors with specificity for a particular chorismate mutase.31-35 This system is especially amenable to computational analysis. The B. subtilis chorismate mutase is monofunctional. Its structure has been solved with the endo-oxabicyclic transitionstate analogue 4 bound and also with the product prephenate 3 bound to 2.2 Å and free enzyme initially to 1.9 Å.36,37 After this work was completed a 1.3-Å crystal structure of the unliganded state became available.38 There seems to be no significant conformational change of the enzyme upon binding.36 The active site and the ligands are highly charged, and mutational studies have shown that the electrostatic properties of the active site residues are critical for catalysis.39-41 The kinetic parameters for the enzyme are Km ) 100 µM and kcat ) 50 s-1 (ref 42). The endo-oxabicyclic transition-state analogue 4 and the product prephenate 3 competitively inhibit the B. subtilis enzyme with Ki’s of 3 µM and 70 µM, respectively.22 The measured kcat/Km is essentially invariant with pH between pH 5 and pH 9, suggesting that the titration states of the important residues are essentially fixed near neutral pH.22 The transition-
882 J. Phys. Chem. B, Vol. 105, No. 4, 2001 state analogue 4 itself is fairly rigid, with few degrees of freedom that may be affected by binding. Given the geometry of the bound complex, the results show that the endo-oxabicyclic analogue 4 matches the electrostatic optima well in most regions, particularly those making hydrogenbonded contacts with the enzyme. There are, however, a few regions of especially large discrepancy, which provide a basis for trial improvements. The suggested improvements include reducing the charge of the C10 carboxylate group and altering the polarity somewhat through substitution of electronegative atoms; computed improvements of up to 2-3 kcal/mol appear available. A unique feature of the crystal structure used here is the presence of 12 separately refined copies of the active site in the asymmetric unit. While not at all a requirement for application of the methodology, multiple structures provided a useful test of robustness. The results also show a very strong similarity among the electrostatic optima computed for the 12 crystallographically independent active sites, which indicates that the algorithm is not overly sensitive to subtle differences in the X-ray coordinates. The remainder of the paper is organized as follows: Section II describes the methods used for these calculations. Results are presented in Section III and discussed in Section IV. II. Methods Structure Preparation. The 2.2-Å X-ray crystal structure of the endo-oxabicyclic transition-state analogue 4 bound to BsCM was the starting point for all calculations (PDB identifier 2CHT).36,37 The enzyme is a pseudo-3-fold symmetric homotrimer with catalytic active site clefts at each of the three subunit-subunit interfaces. The asymmetric unit in the crystal contains four trimeric complexes, which were refined separately in the structure determination. As a result there are twelve versions of the bound-state active site structure. The twelve subunits are labeled A through L, and the twelve active sites are labeled with their constituent subunit pair (e.g., AB, BC, and AC). Except where explicitly noted, all crystallographic water molecules were deleted from the structure, as no water molecules could be identified as consistently participating in ligand binding across all active sites at this resolution. Each trimer was separately extracted from the crystal structure and treated in isolation. Calculations reported were carried out for binding ligand to each active site with the other two active sites empty; trial calculations indicated negligible interaction among the sites. The CHARMM molecular modeling package43 was used with the CHARMM22 all-atom parameter set44 to prepare the structure for computations. Because there was no electron density for the first residue and the last 8-12 residues of each subunit, the N- and C-termini of each chain were modified with acetyl and N-methyl blocking groups, respectively. Effective molecular mechanics parameters for 4 were constructed using the quantum mechanically derived partial atomic charges described below, combined with appropriate nonbond parameters from the CHARMM22 parameter set. These were used to build all hydrogen atoms and any missing heavy atoms into the crystal structure. There are three histidine residues per subunit; all were treated as being neutral. Two are surface exposed and not in a position to make additional protein hydrogen bonds, and the third is partially buried and positioned to accept a hydrogen bond. The electrostatic calculations which followed were carried out using a polar-hydrogen atom representation for the protein; operationally, nonpolar hydrogen atoms were removed from the structure file as the final step of preparation.
Kangas and Tidor General Electrostatic Calculations. Electrostatic binding free energies and potentials were computed in the continuum electrostatic approximation45-47 with a solute dielectric constant of 4, a solvent dielectric constant of 80, a bulk ionic strength of 0.145 M, and a 2.0-Å ion-excluding Stern layer.48,49 The finite-difference linearized Poisson-Boltzmann equation was solved using a modified version of the DELPHI computerprogram package.50-52 The results of 10 translations on a 0.257-Å grid centered on the ligand were averaged in all cases. The 0.257-Å grid spacing was obtained using a three-step focusing procedure on a 129 × 129 × 129 grid, where in the last stage only those partial atomic charges within a cubical region 33.2 Å on a side, centered on the ligand, were treated on the grid. All longer-range interactions were treated with a coarser calculation on a 0.514-Å grid. Other details of the computations were as in our previous work.2,53 For the continuum electrostatic calculations, the PARSE parameter set (partial atomic charges and radii) was employed.54 For 4 and all variants, quantum mechanically derived charges for all atoms were combined with PARSE-like radii (1.4 Å for oxygen, 1.0 Å for hydrogen, 1.7 Å for carbon, and 1.5 Å for nitrogen). The partial atomic charges of the ligands were calculated by first geometry optimizing in the gas phase using the B3LYP/6-31+G* basis set and the Jaguar55 computer program package. Next, molecular electrostatic potentials were computed from a single-point calculation using the Gaussian computer program package and the same method and basis set.56 The partial atomic charges at the atom centers of each ligand were obtained using a two-stage RESP fit to these potentials, constraining the oxygen partial atomic charge of each carboxylate and nitro group to be symmetric.57,58 Charge Optimization. The charge distribution for the analog was optimized using a variable point charge at each of its 26 atomic centers. The bound-state conformations were from the X-ray crystal structure, and the unbound-state conformations of the enzyme and ligand were taken to be the same as the bound state. This procedure is reasonable given the similarity of the protein structures between the bound and unbound states.36 The theory is described in detail in other work.1,2 Here we briefly state what numerical procedures were used to implement that theory. The electrostatic contribution to the binding free energy can be represented as
∆G°es(r;1) ) Qr†RQr + Qr†CQ1 + Q1†LQ1
(1)
where Ql and Qr are the ligand and receptor charge distributions, respectively, and L, R, and C are matrixes that represent the ligand desolvation, receptor desolvation, and ligand-receptor interaction, respectively. The ligand desolvation matrixes, L, were numerically obtained by individually charging each ligand atom center to 1.0e in the bound and unbound states and then taking the difference in the resulting potentials at ligand atom centers. Using the principle of superposition, the potential differences were assembled into the ligand desolvation matrix.2 By reciprocity, this matrix must be symmetric; because the numerical calculations produce small deviations, the matrix was symmetrized by averaging the corresponding off-diagonal elements. The bound-state ligand-receptor interaction vector, C†Qr, was computed by charging up the entire enzyme to its actual charge distribution in the bound state and calculating the potentials at the ligand atom centers. The optimized electrostatic charge distribution for the ligand was then computed from1,2
1 -1 † Qopt 1 ) - L C Qr 2
(2)
Electrostatic Complementarity at Ligand Binding Sites
Figure 3. The electrostatic binding free energy of 4 (unfilled squares), of FO (full optima; unfilled circles), the differences between these values (triangles), and of AO (average optimum; diamonds), for each of the twelve active sites.
The ligand desolvation matrixes, L, were very well behaved: cutoffs based upon singular value decomposition2,59,60 did not need to be applied. An estimate of the error in each eigenvalue, δei, was computed by propagating the standard error in the mean of the matrix elements, Lij, arising from 10 grid translations, through an inner product with each eigenvector, ei, which gives the eigenvalue, e.g., ei ) etiLei. The relative estimated error was small enough (|δei/ei| < 0.20) that no eigenvectors needed to be discarded in the matrix inversion. The LOQO computer program61,62 was employed to obtain solutions to eq 2. The partial atomic charges at the individual atom centers were constrained to be less than or equal to 0.85e in magnitude; however, this constraint was seldom needed. Where noted, other constraints were applied, such as that the total charge on the ligand summed to a given integer. III. Results The electrostatic binding free energy of the endo-oxabicyclic transition-state analogue 4 was computed separately in each of the twelve enzyme active sites of the X-ray crystal structure of B. subtilis chorismate mutase. These free energies are given in Figure 3. Most of the electrostatic binding free energies were unfavorable, ranging from -2.2 to 7.7 kcal/ mol, with an average (( standard deviation) of 2.1 ( 2.7 kcal/ mol. It is interesting to note that in 25% of the X-ray complexes, electrostatic interactions faVored binding. This result contrasts with calculations on a number of complexes, for which the net contribution of electrostatics was unfavorable due to incomplete compensation of large, unfavorable desolvation terms,5-9 and implies that electrostatic interactions are already quite complementary in this highly charged complex. The relatively wide range of computed binding contributions from the twelve individual active sites indicates that substantial energetic differences can result from rather small differences in interpretation of protein structural data. These differences do not appear to result from any limitation in the calculations but rather from geometric differences in the refinements at this resolution. Optimization of the Transition-State Analogue. Optimized charge distributions were computed for each of the twelve enzyme active site-transition-state analogue structures, resulting
J. Phys. Chem. B, Vol. 105, No. 4, 2001 883 in twelve optimized ligand charge distributions and electrostatic binding free energies. These binding free energies provided a lower limit to the electrostatic binding free energy attainable through altering the ligand charge distribution at the atom centers in the specified geometries. As expected, all optimized binding free energies (Figure 3) were much improved over 4 (by 6.2 ( 1.0 kcal/mol). The actual value of the optimized electrostatic binding free energy also varied widely, from -7.9 to 2.3 kcal/ mol. In 11 of 12 cases, it was favorable. It is important to note that, while the electrostatic binding free energy varied widely from site to site, indicative of some dependence on conformation, improVements in the electrostatic binding free energies through optimization were much more consistent, indicative of much less conformational sensitivity. This important result indicates structural uncertainties that may limit calculations of absolute binding affinity need not interfere with calculations of relative binding affinity. An average optimized charge distribution was constructed by averaging over the 12 individual optima for each atom center (denoted AO for average optimum). The average optimum also had a favorable electrostatic binding free energy to all but one active site. As seen in Figure 3, it bound only 1.0 ( 0.3 kcal/ mol worse than the individual optima, still 5.2 ( 0.8 kcal/mol better than 4. In subsequent analysis we use this average optimum for convenience, but similar results were obtained using the actual optimum for each site individually (denoted FO for full optimum). This treatment is justified here because the average optimum bound to each active site practically as well as each individual optimum. The total charge of the average optimum was -1.42e ( 0.17, approximately half a charge unit less negative than 4. In fact, when the total ligand charge was constrained to an integer, 8 of the 12 active sites produced an optimized ligand with a charge of -1e. Constraining the total charge to be -3.0e, -2.0e, -1.0e, 0.0e, or +1.0e, resulted in improvements over 4 in the constrained optimal binding free energy of 0.5 ( 1.6, 5.5 ( 1.1, 5.8 ( 1.0, 2.0 ( 1.9, and -6.8 ( 4.0 kcal/mol, respectively. Total charges for the ligand of -1.0e or -2.0e were much more favorable than other values; furthermore, -1.0e was preferred, based on the largest improvement, and also had the smallest variance by a small margin. Thus, the results suggested the endooxabicyclic transition-state analogue 4 was over-charged for binding to these active sites. Figure 4 shows the difference between the quantum mechanically derived charge distribution for 4 and the computed AO; Figure 5 graphically displays the charges of 4 (in parentheses) and the average optimum and also shows active site residues making hydrogen bonds with the ligand. It can be seen that 65% of the charges for 4 were within 0.25e of those for the average optimum. This indicates that, for the most part, 4 is close to optimum. In regions of the ligand that are hydrogenbonded to protein, most of the partial atomic charges for 4 were especially similar to those of the average optimum. In particular, the C11 carboxylate, the hydroxyl group at position 4, and the ether oxygen (O7), which form hydrogen bonds with the enzyme, have charges quite close to the computed optimum. The area of highest dissimilarity between the charge distributions is the C10 carboxylate region, which does not directly contact the enzyme in any of the twelve active sites in the structure. The computed optimal charge distribution differs from 4 in that much of the negative charge from the C10 carboxylate oxygens is moved onto the more buried carbon atoms. This area may provide an opportunity for engineering an improved binding ligand. Furthermore, the results indicate that the method is fairly
884 J. Phys. Chem. B, Vol. 105, No. 4, 2001
Figure 4. Deviation between optimum and actual partial atomic charges. The deviation in the charge distribution of 4 from that of AO (average optimum) is indicated by the points; the bars represent one standard deviation of the distribution over the twelve active sites. A positive deviation indicates that the average optimized charge is more positive than that of 4.
Figure 5. The transition-state analogue structure together with the charge distributions of 4 and AO. The original analogue 4 charges are shown in parentheses; the average optimized charges are shown without parenthesis. Shown in blue are five of the active site residues that make hydrogen bonds with 4 in many of the twelve active sites. Consensus hydrogen bonds are shown in orange.
robust. For instance, the standard deviation of the optimized charges computed over the twelve active sites is quite small for most atom positions, including the C10 carboxylate oxygens. The exceptions, positions with larger distributions in optimized charge, are frequently less exposed regions of the ligand and correspond to shallow dimensions of the free energy paraboloid. That is, there is a relatively small binding penalty for deviations from optimal charge at these locations.63 Other areas that may also be of interest because they were computed to be sub-optimal are the methylene groups at positions 6 and 9. These are not making interactions with active site residues but may have other indirect effects. It may also be useful to examine other carbon-hydrogen positions for which optimization predicts that the hydrogen should be more negatively charged, as this type of charge perturbation may be introduced with a fluorine substitution.64,65
Kangas and Tidor Modifications to the Transition-State Analog. Differences between the partial atomic charges computed quantum-mechanically for the endo-oxabicyclic transition-state analogue 4 and the optimized charge distribution were used to guide hypotheses for tighter-binding ligands. We looked for modifications that increased the similarity between ligand partial atomic charges and the optimized charge distribution. First, we chose a region of 4 with sub-optimal charge distribution relative to the calculation. Then, we re-optimized only those atom centers in the region of interest, fixing all other charges to their original analog values and constraining the total charge of the ligand to an integer. This approach produces the best local electrostatic modification that could be made to this region, which is generally similar to the region in the full optimization. The reoptimization provides important information by giving an estimate of the maximum electrostatic free energy gain that can be made through local modification; if the gain is small, the region might not be worth modifying. We then sought local modifications that would approximate the locally optimized charge distribution. Examination of Figure 5 revealed that the region near the C10 carboxylate (C1, C10, and the two associated oxygen atoms) was furthest from the full optimum. The partial atomic charges were (0.16e, 0.83e, -0.84e, -0.84e), whereas the charges in the average optimum were (-0.70e, 0.33e, -0.16e, -0.25e). Local optimization of these four atoms, with the total ligand charge constrained to be an integer, resulted in a total ligand charge of -1.0e in 7 of 12 active sites. Averaged across all active sites, local optimization with the total charge constrained to -1.0e (-2.0e) improved the electrostatic binding free energy by 2.3 ( 0.9 (1.7 ( 0.6) kcal/mol. Thus, a local change that reduced the total negative charge from -2e to -1e could enhance binding. The locally optimized charge distribution for -1e was (0.63e, -0.37e, 0.05e, -0.03e). The oxygen charges were reduced to almost zero, C10 became somewhat negative, and C1 became more positive. The dual constraints of total charge of -1e and fixing the charges at the other centers left relatively little freedom for these charges. It should be noted that the binding free energy is most sensitive to the charge magnitudes of the two oxygen atoms, which are reduced from those in 4, similar to the full optimizations but somewhat more so. This, again, is consistent with the picture that the oxygen atoms are overcharged in 4. Modifications to 4 consistent with this scheme include removing the carboxylate and replacing it with a hydrogen atom or methyl group. Alternatively, to maintain some of the polarity of the carboxylate, replacement with a fluorine atom or trifluoromethyl group might be better. Because of our focus on modifying the charge distribution separately from the shape wherever possible, we are most enthusiastic about a carboxylate to nitro substitution 5, with carbon atom C10 replaced with a nitrogen. This substitution changes the total charge of the ligand to -1e and also reduces the charge magnitude on the oxygen atoms, resulting in partial atomic charges of (0.17e, 0.68e, -0.46e, -0.46e). We calculated that 5 should bind 2.2 ( 1.1 kcal/mol better than 4. The charge differences between 4 and 5 extended beyond the local neighborhood of the modification; the charge on atom C4 changed by 0.23e and the charge on C3 changed by about 0.1e toward the full optimum. Allowing only the four local charges to change from those of 4 in the RESP fit to the quantum mechanical potential produced an improvement of 1.4 ( 1.0 kcal/mol. Table 1 describes other variations of 4. The best electrostatic binders 5-7 and 16-18 all included variations of the C10
Electrostatic Complementarity at Ligand Binding Sites TABLE 1: Computed Improvement in Electrostatic Binding Free Energy of Variants Compared to 4a variants of 4
∆∆Ges (kcal/mol) 6.2 ( 1.0 5.2 ( 0.8 0.0 ( 0.0
FO AO 4
fully optimized average optimized original analog
5 6 7 8 9 10 11 12 13 14 15
Single Modifications C10 carboxylate to nitro C10 carboxylate to aldehyde (front) C10 carboxylate to aldehyde (back) pro-R H9 to fluorine pro-S H9 to fluorine H3 to fluorine H2 to fluorine pro-S H6 to fluorine C11 carboxylate to nitro C11 carboxylate to aldehyde (back) pro-R H6 to fluorine
2.2 ( 1.1 2.2 ( 1.2 2.1 ( 1.0 0.4 ( 0.1 0.4 ( 0.1 0.5 ( 0.8 0.1 ( 0.1 0.0 ( 0.2 -3.0 ( 1.6 -4.2 ( 1.9 -1.3 ( 0.3
16 17 18 19 20 21
Multiple Modifications modifications 5 + 8 modifications 5 + 11 modifications 5 + 8 + 11 modifications 8 + 11 modifications 8 + 9 modifications 12 + 15
1.9 ( 1.2 2.1 ( 1.1 1.7 ( 1.2 0.4 ( 0.2 0.4 ( 0.3 -0.2 ( 0.3
a The values are given as the average and standard deviation over the twelve active sites. For carboxylate to aldehyde modifications, the parentheses indicate the oxygen that remains, using the view in Figure 5.
carboxylate region that reduced the ligand charge to -1.0e, either by substituting the nitro group or replacing the carboxylate with an aldehyde. There were also several fluorine variations 8- 9 that were computed to improve the electrostatic binding free energy by about 0.5 kcal/mol each. Some variations that are shown (13-15) resulted in the charges on very sensitive atom positions moving further from the local optimum than the original analog’s charges were, and such variations were computed to degrade electrostatic binding. Robustness of Results. It is important to evaluate the robustness of the calculations presented here. We have examined the electrostatic binding in 12 different X-ray crystal sites and found a significant predicted advantage in each one; thus the predictions are likely to be independent of extremely subtle bound-state conformational detail. Furthermore, most of the variations suggested would have little effect on nonelectrostatic aspects of binding, and one may expect some relaxation of these molecules, relative to the positions in the given X-ray crystal structure, which would further improve their affinity for the enzyme. However, there are other factors that could affect the relative binding free energy of these proposed ligands. Some of these are examined below. First, even though there are 127 residues in each subunit of the homotrimer B. subtilis chorismate mutase,42 the subunits in the X-ray crystal structure are truncated between residues 115 and 120, depending on the subunit, due to a lack of electron density.37 In the original uncomplexed crystal structure of this enzyme,36 all subunits are truncated at residue 115. A recently released 1.3-Å structure of the free trimer shows the C-terminal tail extended away from the structure.38 What conformation these residues adopt in the bound state is unclear. Residues 116127 are Arg, Pro, Asp, Leu, Ser, Leu, Thr, Lys, Asn, Thr, Glu, and Leu. The reported results parallel most other studies and ignore any residues not explicitly in the X-ray crystal structure, assuming them to be disordered and to not significantly affect binding. However, it is possible that these C-terminal residues may to some extent occlude the entrance to the active sites,
J. Phys. Chem. B, Vol. 105, No. 4, 2001 885 essentially causing the bound ligand to be much more buried, or they may interact with the ligand. Because there is no electron density for them, they may not make any important hydrogen bonds but simply reduce the amount of solvent near the active site. We modeled this potential effect by recalculating the results with the addition to the bound state of a 16-Å radius lowdielectric sphere centered on the ligand. This adds an additional volume somewhat larger than that of the missing residues; however, one expects the positions of these residues to be fluctuating and thus they could occupy more space than fixed residues would. This addition represents a limiting case where the C-terminal tail totally covers the bound-state active site, excluding all solvent and fully burying the ligand. The average electrostatic binding advantage of 5 decreased with the addition of this dielectric region to 1.1 ( 3.3 kcal/mol, but the standard deviation increased substantially; this represents an upper limit to the effect of the C-terminal tail. These data suggest that a disordered C-terminal tail might not strongly affect the results presented here. Nevertheless, when some of the variant ligands suggested here are synthesized, they should be assayed against full-length and truncated versions of the protein s the latter presumably corresponding more closely with the calculations. The results may also have a certain dependence on the choice of charge distribution for the ligands and enzyme. To ascertain this dependence, we compared the binding of 4 and 5 using CHARMM PARAM19 partial atomic charges and radii instead of the PARSE partial atomic charges and radii for the enzyme; the charges and radii used for the ligands were unaltered. We find that this affected the relative electrostatic advantage of 5, improving it to 4.0 ( 1.0 kcal/mol. We also examined the effect of the choice of quantum mechanical method and basis set on the ligand charge distributions. Gas-phase minimization of the ligands using the RHF/6-31G* method followed by a two-stage RESP fit to the molecular electrostatic potential produced by this level of theory gave slightly different sets of ligand charge distributions. Use of these charges for 4 and 5 (and the standard PARSE charges and radii for the enzyme), resulted in an advantage of 1.9 ( 1.1 kcal/mol for 5; about 0.4 kcal/mol worse than the advantage obtained at the B3LYP/6-31+G* level. The RHF/6-31G* basis, commonly used in charge fitting less highly charged molecules, does not include the diffuse functions that may be important in accurately describing the electrostatic potential of the divalent anionic ligands. The presence of water molecules bound to the ligand-enzyme complex may have an effect on the overall discrimination between 4 and any variants. The X-ray crystal structure of the bound transition-state analogue-enzyme complex contains 133 water molecules per enzyme molecule;36,37 however, none of these (especially those with low B-factors) make consistent interactions with the ligand across all active sites. In fact, in only 5 of 12 sites is there a bound water molecule within 5 Å of the C10 carboxylate group, and of these, only one has a B-factor less than 20. To examine the maximum effect that bound water molecules may have on the differentiation between the endo-oxabicyclic transition-state analogue 4 and its variants, all bound water molecules within 16 Å of the center of the ligand were included in the electrostatic binding free energy calculations. The resulting discrimination between 4 and 5 worsened only slightly to 2.1 ( 1.2 kcal/mol. Discarding all water molecules with B-factors larger than 20, the discrimination was still 1.9 ( 1.4. Thus, the presence of bound water molecules might slightly reduce the preference for 5 (and perhaps all variations that reduce the charge of the C10 carboxylate).
886 J. Phys. Chem. B, Vol. 105, No. 4, 2001 IV. Discussion and Conclusion Electrostatic optimization calculations carried out on an endooxabicyclic transition-state analogue 4 binding to B. subtilis chorismate mutase provide a detailed view of electrostatic compatibility. The results suggest that hydrogen-bonded polar and charged groups, including the C11 carboxylate, the ether oxygen, and the hydroxyl, are all extremely complementary to the site. A number of regions have been identified that are not as well compensated, the most prominent being the C10 carboxylate. Computations have suggested a nitro substitution at this site, 5, could enhance binding by 2-3 kcal/mol. These conclusions, while sensible in hindsight, are not immediately suggested by a visual inspection of the structure. The chief mechanism by which 5 improved the computed electrostatic binding free energy was by lowering the desolvation penalty of the ligand, rather than by improving ligand-enzyme interactions. For instance, while 5 made poorer interactions with the enzyme by 8.9 ( 2.2 kcal/mol compared to 4, it also paid 11.1 ( 2.0 kcal/mol less in desolvation. The decrease in the desolvation penalty came primarily by reducing the total charge of the ligand by 1 charge unit; because this was in a region not making direct interactions with the active site, the interaction loss was relatively small. Other modifications that similarly reduced the charge were computed to be extremely deleterious to binding. For example, diminishing the charge on the other (C11) carboxylate by changing it to a nitro group (13) or an aldehyde (14) was computed to diminish the electrostatic binding contribution by 3.0 or 4.2 kcal/mol, respectively. The results of our calculations can be compared to a number of careful structure-activity studies carried out using inhibitors to chorismate mutase. The comparison must be indirect, however, because these studies have generally not utilized the enzyme from B. subtilis. Andrews and co-workers used a series of substituted adamantanes to probe inhibition of chorismate mutase activity in the bifunctional enzyme from Escherichia coli K-12.30 Results showed the importance of a negatively charged group, as adamantane-1-carboxylic acid was about 6-fold more potent than 1-hydroxyadamane (neutral) and 200fold more potent than 1-aminoadamantane (positively charged). The role of an appropriately placed hydroxyl group (analogous to position 4 in 4) enhanced inhibition by 8-fold. Most interestingly, however, addition of a second carboxylic acid was unfavorable by between 5- and 40-fold. The former is for adamantane-1,3-diacetic acid relative to adamantane-1-acetic acid, and the latter is for adamantane-1,3-dicarboxylic acid relative to adamantane-1-carboxylic acid. These results correlate extremely well with the calculations reported here, in which the C11 carboxylate and hydroxyl appear close to optimal, yet the C10 carboxylate appears overcharged. The loss in binding experimentally was larger for carboxylic acid groups than for acetic acid groups, which could be due to a greater desolvation penalty paid by the shorter functional group. Bartlett and co-workers have also carried out an elegant structure-activity study of the bicyclic series of inhibitors that includes 4 and tested them against the bifunctional enzyme from E. coli.32 The role of the C11 carboxylate, in particular was probed. Altering the stereochemistry at C11, giving the exooxabicyclic analogue, reduced binding about 700-fold. Moreover, several methods of changing the hybridization of C8 from sp3 to sp2 were attempted. Introduction of a double bond between C8 and C9 reduced binding by about 90-fold and was thought to distort the ring away from its structure in the transition state. Attempts to make the oxabicyclic nitronate were unsuccessful due to the lability of R-nitro ethers. However, the oximinolac-
Kangas and Tidor tone was synthesized (which was a poor inhibitor), and the nitronate was introduced into the carbabicycle (where it diminished binding by about 60-fold). Moreover, the endooxabicycle 4 was 250-fold more potent than the endo-carbabicycle. These results point to the importance of the C11 carboxylate and O7 ether, consistent with the results of our calculations (albeit in the active site from a different species); however, a corresponding study of the C10 carboxylate has not been reported. Given the difference in geometry around C11 between the endo-oxabicyclic analogue 4 and the presumed transition state 2, one puzzle that remains is why the endooxabicyclic compound is difficult to improve upon experimentally by variations in the neighborhood of C11 and why the calculations find the charge distribution in this region so close to optimal. Our results are relevant to mechanistic studies of chorismate mutase. Besides suggesting a specific compound, 5, that might bind better to BsCM, the calculated optimized charge distribution represents, in some sense, what the enzyme is designed to bind tightest. This molecular complement can be compared to species present along the catalytic pathway, including the transition state. When the total charge was constrained to -2e, the average optimized ligand charge distribution was highly polarized, with a charge of -1.47e on the enolpyruvyl-derived moiety and -0.53e on the cyclohexadienyl-derived moiety. This compares with -1.15e and -0.85e for each of the two moieties, respectively, in calculations of the transition-state analogue 4; the values are -1.34e and -0.66e for a model of the gas-phase transition state computed at the B3LYP/6-31G* level of theory and fit with RESP charges.66 These results suggest that the transition state is more highly polarized than the analogue, which is consistent with experimental and theoretical studies indicating that C-C bond cleavage precedes C-O bond formation for BsCM29 as well as for the uncatalyzed reaction.17,27 Our results further indicate that the enzyme is designed to bind a molecule at least as polarized as the transition state. One difficulty with catalytic antibodies produced using the analogue could be that they do not adequately stabilize the actual transition state due to its additional polarity over the analogue. While this suggestion has been raised before,67 this is the first demonstration that an enzyme that catalyzes the reaction is more complementary to the extra polarization only achieved in the transition state. It will be important to carry out corresponding calculations on the catalytic antibody for a complete comparison. The substitution of a nitro group for the C10 carboxylate suggested here may not only produce an enhanced inhibitor of B. subtilis chorismate mutase, but it may also give a more specific inhibitor. The endo-oxabicyclic transition-state analogue 4 binds well to a number of different chorismate mutase enzymes (Ki of 3 and 2 µM for B. subtilis and E. coli, respectively).22,68 Until recently, relatively specific inhibitors were not available. Ganem and co-workers searched a combinatorial library of compounds and found that S-(-)-dinitrobiphenic acid is a selective inhibitor for E. coli chorismate mutase (Ki ) 13 µM) but is inactive against the enzymes from B. subtilis and S. cereVisiae.35 Because E. coli and Saccharomyces cereVisiae chorismate mutases make extensive hydrogen bonds with the C10 carboxylate of the analogue, the nitro substitution will might decrease affinity.69,70 This suggests the nitro derivative 5 could be selective for B. subtilis chorismate mutase over E. coli or S. cereVisiae, or both. The work described here represents a novel approach that could have broad applicability to structure-based design studies. Starting with the structure of a trial ligand bound to a molecular
Electrostatic Complementarity at Ligand Binding Sites target, calculations identified regions of good and regions of poor electrostatic complementarity. While in hindsight these regions make sense, they were not obvious originally. In particular, the computations suggest that the entire hydrogenbonded interface in the current structure is quite complementary. The implication is that modifications to this interface, which is where one would normally introduce changes, might not be the most productive. Rather, the computations identify the other face of the ligand as the most underoptimized, which initially might seem counterintuitive. One use of this methodology could be simply to identify underoptimized regions and then to focus combinatorial synthetic methods to attempt to isolate enhanced binders. A limitation of the calculations is that the optimized charge distribution is linked to the ligand shape for which it was computed. As the ligand shape is varied, the optimized charge distribution can change. Here the principal substitution suggested involves no change in shape. In other cases, it could be necessary to recompute the optimum as the geometry is varied. Fortunately, the computations are not prohibitively expensive and such re-optimization is generally possible. The current results also have important implications for the charge-optimization approach. The methodology is now sufficiently developed that this is the first example of a suggestion for a small-molecule analogue of an existing ligand that has resulted from the approach. Moreover, the computed binding free energy improvement is significant. Some studies have indicated that continuum contributions may overestimate binding improvements,71,72 and it may be possible that such effects will be present here as well, through other studies have shown much closer agreement between experiment and theory.73 The average over 12 active sites is 2.2 kcal/mol. This is not a trivial value. Although the method only requires one active site, that similar results were obtained in all 12 active sites indicates its robustness, particularly with respect to small variations is structure refinement. The additional information obtained with multiple active sites can be used to estimate a higher-order binding enhancement change, because the set of active sites can be treated as an ensemble. Modifying the ligand can be viewed as causing a shift in the population in the final ensemble toward those members that bind preferentially the perturbed over the original ligand. When this is done here, the computed binding free energy improvement for 5 increases to beyond 3 kcal/mol. Further research is required before charge optimization becomes a fully general tool for de novo structure-based design. For example, it would need to be shown that calculations carried out on unliganded active sites can lead to suggestions for tightbinding small-molecule inhibitors. Current work in that direction is being carried out. Nevertheless, the results presented here show that charge optimization is ready to be used in the important but more limited case in which a structure with bound ligand exists and one wishes to assess regions of sub-optimal complementarity to attempt further synthesis, perhaps with combinatorial methods. Acknowledgment. We thank Barry Honig for making the DELPHI computer program available, Martin Karplus for CHARMM, and Robert J. Vanderbei for LOQO. We are grateful to members of our research group, particularly Brian H. Sherman and Christopher D. Snow, as well as to Donald Hilvert, Stephen L. Buchwald, David W. Christianson, Rick L. Danheiser, Peter S. Kim, and JoAnne Stubbe for helpful discussions. This work was partially supported by the National Institutes of Health (GM55758 and GM56552). E.K. was supported as a fellow by the National Science Foundation and by IBM.
J. Phys. Chem. B, Vol. 105, No. 4, 2001 887 References and Notes (1) Lee, L.-P.; Tidor, B. J. Chem. Phys. 1997, 106, 8681. (2) Kangas, E.; Tidor, B. J. Chem. Phys. 1998, 109, 7522. (3) Chong, L. T.; Dempster, S. E.; Hendsch, Z. S.; Lee, L.-P.; Tidor, B. Protein Sci. 1998, 7, 206. (4) Lee, L.-P.; Tidor, B. Nature Struct. Biol. 2001, 8, 73. (5) Misra, V. K.; Sharp, K. A.; Friedman, R. A.; Honig, B. J. Mol. Biol. 1994, 238, 245. (6) Misra, V. K.; Hecht, J. L.; Sharp, K. A.; Friedman, R. A.; Honig, B. J. Mol. Biol. 1994, 238, 264. (7) Sharp, K. A. Biophys. Chem. 1996, 61, 37. (8) Shen, J.; Wendoloski, J. J. Comput. Chem. 1996, 17, 350. (9) Novotny, J.; Bruccoleri, R. E.; Davis, M.; Sharp, K. A. J. Mol. Biol. 1997, 268, 401. (10) Kangas, E.; Tidor, B. Phys. ReV. E 1999, 59, 5958. (11) Bartlett, P. A.; Johnson, C. R. J. Am. Chem. Soc. 1985, 107, 7792. (12) Smith, W. W.; Bartlett, P. A. J. Org. Chem. 1993, 58, 7308. (13) Haslam, E. The Shikimate Pathway; Wiley: New York, 1974. (14) Hilvert, D.; Carpenter, S. H.; Nared, K. D.; Auditor, M.-T. M. Proc. Natl. Acad. Sci. U.S.A. 1988, 85, 4953. (15) Jackson, D. Y.; Jacobs, J. W.; Sugasawara, R.; Reich, S. H.; Bartlett, P. A.; Schultz, P. G. J. Am. Chem. Soc. 1988, 110, 4841. (16) Jackson, D. Y.; Liang, M. N.; Bartlett, P. A.; Schultz, P. G. Angew. Chem., Int. Ed. Engl. 1992, 31, 182. (17) Addadi, L.; Jaffe, E. K.; Knowles, J. R. Biochemistry 1983, 22, 4494. (18) Sogo, S. G.; Widlanski, T. S.; Hoare, J. H.; Grimshaw, C. E.; Berchtold, G. A.; Knowles, J. R. J. Am. Chem. Soc. 1984, 106, 2701. (19) Copley, S. D.; Knowles, J. R. J. Am. Chem. Soc. 1985, 107, 5306. (20) Guilford, W. J.; Copley, S. D.; Knowles, J. R. J. Am. Chem. Soc. 1987, 109, 5013. (21) Copley, S. D.; Knowles, J. R. J. Am. Chem. Soc. 1987, 109, 5008. (22) Gray, J. V.; Eren, D.; Knowles, J. R. Biochemistry 1990, 29, 8872. (23) Wiest, O.; Houk, K. N. J. Org. Chem. 1994, 59, 7582. (24) Wiest, O.; Houk, K. N. J. Am. Chem. Soc. 1995, 117, 11628. (25) Davidson, M. M.; Gould, I. R.; Hillier, I. H. J. Chem. Soc., Chem. Commun. 1995, 63. (26) Lyne, P. D.; Mulholland, A. J.; Richards, W. G. J. Am. Chem. Soc. 1995, 117, 11345. (27) Ganem, B. Angew. Chem., Int. Ed. Engl. 1996, 35, 936. (28) Kast, P.; Tewari, Y. B.; Wiest, O.; Hilvert, D.; Houk, K. N.; Goldberg, R. N. J. Phys. Chem. B 1997, 101, 10976. (29) Gustin, D. J.; Mattei, P.; Kast, P.; Wiest, O.; Lee, L.; Cleland, W. W.; Hilvert, D. J. Am. Chem. Soc. 1999, 121, 1756. (30) Andrews, P. R.; Cain, E. N.; Rizzardo, E.; Smith, G. D. Biochemistry 1977, 16, 4848. (31) Clarke, T.; Stewart, J. D.; Ganem, B. Tetrahedron Lett. 1987, 28, 6253. (32) Bartlett, P. A.; Nakagawa, Y.; Johnson, C. R.; Reich, S. H.; Luis, A. J. Org. Chem. 1988, 53, 3195. (33) Clarke, T.; Stewart, J. D.; Ganem, B. Tetrahedron 1990, 46, 731. (34) Wood, H. B.; Buser, H.-P.; Ganem, B. J. Org. Chem. 1992, 57, 178. (35) Husain, A.; Galopin, C. C.; Zhang, S.; Pohnert, G.; Ganem, B. J. Am. Chem. Soc. 1999, 121, 2647. (36) Chook, Y. M.; Ke, H.; Lipscomb, W. N. Proc. Natl. Acad. Sci. U.S.A. 1993, 90, 8600. (37) Chook, Y. M.; Gray, J. V.; Ke, H.; Lipscomb, W. N. J. Mol. Biol. 1994, 240, 476. (38) Ladner, J. E.; Reddy, P.; Davis, A.; Tordova, M.; Howard, A. J.; Gilliland, G. L. Acta Crystallogr., Ser. D 2000, 56, 673. (39) Kast, P.; Hartgerink, J. D.; Asif-Ullah, M.; Hilvert, D. J. Am. Chem. Soc. 1996, 118, 3069. (40) Kast, P.; Asif-Ullah, M.; Jiang, N.; Hilvert, D. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 5043. (41) Cload, S. T.; Liu, D. R.; Pastor, R. M.; Schultz, P. G. J. Am. Chem. Soc. 1996, 118, 1787. (42) Gray, J. V.; Golinelli-Pimpaneau, B.; Knowles, J. R. Biochemistry 1990, 29, 376. (43) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187. (44) MacKerell, A. D., Jr.; Bashford, D.; Bellott, M.; Dunbrack, R. L., Jr.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; Joseph-McCarthy, D.; Kuchnir, L.; Kuczera, K.; Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom, B.; Reiher, W. E., III; Roux, B.; Schlenkrich, M.; Smith, J. C.; Stote, R.; Straub, J.; Watanabe, M.; Wiorkiewicz-Kuczera, J.; Yin, D.; Karplus, M. J. Phys. Chem. B 1998, 102, 3586. (45) Honig, B.; Nicholls, A. Science (Washington, D. C.) 1995, 268, 1144. (46) Sharp, K. A.; Honig, B. J. Phys. Chem. 1990, 94, 7684. (47) Davis, M. E.; McCammon, J. A. Chem. ReV. 1990, 90, 509.
888 J. Phys. Chem. B, Vol. 105, No. 4, 2001 (48) Bockris, J. O’M.; Reddy, A. K. N. Modern Electrochemistry; Plenum: New York, 1973. (49) Gilson, M. K.; Honig, B. H. Nature (London) 1987, 330, 84. (50) Gilson, M. K.; Sharp, K. A.; Honig, B. H. J. Comput. Chem. 1988, 9, 327. (51) Gilson, M. K.; Honig, B. Proteins: Struct., Funct., Genet. 1988, 4, 7. (52) Sharp, K. A.; Honig, B. Annu. ReV. Biophys. Biophys. Chem. 1990, 19, 301. (53) Hendsch, Z. S.; Tidor, B. Protein Sci. 1999, 8, 1381. (54) Sitkoff, D.; Sharp, K. A.; Honig, B. J. Phys. Chem. 1994, 98, 1978. (55) Jaguar 3.5; Schro¨dinger, Inc.: Portland, OR, 1998. (56) 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.; Millam, 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.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Cioslowski, J.; Ortiz, J. V.; 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.; Gonzalez, C.; Challacombe, M.; Goll, P. M. W.; Johnson, B. G.; Chen, W.; Wong, M. W.; Andres, J. L.; Head-Gordon, M.; Replogle, E. S.; Pople, J. A. Gaussian 98 (ReVision A.1); Gaussian, Inc.: Pittsburgh, PA, 1998. (57) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. J. Phys. Chem. 1993, 97, 10269. (58) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Kollman, P. A. J. Am. Chem. Soc. 1993, 115, 9620.
Kangas and Tidor (59) Press, W. H.; Flannery, B. P.; Teukolsky, S. A.; Vetterling, W. T. Numerical Recipes: The Art of Scientific Computing; Cambridge University Press: Cambridge, 1986. (60) Strang, G. Introduction to Linear Algebra; Wellesley-Cambridge Press: Wellesley, Massachusetts, 1993. (61) Vanderbei, R. J. Optim. Methods Software 1999, 12, 451. (62) Vanderbei, R. J. Optim. Methods Software 1999, 12, 485. (63) Lee, L.-P.; Tidor, B. Protein Sci. 2001, in press. (64) Gao, J. M.; Qiao, S.; Whitesides, G. M. J. Med. Chem. 1995, 38, 2292. (65) Doyon, J. B.; Jain, A. Org. Lett. 1999, 1, 183. (66) These values agree quite well with CHELPG-derived charges from ref 23 for the transition state and are somewhat overpolarized relative to values given for the analog. (67) Haynes, M. R.; Stura, E. A.; Hilvert, D.; Wilson, I. A. Science (Washington, D.C.) 1994, 263, 646. (68) Lee, A. Y.; Stewart, J. D.; Clardy, J.; Ganem, B. Chem. Biol. 1995, 2, 195. (69) Lee, A. Y.; Karplus, P. A.; Ganem, B.; Clardy, J. J. Am. Chem. Soc. 1995, 117, 3627. (70) Stra¨ter, N.; Schnappauf, G.; Braus, G.; Lipscomb, W. N. Structure 1997, 5, 1437. (71) Hendsch, Z. S.; Jonsson, T.; Sauer, R. T.; Tidor, B. Biochemistry 1996, 35, 7621. (72) Spector, S.; Wang, M.; Carp, S. A.; Robblee, J.; Hendsch, Z. S.; Fairman, R.; Tidor, B.; Raleigh, D. P. Biochemistry 2000, 39, 872. (73) Chong, L. T.; Hendsch, Z. S.; Tidor, B. In preparation.
