Published on Web 09/01/2004
Computational Titration Analysis of a Multiprotic HIV-1 Protease-Ligand Complex Francesca Spyrakis,† Micaela Fornabaio,† Pietro Cozzini,‡ Andrea Mozzarelli,*,† Donald J. Abraham,§ and Glen E. Kellogg*,§ Department of Biochemistry and Molecular Biology, UniVersity of Parma, 43100, Italy, Department of General and Inorganic Chemistry, UniVersity of Parma, 43100, Italy, and Department of Medicinal Chemistry and Institute for Structural Biology & Drug DiscoVery, Virginia Commonwealth UniVersity, Richmond, Virginia, 23298-0540 Received June 10, 2004; E-mail:
[email protected];
[email protected] It has generally proven difficult to correlate solution binding data and crystallographic structural data with an in-depth understanding of the atomic-level chemical details of binding. This is especially true in cases where the binding site and/or the ligand are multiprotic, i.e., possessing multiple ionizable groups. While ionization states are governed by the pH of the media and the pKa’s of the ionizable moieties, localized effects can significantly impact ionization, and each ionizable group is not independent of all others in the complex. The desired understanding is further obscured by experimental constraints: (1) there are relatively few cases where binding as a function of pH has been reported for structurally well-characterized protein-ligand complexes, (2) the pH at which solution binding is measured may not be the optimum pH for binding, and (3) the pH of crystallization may not be the same as that for solution measurements. Thus, the linkage between solution data and crystal structure may be fairly tenuous, and predictions, e.g., in virtual screening, may be of questionable value. Using our free energy scoring tool HINT, which has been validated for a number of biomolecular systems where protonation states were not explicitly considered,1 we recently explored ligands bound to influenza virus neuraminidase and proposed a new calculation paradigm that we termed “computational titration”.2 We were able to model a multiplicity of protonation states for the complexes and propose families of “isoprotic” structures. Because this new paradigm requires building and optimization of tens to thousands of models with specific protonation for each proteinligand complex, some amount of manual filtering was necessary and models that would obviously be of high energy were excluded. A more accurate method of energy evaluation would also include a statistical mechanics evaluation of all models. With these considerations we have implemented a new computational tool that in silico builds and energy scores molecular models of all protonation states for a ligand-bound active site, including the energetic effects of bridging water molecules.3 To evaluate this tool, we have modeled the pH-dependent inhibition of HIV-1 protease by Glu-Asp-Leu, as recorded experimentally by Louis et al. between ca. pH 3.0 and 5.0.4 A 2.0 Å resolution structure of the same complex was obtained on crystals grown at pH 4.2.4 Figure 1 illustrates the environment of the ligand at its binding site nominally at “pH 7”. There are four Asp (25R, 29R, 30R, 25β) residues interacting with the peptide, and the peptide itself has three carboxylate groups and an amine from the side chains and termini. Thus, there are eight potentially ionizable moieties at the binding site, although one, the glutamate carboxylate of the peptide, is positioned largely out of the cavity and is directed into the solvent. Two potentially significant3 water molecules (301 and 313) are also present at the active site. The close interactions between peptide † ‡ §
Table 1. Computational Titration Results for Glu-Asp-Leu in HIV-1 model charge
No. modelsa
min scoreb,c
max scoreb,c
avg scoreb,c
Boltzmannb,c,d
-7 -6 -5 -4 -3 -2 -1 0 1
1 15 98 364 840 1232 1120 576 128
6.04 1.11 -2.52 -5.80 -8.16 -8.47 -8.57 -8.37 -7.94
6.04 6.53 6.38 5.51 4.77 3.51 0.80 -1.23 -3.24
6.04 3.91 2.00 0.30 -1.19 -2.50 -3.62 -4.58 -5.39
6.04 2.35 -0.44 -2.84 -4.61 -5.35 -5.79 -6.07 -6.23
a Each carboxylate has three possible states: deprotonated (1) and protonated at either carboxy oxygen (2). b Scores converted to free energy with eq 2.3 c In kilocalories per mole. d Boltzmann-weighted free energy of ligand binding.
and protein (Figure 1) suggest that at least some of the protein and/or ligand carboxylates must be protonated for viable ligand binding. In particular,3-5 the Asp25R/Asp25β region is normally protonated in some manner. The molecular models were prepared using the Sybyl modeling program as described previously.1c,2,3 First, the coordinates for the complex were retrieved from the PDB (accession code 1a30), protons were added to all atoms (protein, peptide, and water), and the resulting model was subjected to an energy minimization that allowed only the (nonexperimentally determined) hydrogens to relax. Finally, the protein, peptide, and water were separated for computational titration. Analysis of this molecular model indicated 4374 unique protonation models. These models (Table 1) range from the most basic (all eight sites deprotonated, charge -7) to the most acidic (all eight sites protonated, charge +1). The computational titration procedure built molecular models for each case, optimized the positions of polar hydrogens (-OH, -NH2, and -SH), and reoriented the two water molecules as necessary with exhaustive torsional rotation algorithms6 to improve protein-ligand-water hydrogen bonding. The reported HINT1 score (HTOTAL) for each protonation model is:
HTOTAL ) Hprotein-ligand + Hligand-water
(1)
Biochemistry and Molecular Biology, University of Parma. General and Inorganic Chemistry, University of Parma. Virginia Commonwealth University.
where the two terms are from protein-ligand and ligand-water interaction scores.3 It is important to note that all of the 4374 models differ only in the presence or placement of protons. Thus, any of them would still fit within the original electron density envelopes of the crystallographic data. We might term these models to be “isocrystallographic”. This is not to say, however, that the converse is necessarily true: some of these proton models likely correspond to crystallographically distinct structures, but analysis is beyond the scope of the present work. However, we would expect to see an indication of potential structural rearrangement manifested with highly unfavorable free energy binding scores.
J. AM. CHEM. SOC. 2004, 126, 11764-11765
10.1021/ja0465754 CCC: $27.50 © 2004 American Chemical Society
11764
9
COMMUNICATIONS
Figure 1. Stereoview of Glu-Asp-Leu (rendered ball and stick) bound in HIV-1 protease from 1a30.4 Two (301 and 313)3 bridging water molecules (rendered stick) are shown. Model is displayed at nominal pH 7, where all carboxylates are ionized and the primary amine is protonated. Distances (Å) between atoms for principal hydrogen bonds and unphysical interactions involving carboxylates, carbonyls, and/or water molecules are indicated.
pH-dependent binding.4 The calculated and experimental titration curves were aligned (Figure 2) by scaling and shifting the pH axis of the experimental binding data to this presumed parallel relationship. Other alignment criteria and/or assumptions are possible. Here the difference between the measured and predicted ∆Gbinding is about 0.6 kcal mol-1 throughout the pH range. The crystal model thus appears to correspond to a site charge of -3 to -2. Among the 2072 models of this site charge range, many of the more energetically reasonable ones have protonated two sites of Asp25R, Asp25β, or the C-terminal carboxylate of the peptide, either Asp30R or the Asp carboxylate of the peptide, and the N-terminal amine of the peptide. Protonation of the peptide Glu carboxylate is largely irrelevant as it interacts mostly with solvent, and may be the fifth added proton for -2 charge models. While this calculation should be applicable to nearly any protein-ligand complex with structural data, we are quite interested in further calibrations with pH-dependent binding data. Acknowledgment. We thank the Italian Ministry of Instruction, University and Research, and Virginia Commonwealth University for support. Figure 2. Composite of calculated and experimental titration results for Glu-Asp-Leu bound in HIV-1 protease. Louis et al.4 experimental data have been scaled and shifted on the pH/site charge axis to match Boltzmann data. The arrow indicates pH of crystallization for 1a30.4
Figure 2 shows the results of the titration analysis (also Table 1). The scores for each site/model charge vary over a wide range, indicating both favorable and (highly) unfavorable protonation models within each isoprotic set. For each, an arithmetic and Boltzmann average of the scores converted to free energy with3
∆Gbinding ) -0.0016HTOTAL - 4.873
(2)
is shown. Equation 2, derived from a set of 23 HIV-1 proteaseligand complexes (not including 1a30), has a standard error of (0.84 kcal mol-1.3 The Boltzmann curve shows a particularly interesting trend: an initially rapid decrease in free energy for the first few protons added to the system as the highly negative active site is neutralized, followed by a leveling off as the site appears to compensate by adopting low-energy torsional conformations for the resulting acids. This behavior parallels nicely the experimental
Supporting Information Available: Detailed description and HINT scores of protonation-specific molecular models. This material is available free of charge via the Internet at http://pubs.acs.org. References (1) (a) Wireko, F. C.; Kellogg, G. E.; Abraham, D. J., J. Med. Chem. 1991, 34, 758-767. (b) Kellogg, G. E.; Abraham, D. J. Eur. J. Med. Chem. 2000, 35, 651-661. (c) Cozzini, P.; Fornabaio, M.; Marabotti, A.; Abraham, D. J.; Kellogg, G. E.; Mozzarelli, A. J. Med. Chem. 2002, 45, 2469-2483. (2) Fornabaio, M.; Cozzini, P.; Mozzarelli, A.; Abraham, D. J.; Kellogg, G. E. J. Med. Chem. 2003, 46, 4487-4500. (3) Fornabaio, M.; Spyrakis, F.; Mozzarelli, A.; Cozzini, P.; Abraham, D. J.; Kellogg, G. E. J. Med. Chem. 2004, 47, 4507-4516. (4) Louis, J. M.; Dyda, F.; Nashed, N. T.; Kimmel, A. R.; Davies, D. R. Biochemistry 1998, 37, 2105-2110. (5) (a) Smith, R.; Brerenton, I. M.; Chai, R. Y.; Kent, S. B. H., Nat. Struct. Biol. 1996, 3, 946-950. (b) Wang, Y. X.; Freedberg, D. I.; Yamazaki, T.; Wingfield, P. T.; Stahl, S. J.; Kaufman, J. D.; Kiso, Y.; Torchia, D. A. Biochemistry 1996, 35, 9945-9950. (6) Kellogg, G. E.; Chen, D. L. Chem. BiodiVersity 2004, 1, 98-105.
JA0465754 J. AM. CHEM. SOC.
9
VOL. 126, NO. 38, 2004 11765
