PROTEINS Structure, Function, and Genetics 25:300-314 (1996) Identification and Analysis of Long-Range Electrostatic Effects in Proteins by Computer Mode1ing:Aspartate Transcarbamylase Himanshu Oberoi,*Jaishree Trikha; Xiaoling Yuan,*and Norma M. Allewell Department of Biochemistry, University of Minnesota, St. Paul, Minnesota 55108 ABSTRACT While ion pairs are readily identified in crystal structures, longer range electrostatic interactions cannot be identified from the three-dimensional structure alone. These interactions are likely to be important in large, multisubunit proteins that are regulated by allosteric interactions. In this paper, we show that these interactions are readily detected by electrostatic modeling, using, as an example, unliganded Escherichia coli aspartate transcarbamylase, a widely studied allosteric enzyme with 12 subunits and a molecular weight of 310 kD. The Born, dipolar, and sitesite interaction terms of the free energy of protonation of the 810 titratable sites in the holoenzyme were calculated using the multigrid solution of the nonlinear Poisson-Boltzmann equation. Calculated titration curves are in good agreement with experimental titration curves, and the structural asymmetry observed in the crystal structure is readily apparent in the calculated free energies and pKl,z values. Most of the residues with PK,,~values that differ substantially from those of model compounds are buried in the low dielectric medium of the protein, particularly at the intersubunit interfaces. The dependence of the site-siteinteraction free energies on distance is complex, with a steep dependence at distances less than 5 A and a more shallow dependence at longer distances. Interactions over distances of 6 to 15 A require a bridging residue and are often not apparent in the structure. The network of interactions between ionizable groups extends across and between subunits and provides a potential mechanism for transmitting long-range structural effects and allosteric signals. lar distances. The transmission of the signal is often linked to changes in tertiary and quaternary structure; however, the forces that drive these structural changes and the mechanism of their propagation through the molecule are poorly understood. A role for electrostatic forces can readily be envisioned, because they are long range within the low dielectric environment of the protein and are critical to the overall structure and solvation of macromolecules. Moreover, most cellular regulatory mechanisms depend on changes in charge, for example, allosteric regulation of carbohydrate metabolism by nucleotides and phosphorylated sugars and the regulation of cell signaling pathways by protein phosphorylation. 1-3 Understanding the role of electrostatic effects in these processes requires the ability to investigate the energetics of changes in the protonation state of individual groups and charge-charge interactions in detail. Much of the information that one would like to have is not directly accessible experimentally, particularly for large macromolecular assemblies. As theoretical methods have improved and computational power has increased, computer modeling has become an increasingly attractive alternative. Approaches based on continuum models have been used to model enzyme catalysis, ligand binding, protein-protein association, and pK shift^.^-^ However, methods that are sufficiently rapid to model large macromolecular assemblies and to make predictions that can be tested by experiment have not been available until recently. We and others have recently developed the multigrid method for solving the nonlinear Poisson-Boltzmann equation (MGPB)and have shown that it produces a substantial increase in the speed and accuracy with which electrostatic interactions in 0 1996 Wiley-Liss, Inc. Key words: electrostatic modeling, PoissonBoltzmann equation, finite difference, multigrid, allosteric regulation INTRODUCTION Allosteric regulation of protein function requires that signals be transmitted over long intramolecu0 1996 WILEY-LISS, INC. Received July 7, 1995; revision accepted January 11, 1996. Address reprint requests to Norma M. Allewell, Department of Biochemistry, University of Minnesota, St. Paul, MN 55108. *Current address: Belmont Research, Inc., 84 Sherman St., Cambridge, MA 02140. 'Current address: Department of Biological Chemistry and Molecular Pharmacology, Harvard Medical School, 240 Longwood Ave., Boston, MA 62115. *Current address: Incstar, 1990 Industry Blvd., Stillwater, MN 55082. LONG-RANGE ELECTROSTATIC EFFECTS Fig. 1. Alpha carbon trace of the holoenzymeviewed down the threefold axis. The two catalytic trirners are in the center of the figure; the three regulatory dimers are on the periphery. The two catalytic trimers and the r chains in each regulatory dirner are related by the approximate non-crystallographic twofold axes. Produced with DAMPS (H. Oberoi and N. Allewell). proteins can be r n ~ d e l e d .The ~ , ~concurrent development of new approaches for calculating individual site titration curves1' and treating clusters of charges (Oberoi and Allewell, in preparation) makes it possible to calculate pK,,, values and model the titration behavior of molecules with several hundred ionizable groups in reasonable time on commonly available workstations. In this paper we describe the application of MGPB to unliganded aspartate transcarbamylase (ATCase) from Escherichiu coli. ATCase is a challenging test case since it is a large molecule with 2,778 residues and a molecular weight of 310 kD. There are 24 amino and carboxyl terminal groups and 786 titratable side chains, more than have been modeled for any other system. This widely studied allosteric enzyme is comprised of 12 polypeptide chains: six c chains, organized as two catalytic trimers, and six r chains, organized as three regulatory dimers, with approximate D, symmetry (Fig. 1; reviewed in Refs. 11-14). There are 85 ionizable residues on each c chain and 46 on each r chain. Both its substrates, L-aspartate and carbamyl phosphate, and the nucleotides that regulate its affinity for L-aspartate are polyanions, and a substantial body of experimental results indicates that electrostatic effects are likely to be central in catalysis and r e g u l a t i ~ n . l ~ - ~ * We have previously used the modified TanfordKirkwood model to investigate the effects of assembly, conformation, and ligation on the electrostatic properties of These calculations indicated that long-range perturbations of the electro- 301 static properties of the protein result from substrate and nucleotide binding and assembly of the enzyme from its subunits and suggested that a more rigorous analysis was warranted. The results obtained with the MGPB approach, which are in good agreement with experimental titration curves and generally yield individual site pKll2 values in the pH range 2-14, confirm these conclusions. Groups whose pK,,, values differ most from those of model compounds are a t the intersubunit interfaces where they are buried in the low dielectric of the protein. The asymmetry observed in the structure around the approximate non-crystallographic twofold axesz7 is readily apparent in the calculated electrostatic free energy terms. The sitesite interaction energy has a complex dependence on the distance between interacting sites, with a number of interactions between groups separated by distances of up to 11A having site-site interaction free energies greater than 2 kcal . mol-l. Three types of electrostatic interactions are identified 1)salt linkages; 2) interactions involving a single bridging group; and 3) networks of interactions extending over several residues and across interfaces for distances as large as 35 A. The longer range interactions are likely to be important in site-site communication of structural and functional effects. MATERIALS AND METHODS The Model The initial coordinates were those for unliganded ATCase crystallized in the space group P321 at pH 5.8 (Brookhaven Protein Data Bank entry 6AT1)28 with polar hydrogen atoms added with the HBUILD option of XPLOR." To idealize the geometry and minimize the effects of crystal contacts, the 5 ps structure from a molecular dynamics simulation was energy minimized with 120 cycles of conjugate gradient minimization with a 20 kcal . mo1-l constraint on the C" position^.^' The electrostatic terms on the ionizable residues were switched off to prevent formation of incorrect hydrogen bonds during minimization. Molecular dynamics calculations and energy minimization were carried out with X-PLOR. The dielectric constant of the interior of the protein was generally set at 4, while the solvent surrounding the protein was treated as a continuum with a dielectric constant of 80. Other calculations with a dielectric of 20 yielded interaction energies and pK,,, values for residues in the interior of the protein that were on the order of 5 1 0 % smaller. Partial charges and van der Waals radii were taken from the CHARMM19 force field.31 Ions in the solvent were assumed to have a Boltzmann distribution and to obey the Debye-Huckel law.,' The ion exclusion shell was assigned a thickness of 2 A (approximately the radius of the hydrated sodium ion). 302 H. OBEROI ET AL. All calculations were carried out at an ionic strength of 0.2 M and a temperature of 298°K. For the calculations of W,, the entire molecule was mapped on to a three-dimensional grid, with a mesh spacing of 0.5-1.4 A, except for residues <20 A from a titrating site where grids with a spacing of 0.25 A centered on each ionizable residue and including all residues enclosed within a box 40 A on each side were used. AAGB,,, and AAGdipolewere calculated only for the asymmetric unit, consisting of one c chain from each catalytic trimer (c, and c,) and a bridging r subunit dimer (rl and r,); values for the rest of the molecule were generated by making use of the threefold molecular symmetry. Calculations carried out on the whole molecule yielded very similar results. Poisson-BoltzmannEquation Electrostatic potentials were calculated with the MGPB approach' to solving the finite difference form of the nonlinear Poisson-Boltzmann equation Vo[(~(r)V+(r)] - k2(r)sinh[+(r)] + 4~rp(r)= 0 (1) where r is a vector between a point and the origin, +(I)is the electrostatic potential, E(r) is the dielectric constant, k is the modified Debye-Huckel parameter, and p(r) is the charge density function.33 The modified Debye-Huckel parameter is given by (2) where 1 / ~the , Debye screening di~tance,~' is independent of the dielectric constant. The finite difference form of the nonlinear Poisson-Boltzmann is given by: k = K & +ti"" = Cei+i +4 ~ d h C e i + + K~P[I + 433 + $251 +....+3(2n + I)!] each level are carried through only a few iterations. For intergrid transfer, we have used the method of local inversion,38 which interpolates from a coarse to fine grid by locally inverting the discrete operator. The interpolation scheme is that described by Dougla~.~~.~~ Titrating Charges Titration curves were calculated by dividing the charges on the protein into titrating charges, corresponding to ionizable groups on the protein, and nontitrating or background charges that are taken from the CHARMM19 parameter set.3' The charge was assigned to a single atom, based on solvent accessibility41and proximity of other ionizable groups, i.e., to the most exposed atom, or, when the solvent accessibilities were equal, to the atom that is closest t o another titrating atom.' Charges on ionizable residues were adjusted to sum to zero in the uncharged state. Free Energies of Protonation and Titration Curves The method used to calculate the electrostatic contribution to the titration behavior of proteins has been described in Bashford and Gerwert4' and is based on the thermodynamic cycle described in Jorg e n ~ o nThe . ~ ~finite difference procedure is used to obtain the three terms that contribute to the free energy of protonating each ionizable group: AAG,,, the interaction term that results from the polarization of the surrounding medium by the charge; AAGdipole,the term that results from the interaction of the titrating charge with the nontitrating (fixed)charges on the protein; and W,, the term that results from the interaction of the titrating charge with other titrating charges. 1 A A G o r n = - QP[+Pprotein - +Pmodel] 2 (3) +i 1 - - QU[+%rotein- +Umodel] where the nonlinear term is represented as a series, h is the grid spacing in A, +o is the potential a t the central grid point, q, is the charge at this grid point, and the summations are over the six neighboring grid points (i = 1-6).34-36 Multigrid The MGPB method maps the molecule onto nested grids of successively smaller spacing and calculates potentials by iterating through the grids, in both directions, until convergence is achieved.37 A few iterations of a smoother (Gauss-Seidel) are applied a t each level before transferring the correction to a n adjacent level through a restriction (fine to coarse) or prolongation (coarse to fine) operator.8 The multigrid approach reduces the time required for a calculation substantially, since most of the iterations are carried out on coarse grids and the calculations at 2 (4) A A G =~Cqi[+Pprotein ~ ~ ~ ~ ~- +%rotein] i - xqi[+Pmodel - +Umodel] (5) I W, = [QP - Qy](+Pprotein, i - +"protein, i) (6) The model compound is the residue in the solvent in the conformation in which it exists in the protein. Qp and Q" are the charges on the titrating group when the titrating group is protonated and unprotonated, respectively. Qp has a value of 1for cationic groups and 0 for anionic groups, while Q" has values of 0 and - 1for cationic and anionic groups, respectively. qiis the static charge on a nontitrating atom; +p and are the potentials calculated with the titrating +" 303 LONG-RANGE ELECTROSTATIC EFFECTS groups protonated and unprotonated, respectively. Positive values AAG,,,, and AAGdipole correspond to positive values of the difference in the free energy of protonation in the protein and the model compound. The pK of site i is given by where Q$' is the charge of site j in the unprotonated state, Oi is the fractional protonation of site i, and pK,, is the pK, of the residue in the protein with all other sites unprotonated -140 5 Two approaches can be used to calculate Oi. Tanford and R o ~ b yuse ~ ~the Henderson-Hasselbach equation to calculate a set of self-consistent partial charges, while Bashford and Karplus" take a Boltzmann-weighted sum over all possible protonation states at each pH. We calculated individual site titration curves using both the Tanford-Roxby method44 and a hybrid approach45 that uses the Tanford-Roxby equation for residues that have sitesite interactions lower than 0.9 pH units, and the Boltzmann-weighted" summation with MonteCarlo sampling for others.46 Although both approaches give similar results when site-site interactions are weak, the Tanford-Roxby equation often fails to converge when site-site interactions are strong. The statistical thermodynamic approach, while slower, converges under these conditions. Hence the hybrid approach is most efficient. A calculation for the holoenzyme requires approximately 90 hours on an IBM RISC 6000 workstation. Calculating pK,,, values with the Tanford-Roxby method requires about 2 hours; the hybrid scheme requires about 7 hours. The time required increases with the number of strongly interacting groups. Experimental Titration Curves ATCase was isolated from E . coli strain EK1104 transformed with plasmid pEK2 (a gift from E. Kantrowitz, Boston College) as previously described.47Purity, assessed by nondenaturing polyacrylamide gel electroph~resis~~, was greater than 95%)with only c6r4,cg, and aggregates as contaminants. Isoionic, salt-free stock solutions of the protein were prepared by gel filtration and dialyzed against 0.2 M KC1 overnight at 4°C. Protein concentrations were determined spectrophotometrically a t 280 nm using an extinction coefficient of 0.59 mL . (mg . ~ m 1 - l . ~ ~ Potentiometric titrations were carried out on a Metrohm E425 pH-stat at 25°C and a pH691 pH meter calibrated a t pH 4 and pH 10 and 25°C before each experiment. Solutions were degassed for 20 6 7 8 PH 9 10 Fig. 2. Comparison of calculated and experimental titration back titration; -, calculated titracurves. 0, forward titration; 0, tion curves at a grid spacing of 0.6 A; - - - -, titration curves calculated at a grid spacing of 0.3 A. The error bars correspond to the forward titration, which was the average of four experiments. minutes before the experiment and the titration vessel was flushed with nitrogen gas (ultrapure carrier grade) during the experiment. The titrant was standardized at 100 mM HC1 for forward titrations and 100 mM KOH calibrated with potassium acidic phthalate for back titrations. Titration curves were obtained by the blank subtraction method of Tanford,5owith the blank solution consisting of 5.0 ml of the 0.2 M KC1 solution used to dialyze the protein. Protein concentrations ranged from 5 to 10 mg . m1-l; results were independent of protein concentration, although these experiments were not as reproducible as those carried out in this laboratory with monomeric proteins such as ribonuclease A and lysozyme. Error estimates are shown in Figure 2. In most experiments, the protein was titrated from pH 6.5-10.5 and back to 6.5. These limits were chosen because the protein precipitates below 6.5 and loses its cooperativity (reversibly) above pH 10.X51There is no dissociation within this pH range, particularly at the protein concentrations used. RESULTS Comparison of X-Ray and Energy Minimized Structures Differences in the deposited coordinates, the structure from the 5 ps molecular dynamics trajectory, and the energy-minimized structure are small, with root mean square (rms) differences of 0.19 and 0.23 A for backbone and side chain atomic positions, respectively. Thus the results of the electrostatic calculations for all three structures would be expected to be, and, in fact, are very similar, although the energy minimization did eliminate several unfavorable site-site interactions. All the re- 304 H.OBEROI ET AL. sults presented below are for the energy-minimized structure. Comparison With Experimental Data Experimental and calculated titration curves are compared in Figure 2. The enzyme does not dissociate over the pH range examined, although its allosteric properties and hence subunit interactions are a sensitive function of pH (cf. Refs. 15, 16, and 19). These curves were obtained with the same sample of protein and there is a small amount of hysteresis in the experimental curves so that the forward and back titrations differ by 4-10 protons. The basis of this hysteresis warrants further exploration in the future. The forward titration gives better agreement with the titration curve calculated with a grid spacing of 0.3 A, while the back titration agrees better with the curve calculated with a grid spacing of 0.6 b. Both the forward titration and the curve calculated with the finer mesh would be expected to be more accurate; however, a series of calculations carried out with several grid spacings between 0.3 b and 0.6 b indicated that there is not a simple correlation between grid spacing and agreement with experiment. This is probably because the calculations are very sensitive to the relative positions of ionizable atoms and grid points and this is not a simple function of the grid spacing. Both calculated curves diverge from the experimental curves a t high pH, where the average solution structure would not be expected to be the same as the crystal structure, since the enzyme is non-cooperative at pH values above The agreement of the experimental and titration curves indicates that the calculations are able to reproduce the overall titration behavior of the protein with reasonable accuracy. Exact agreement would not be expected since the sensitivity of the catalytic and allosteric properties of the enzyme to pH indicates that the structure is sensitive to pH. The only other experimental result with which the calculations can be compared is the pK, of His 134 in the c chain, which has been reported to be less than 6 in the unliganded catalytic subunit and ligated with N-phosphonacetyl-L-Asp (PALA), based on the observation that the proton resonance to which it was assigned does not titrate over the pH range 6-8.5.5z The calculated pK,,, values are 5.8 and 6.4for the two c chains in the asymmetric unit. Villoutreix et al.53reported a value of -3.5 pH units for this residue in the Bacillus subtilis ATCase structure obtained by structure-based homology modeling of a n ATCase catalytic trimer ligated with PALA. Potential Surfaces Figure 3 shows the electrostatic potential surface at a contour level of 1.5kT/e for a pair of c:r chains within the molecule at pH 2,6, and 11. Although the net charge on the molecule is negative, substrate and nucleoside triphosphate binding sites are positive, as reported previouslyz6 and remain positive over the entire pH range of the titration curve. Since both of the substrates (carbamyl phosphate and L-Asp) are negatively charged, as are the regulatory nucleoside triphosphates, the regions of positive potential will assist in docking and electrically polarizing these ligands a t their binding sites, while the negative potential around other regions of the molecule will tend to prevent nonproductive binding a t other sites. Contributions of AAG,,, and AAG&ipole AAGBorn and AAGdipole for all ionizable residues of the cl and r l chains are plotted in Figure 4.The form of Equations 4 and 5 is such that positive AAGB,,, terms indicate that the solvent accessibility of the titrating group is lower in the protein than in the model compound while positive AAGdipole terms indicate an unfavorable interaction between the protonated species and the static charges of the protein relative to the model compound in solution, for both cationic and anionic groups. Most of the large AAGdipole terms are negative and are paired with positive AAGBo,, terms. The large positive AAGB,,, terms indicate that the solvent accessibility of these anionic groups is low, as expected. The large negative AAGdipol, terms indicate that they are stabilized by the static charges of the protein. As a result of this compensation, the pKadj values of 500 of the 810 ionizable groups differ by less than 3 pH units from the values for the model compounds. Most of the large negative AAGdipole terms (and hence large positive AAGB,,, terms) correspond to anionic groups. The most likely reason for anionic groups having larger AAGdiPle and AAGBorn terms is that they are selectively positioned in the protein so that they can interact favorably with the dipoles of the secondary structural elements. In addition, because the side chains of Glu and Asp are relatively short, they will tend to be closer to the surface of the protein and hence more affected by the dielectric discontinuity. Long-Range Interactions While intuitively one might expect the relationship between the magnitudes of the W, and distance to be hyperbolic, effects such as desolvation and a non-uniform dielectric may lead to a more complex dependence in complex systems.54As shown in Figure 5a, the dependence of W, on distance is roughly biphasic, with a steep dependence over short distances and a shallow dependence over long distances. However, there is a great deal of scatter across the entire range. The calculated points fall below the curve calculated using Coulomb’s law and a dielectric constant of 4, while curves calculated with dielectric constants 25-40 pass through the points. The damping of site-site interactions in the LONG-RANGE ELECTROSTATIC EFFECTS 305 Fig. 3. Potential surfaces around a c:r pair within the molecule at pH 2 (a), 6 (b),and 11 (c). The catalytic chain is in gold; the regulatory chain in green. Blue, positive potentials; red, negative potentials. Potential surfaces were calculated at a contour level of 1.5 kT/e and are shown for only one c and one r chain for clarity. Produced with DAMPS (H. Oberoi and N. Allewell). protein is the result of several factors including interactions of titrating charges with static charges on the protein, screening of interacting charges by the titrating charges that lie between them, and the irregular dielectric boundary. These results support the suggestion that use of a dielectric constant greater than 20 may improve agreement with experiment. 55 All of the interactions with W, values 2 2 kcal . mol-' are listed in Table I, along with distances between the atoms to which charges were assigned and solvent accessibilities. The error in W, is estimated as being in the range of about 1 kcal . mol-l, based on the error in the calculated potential^.^^ Most of the residues involved in these interactions are clustered in four regions of the molecule: in a zone that extends from the active site to the cl:c2 interface, a t the c1:rl and cl:c4 interfaces, and within the r chain. A number of these residues are involved in site-site interactions with several other residues, often over distances greater than 10 A. These strong, long-range electrostatic interactions are of particular interest because they provide a potential mechanism for long-range structural interactions and functional communication within and between subunits. In each c:r pair there are 31 site-site interactions over distances greater than 5 A that result in W, values greater than 2 kcal mol-l. Nine of these strong, long-range interactions (W, 2 2 kcal . mol-l; distance 2 5A) involve ion pairs or three-center clusters (Fig. 6a, i and ii) that are not involved in other strong, long-range interactions. The remainder are part of a network involving 15 residues whose connectivity is mapped in Figure 6a, iii). This network includes Arg c54, Arg c105, and His c134 a t the active site, Glu c109 at the cl-c4 interface, Glu r142 a t the cl-rl interface, and Tyr 98 in c2. The interactions in this network involve a bridging residue; for example, His c134 acts as a bridge between Arg c296 and Asp c129; Asp c129 bridges Glu c109 and Arg c167 in the catalytic chain and Glu r142 in the regulatory chain. As shown in Figure 6b, many of these interactions are complex; for example, Asp c129 interacts with Glu r142 in the regulatory chain and with the active site both directly (with His c134, Arg c105) and indirectly through bridging interactions involving Glu c109 and Arg c167. In many cases, there are multiple connectivities between sites; for example, between His c134 and Asp c129. The extent to which different networks of interactions are used may be influenced by variables such as pH and state of ligation. The site-site interactions of His c134 with neighboring residues are of particular interest since this residue may be involved in the catalytic mechanism.13 His c134 is involved in short, medium, and long interactions with other residues. The shortrange interactions include not only the interaction 3 306 20 H.OBEROI ET AL. I c chain M I I Residue Fig. 4. AAG,,,, and r l chains. (shaded bars) and AAGdipale(hatched bars) for the ionizable residues in the c1 with Arg c105 identified from the crystal structure, but also strong interactions with Tyr c98 (in an adjacent chain of the same trimer), Asp c129, Asp c141, Arg c167, Tyr c226, and Arg c296. If all these residues were fully charged, the net effect of these interactions would be to decrease the pK,, of His c134 by 18 pH units. However, the interactions of these residues with other titrating groups results in a net shift of +10.6 pH units, bringing the calculated pK,,, to 5.8. We have focused on a subset of site-site interactions that have energies greater'than 2 kcal . mol-' and distances greater than 5 A and hence are likely to be involved in long-range structural and functional efforts. However, many other interactions are likely to have structural andlor functional significance. For example, Asp 236 (c4), Glu 109 (cl),Asp 129 (cl),His 134 (cl),and Tyr 98 (c2) trace a pathway with interactions greater than 2 kcale mol-' that links all six c chains in the two catalytic trimers. In addition, an interaction between Glu c109 and Lys r143 with a calculated W, of 6.8 kcal . mol-' links this pathway with the regulatory subunit. Hence the complete network of electrostatic interactions is considerably more com- plex than the restricted set of interactions diagrammed in Figure 6. Asymmetry While the molecular symmetry of ATCase about the crystallographic threefold axis in the P321 crystal form must be exact, there are small but significant differences between the chains related by the approximate non-crystallographic twofold axesz7 (reviewed in Refs. 13 and 14). This asymmetry may be functionally significant since both carbamyl phosphate and nucleoside triphosphates bind with negative c o ~ p e r a t i v i t y . ~ The ~ , ~magnitude ~ of the asymmetry in the electrostatic properties of the two halves of the molecule induced by the structural asymmetry can be assessed from the calculated pK,,, values listed in Table I1 for residues a t the subunit interfaces and active sites of the c l , c6, and rl, r6 chains (the asymmetric unit). Lys c31, Lys c40, Arg c56, Arg c65, Lys c84, Tyr c98, Arg c167, Lys c178, Lys c244, Asp c278, and Asp r39, Asp r55, Lys r139 have pK,,, values that differ by at least 2 pH units. Residues a t the cl:c2 interface have the highest asymmetry with 6 of the 14 ionizable resi- 307 LONG-RANGE ELECTROSTATIC EFFECTS 0 I I I 5 10 15 Distance (A) a. 150 5 100 C Q) 3 CT F LL 50 0 b. Distance (A) Fig. 5. a: Dependence of absolute values of W, (kcal . mol-’) on distance (A). The site-site interactions are shown as open squares. The ewes were calculated with the equation W, = Ckr, where C = 331 and r is the distance in A, and dielectric constants (E) of 4 (-), 25 (- - - -), and 40 (. - . - . -). b: Histograms showing the distribution of distances between interacting sites for various ranges of W,. dues a t this interface having electrostatic terms that differ between the two halves of the molecule by more than 2 pK units. On the other hand, many other residues show very little asymmetry. These calculations therefore allow the ionizable groups a t the active site and subunit interfaces t o be divided into two groups: those that show significant electrostatic asymmetry and are likely to be involved in the negatively cooperative binding of carbamyl phosphate and nucleoside triphosphates, and those that TABLE I: Site-Site Interactions (WJ2 2 kcal-mol-' SA SA Residue (A2) Arg 17 cl 36.2 Glu 37 c l Lys 40 cl Glu 50 c l Arg 54 c l 43.1 9.5 20.5 7.4 His 64 cl Arg 65 cl Glu 86 cl 2.5 8.4 0.3 Asp 90 c l Tyr 98 cl 3.0 1.0 Asp 100 cl Arg 105 c l 21.0 18.9 Glu 109 el 3.4 Arg 113 cl 22.3 Glu 117 cl 12.0 Asp 129 cl 0.0 His 134 c l 0.1 Asp 141 c l 0.0 Glu 149 cl Lys 164 cl 1.8 2.5 Tyr 165 c l 0.0 His 170 c l Tyr 185 c l 0.4 0.1 Glu 216 cl 10.4 Glu 239 c l 7.8 Tyr 240 cl Lys 258 c l His 265 c l 1.5 12.6 0.0 Arg 269 c l 1.4 Glu 272 cl 0.0 Asp 278 cl 7.9 Arg 296 c l 2.9 Arg 41 r l Lys 56 r l 17.8 2.9 Arg 102 r l Arg 128 r l Tyr 140 r l Lys 143 r l 58.8 41.4 0.8 0.0 Residue Asp 153 cl Asp 180 cl Lys 40 c2 Glu 37 c3 Arg 167 c l Glu 86 c2 Tyr 98 c2 His 64 c2 Asp 100 c3 Asp 90 c l Arg 54 c3 Asp 278 c3 Arg 54 c3 Asp 141 c3 His 134 c3 Lys 42 cl Glu 50 cl Arg 167 c l Arg 296 c l Lys 143 r l Asp 236 c4 Asp 129 cl Asp 14 cl Glu 142 r l Lys 139 rl Tyr 140 r l His 134 cl Arg 167 cl Glu 142 r l Arg 105 c l Asp 141 c l Arg 167 cl Tyr 226 cl Arg 296 cl Tyr 98 c2 Tyr 226 cl Tyr 98 c2 Asp 223 cl Tyr 165 c l Tyr 165 c4 Glu 239 c4 Arg 234 el Tyr 165 c4 Glu 239 c4 Tyr 197 cl Glu 221 cl Arg 183 cl His 156 cl His 255 c l Asp 253 cl Tyr 165 c4 Lys 164 c4 Asp 271 cl Asp 223 c l Tyr 226 cl Glu 272 c l Asp 278 e l Tyr 285 cl Asp 90 c2 Tyr 240 c l Asp 271 el Tyr 285 cl Asp 90 c2 Tyr 98 c2 Tyr 226 cl Asp 141 cl Glu 62 r l Gly 8 r l Asp 19 rl Asp 104 r l Glu 144 r l Lys 139 r l Asp 236 cl Asv 236 c4 (A2) 12.5 46.3 9.5 43.1 11.8 0.3 1.0 2.5 21.0 3.0 7.4 7.9 7.4 0.0 0.1 24.5 20.5 11.8 2.9 0.0 8.0 0.0 23.8 0.4 13.8 0.8 0.1 11.8 0.4 18.9 0.0 11.8 0.0 2.9 1.0 0.0 1.0 5.9 0.0 0.7 8.0 28.1 0.7 8.0 8.8 37.7 55.8 4.0 6.0 15.1 0.7 1.4 3.8 5.9 0.0 0.0 7.9 0.2 3.0 1.5 3.8 0.2 3.0 1.0 0.0 0.0 11.1 0.0 20.0 7.8 12.1 13.8 5.4 8.0 w, Distance (kcal . mol-l) (A) 2.0 2.9 2.6 2.7 2.7 7.1 2.8 2.0 6.1 4.8 7.2 2.1 2.9 2.2 2.3 2.6 3.1 4.4 2.3 6.8 2.3 2.3 4.1 3.2 5.1 8.2 2.7 3.6 5.0 5.0 2.8 3.2 2.4 5.7 2.3 17.5 2.2 3.1 7.2 2.0 5.8 2.4 2.1 5.4 2.8 3.0 3.0 2.8 2.6 2.0 6.8 4.1 4.7 6.4 2.8 8.1 8.8 4.3 5.5 2.0 3.4 11.0 2.1 4.2 2.2 3.1 4.5 6.6 7.6 2.4 3.4 7.6 3.8 4.2 3.3 5.4 3.5 3.5 4.8 3.1 7.7 10.0 3.3 5.1 3.1 5.7 7.7 10.7 11.4 4.9 4.7 4.5 7.7 3.7 4.7 8.4 3.7 4.7 4.0 3.2 7.3 5.3 6.2 4.6 9.2 5.9 9.3 6.0 11.4 3.5 10.7 7.8 3.2 8.2 3.1 4.5 8.0 4.5 5.0 3.2 4.8 4.2 4.8 4.6 2.9 3.9 3.1 2.9 7.2 7.0 3.0 5.3 4.6 5.5 3.8 4.5 5.7 6.1 9.8 7.8 4.8 6.6 7.6 3.5 3.3 3.2 4.9 3.6 All the interactions that fulfill this criterion for a given residue in column 1 are listed in column 2. Distances are between ionizable atoms (see methods section). SA, solvent accessibility for residue. 309 LONG-RANGE ELECTROSTATIC EFFECTS cl Arg 17 cl Glu 149 cl Arg 269 c1 His 64 - cl Arg 54 cl Asp 180 cl Asp 223 0) cl Tyr 285 C2 His 64 I c2 Tyr Tir 98 cl Arg 296 cl Asp 141 I Cl Arg 105 -cl cl Lys 164 cl Tyr 165 r l Gly0 His 134 cl Arg 167 r l Asp 19 (li) PI cl A C l Tyr 226 12-rI cl G u 109 cl His 265 c1 Glu 272 - Cl Tyr 240 a. (iii) Fig. 6. (Legend on p. 310.) Glu 142 310 H. OBEROI ET AL. do not show electrostatic asymmetry and are not likely to be directly involved. The possibility that this asymmetry is the result of grid artifacts is ruled out in three ways. Each electrostatic term is calculated three times for the whole molecule, with the titrating group in a different position in the grid each time. Differences in the free energy terms between symmetry-related chains (e.g., c l and c2) were -1% in these calculations. Similarly, calculations in which the position of the molecule in the grid was varied differed by 2-4%. Finally, the structural variability generated by the molecular dynamics run results in differences in the Born and dipole terms and calculated pK,,, values that generally do not exceed one pH unit. In contrast, free energy terms and pK,,, values for the same residue in the two halves of the molecule related by the non-crystallographic twofold axis often differ by several pH units. Effects of Mutations These calculations also make it possible to analyze quantitatively the electrostatic effects in various single site mutations that have been shown to affect function. Two examples are presented below for illustrative purposes. Mutating Asp cl00 to Ala or Asn has been shown to decrease the stability of the catalytic trimer of the enzyme and increase the affinity of the holoenzyme for L - A s ~ . ~This ' was attributed to the elimination of a hydrogen bond to Arg c65 that was identified in the crystal structure. We calculate a value of 6 kcal . mol-' for the W, term for this interaction. However, the calculations also indicate that there is a second interaction with Lys c42 (W, = 2.5 kcal . mol-l, distance = 4.9 A) and suggest that it would be of interest to explore the functional consequences of mutating Lys c42. Arg c54, whose side chain interacts with Glu c86 and whose backbone carbonyl group interacts with the OH of Tyr c98 a t the cl:c2 interface, has also been implicated in the catalytic mechanism, since mutation to Ala results in a significant decrease in activity." The calculated pK,,, for Arg c54 in cl is 6.4 (5.8) and the magnitude of the Arg c54-Glu c86 electrostatic interaction is 7.1 kcal . mol-l. However, the interaction of its sidechain with Glu c86 is absent at the comparable c4:c5 interface in the second trimer, suggesting that the structural basis of its role in catalysis may be more complex. These Fig. 6. a: Connectivities between sites separated by distances of greater than 5 A with W,, terms greater than 2.0 kcal . mol-'. i: Pairwise interactions. ii: Bridging interactions. iii: Long-range interactions. Residues that are involved in long-range interactions also form paiwise interactions with neighboring residues. b: Mapping of residues in a on the molecule. Pairwise interactions (a, i) are shown with dashed lines; three-center interactions (a, ii) with solid lines and circles; and interactions in the network (a, iii) with solid lines. differences in site-site interactions between subunits related by the non-crystallographic symmetry axis may be the cause of the differences in the Born and dipole terms for Glu c86. DISCUSSION ATCase is the largest protein to which the finite difference method has been applied and one of the few for which a complete set of pK,,, values have been calculated. The multigrid method developed in this laboratory allows calculations t o be carried out several orders of magnitude more rapidly than is possible with other methods. Comparison of calculated and experimental titration curves provides a test of the validity of the calculations, and for larger molecules is often the only experimentally accessible quantitative comparison available. The results presented here illustrate the potential of calculations of this type for identifying long-range electrostatic effects that cannot be identified directly from the crystal structure and for analyzing the effects of single site mutants. Calculations based on structures in different states of ligation and different quaternary states are under way, and comparison of these with the results for the unliganded structure will allow us to assess the role of these effects in the catalytic and allosteric mechanisms. AAG,,, and hAGdipolemay be considered to be the static components of electrostatic effects because they depend largely on the local environment of individual ionizable residues rather than their dynamic interactions with other ionizable groups. Hence AAGBornand AAGdipolecan be used to examine the electrostatic consequences of specific structural features. AAG,,, is a measure of the change in the solvation of the charge of the ionizable group that results from incorporating it into the folded protein while AAGdipoleis a measure of the strength of the interaction of the ionizable group with the static (nontitrating) charges in the protein. The relative magnitudes of these terms and their signs yield information about both stability and function. For example, a buried charge with an energetically unfavorable AAG,,, that would normally destabilize the protein and have a highly perturbed pK value may have a AAGdipolethat compensates for the burial so that the group neither destabilizes the protein nor has an altered pK,,,. While the AAGB, and AAGdipoleterms are primarily a measure of the local environment, the sitesite (W,) interaction terms often depend on longrange interactions that are linked to changes in tertiary and quaternary structure. Compensatory changes in the apparent pK,,, of residues such as those illustrated by His c134 and Asp el29 require cooperative changes in the titration behavior of distant residues, which in turn may induce structural changes over several residues. Such cooperative interactions are likely to be important in the trans- 311 LONG-RANGE ELECTROSTATIC EFFECTS TABLE 11. (Continued) TABLE 11. AAG,, and AAGaipOlein Units of pH and Calculated P K ~Values ,~ for Residues at the Active Site and Subunit Interfaces* AAGBorn Interdomain interface, c chain Arg c 17 1.3 Lys c 31 -6.4 -3.0 Arg c 54 -4.3 Arg c 105 -4.1 Asp c 129 10.1 Asp c 141 11.3 9.8 Asp c 153 4.1 Arg c 167 -4.7 -4.1 Lys c 178 -3.2 -2.5 Asp c 180 1.3 Asp c 223 5.6 Tyr c 226 14.1 13.1 Lys c 258 -3.0 Lys c 262 -4.2 His c 265 -8.0 Arg c 269 -4.8 Asp c 278 6.2 5.7 His c 282 -4.6 Tyr c 285 12.2 Interdomain interface, r chain 5.7 Tyr r 77 Arg r 102 -0.8 4.2 Asp r 104 cl:r4 interface -2.4 Lys c 244 -0.1 -7.7 Lys r 143 Glu r 144 3.4 rl:r6 interface 0.6 Glu r 10 0.3 Asp r 39 4.7 -2.4 Arg r 41 -1.0 Arg r 55 -1.4 cl:c4 interface -4.2 Lys c 164 7.5 Tyr c 165 3.0 Asp c 236 4.4 Glu c 239 cl:c2 interface 1.6 Glu c 37 0.8 Lys c 40 -0.8 -1.3 His c 41 -4.3 Arg c 54 Arg c 56 -3.0 -4.0 His c 64 -5.7 -3.8 Arg c 65 -2.9 AAGdipole PK1/2 1.3 0.8 1.2 -2.2 1.0 -20.0 -18.8 -7.9 -1.0 -0.3 1.1 1.1 2.3 0.7 -2.7 -24.8 -18.7 0.5 2.4 -10.8 0.8 -1.8 -3.7 1.0 -14.9 12.7 4.5 8.5 6.4 10.3 7.2 6.4 6.7 7.1 8.1 10.3 8.2 10.3 5.6 6.6 -8.2 -0.4 -0.5 7.7 11.5 7.7 4.6 8.1 8.5 10.2 9.2 5.8 2.3 7.8 -1.2 -1.1 -12.0 -2.6 7.0 9.4 -0.8 0.4 0.3 0.4 0.8 2.5 3.8 4.8 8.9 10.7 12.6 -6.1 -10.3 -6.6 -8.6 0.7 -0.2 0.0 0.4 -2.2 3.9 -3.0 1.6 4.1 4.1 5.4 8.1 0.6 4.9 4.6 9.4 5.2 6.4 11.8 1.4 13.3 (Continued) Lys c 84 Glu c 86 Asp c 90 Tyr c 98 Asp c 100 Lys c 232 Arg c 269 c1:rl interface Glu c 109 Arg c 113 Glu c 117 Glu c 204 Glu r 119 Arg r 130 Lys r 139 Tyr r 140 Glu r 142 Lys r 143 Active site Arg c 54 Lys c 84 Arg c 105 His c 134 Arg c 167 Arg c 229 Glu c 233 -1.9 -0.1 7.9 5.2 4.7 7.5 4.0 -1.0 -2.3 -4.8 -2.2 -1.6 -4.4 -2.5 -1.6 -9.0 0.3 -2.2 -1.6 0.8 6.3 -3.0 4.3 0.7 6.3 -0.6 -3.0 -1.9 8.8 6.3 -7.7 -16.9 -1.5 -9.6 -0.6 -10.8 0.4 -5.5 -3.4 -16.0 -11.1 -12.0 -4.3 -1.9 0.1 -4.1 -7.3 -4.7 -4.1 -2.3 2 -2.2 -2.2 -1.6 1.0 -3.8 -0.3 1.1 -1.4 -3.6 6.3 8.9 8.6 6.9 7.0 9.1 8.0 7.3 9.3 10.2 8.7 4.5 12.0 2.2 5.2 3.8 6.4 6.3 8.9 10.3 6.4 8.1 10.3 9.3 2.8 *1 kcal = 1.4pH units. Errors in free energies are on the order of 1 kcal . m-'. Values for the two halves of the molecule related by the approximate non-crystallographic twofold axes are shown when values of AAG,,,,, AAGdipoler and pK,,, differ by more than 2 pH units. In these cases, the upper value corresponds to cl or rl and the lower value to c4 or r6. Calculated pK,,, values outside the range 2-14 are not listed. mission of signals over long distances in ATCase and other multisubunit proteins. While pK,,, values are very sensitive to how charges are assigned and small changes in side chain positions, W, values are less sensitive t o charge assignment and hence are a more reliable parameter for analysis. Although the inherent asymmetry of the structure would suggest that its electrostatic properties would be asymmetric, the magnitude of the asymmetry in the calculated electrostatic free energy terms is striking. These differences persist and are even more striking in the PALA-liganded structure (Oberoi et al., in preparation). Differences in the calculated pK,,, values of individual residues in cl-c6 or rl-r6 related by the approximate non-crystallographic twofold axis may result from differences in AAGBorn,AAGdipole(and therefore pK,,) or the sitesite interaction (W,) term. The fact that this electrostatic asymmetry is demonstrable in two struc- 312 H. OBEROI ET AL. tures with very different quaternary states and charge configurations attests to the robustness and reliability of the calculations. It also suggests that comparisons of the electrostatic properties of the liganded and unliganded structure may provide additional insights into the functional significance of the asymmetry. Although subject to the limitations of continuum models,4.6,60grid methods allow considerable refinement of the physical model over previous approaches. Several additional steps can be taken to reduce the error further. The effect of varying the dielectric constant either throughout the grid or within regions of it remains to be investigated. While results obtained with the modified Tanford-Kirkwood model suggested that the exact value of the dielectric constant within the protein did not have much effect:’ evidence that a value of 20, rather than the generally accepted value of 4, yields better agreement with experiment has recently been presented.55 Since direct experimental determinations and simulations yield values near 4,60 the higher value may compensate for neglecting structural features such as bound solvent, protein motion, and ion binding. It is also likely that the increased screening due to the higher dielectric constant may compensate for side chain mobility and resonance, since neglecting these effects would effectively exaggerate the magnitude of the interaction terms and the calculated pK,,, values. In other studies the effect of changing the protein dielectric constant to 10 was shown to have a modest effect on calculated pK shifts while a dielectric constant of 78 produced a reduction of nearly 50%in the calculated pK shifts.62Such effects may in part be responsible for the discrepancy between the calculated and experimental titration curves.63 Standard continuum models generally assign a high dielectric constant to all regions of the molecule accessible to a spherical probe. Since this includes regions internal to the molecule that may contain bound water, the dielectric constant in these regions may be lower. Developing methods to include water in the calculations, as was done in Yang et al.45will improve the accuracy of the calculations. We find that including bound solvent in the internal cavity of adipocyte lipid binding protein produces significant changes in calculated pKlI2 values (Oberoi et al., in preparation). A second significant source of error is the neglect of side chain motion. The first systematic attempt to examine the consequence of side chain motion using a continuum model was by Wendoloski and Matthew,64 who used the modified Tanford-Kirkwood model and the Tanford-Roxby equation to model pH dependence. Differences in calculated pK,,, values among several tuna cytochrome c structures from a molecular dynamics run were of the same order of magnitude as those between finite-difference and Tanford-Kirkwood calculations. You and B a ~ h f o r d ~ ~ have recently used a finite difference solution of the Poisson-Boltzmann solution to study the effects of conformational flexibility. We have shown for hen egg white lysozyme that crystal structures derived by different refinement protocols also have calculated pK,,, values that differ by > 10 pH units.’ The substantial differences in the results reported here for the two halves of ATCase are a third indication of the sensitivity of these calculations to side chain positions. Although computationally intensive, the incorporation of side chain motion is in principle straightforward. A third issue is the method used to assign charges to individual atoms. We obtained good agreement with experiment when charges were assigned to a single atom on the basis of solvent accessibility and proximity of neighboring charges, Although this approach has been reasonably successful with the molecules to which we have applied it to date, methods with stronger theoretical justification that yield equally good agreement with experiment are clearly desirable. The results presented here suggest a role for three classes of electrostatic interactions in intramolecular communication. One class of interactions, over distances of 3-5 A, involves pairs of residues that interact directly and is generally discernible from direct observation of the crystal structure. A second class of interactions requires a “bridging” residue and therefore involves three residues and extends over distances of 6-15 A. These interactions are not readily apparent from the structure since at least one pair of residues may be 6-8 apart. The complete set of ion pair and bridging interactions forms a third type of pathway, a network that extends across and between subunits and provides a potential mechanism for generating long-range structural changes and transmitting allosteric effects. When calculations for various liganded structures become available for comparison with these calculations, it should be possible to assess the role of these longrange interactions in the catalytic and allosteric mechanisms of the enzyme in greater detail. The actual network of long- and short-range interactions within the molecule is complex, with multiple pathways between many ionizable groups. The extent to which any given pathway is used will depend on pH, since the magnitudes of site-site interactions vary with pH, and also with the state of ligation of the protein. In addition, the complexity of the network suggests that communication between different parts of the molecule may be plastic, so that a mutation in one pathway may be compensated for by communication through a different pathway. This would in turn mean that failure of a mutation to alter function cannot be interpreted to mean that that residue does not have a functional role. LONG-RANGE ELECTROSTATIC EFFECTS Networks of interactions over long intramolecular distances such as we have described here may serve the dual purpose of maintaining subunit interactions and channeling allosteric and cooperative signals. To validate and understand the role of the types of interactions illustrated in this paper will require the design of new sets of single site mutants. Most single site mutations to date have been designed to investigate the roles of residues that are either at the active or effector site or are involved in subunit-subunit interactions. Electrostatic modeling has the potential to identify residues that may be equally important from a structural and/or functional point of view, but whose role is less obvious. ACKNOWLEDGMENTS This work was supported by NIH grant DK-17335 and the University of Minnesota. The facilities of the Minnesota Supercomputer Institute were used for some calculations. We thank Drs. Vince J. LiCata and Victor Bloomfield for critical discussion, and Hiroki Morizono and Dr. Meena Hariharan for help with the figures. REFERENCES 1. Stone, R.L., Darion, J.E. Protein-tyrosine phosphatases. J. Biol. Chem. 26931323-31326,1994, 2. Samuelson, L.E., Klausner, R.D. Tyrosine kinases and tyrosine-based activation motifs. Current research on activation via the T cell antigen receptor. J. Biol. Chem. 267: 25913-25916,1992. 3. Johnson L.N., Barford, D. The effects of phosphorylation on the structure and function of proteins. Annu. Rev. Biophys. Biomol. Struct. 22:199-298, 1993. 4. Warshel A., Russell, S.T. Calculation of electrostatic interactions in biological systems and in solution. Q. Rev. Biophys. 17:283-422, 1984. 5. Matthew, J.B., Gurd, F.R.N., Garcia-Moreno, B., Flanagan, M.A., March, K.L., Shire, S.J. pH-dependent processes in proteins. CRC Crit. Rev. Biochem. 18:91-197, 1985. 6. Sharp K.A., Honig, B.H. Electrostatic interactions in macromolecules: Theory and applications. Annu. Rev. Biophys. Biophys. Chem., 19:301-332,1990. 7. Allewell, N.M., Oberoi, H. Electrostatic effects in protein folding and function. Methods Enzymol. 202:3-19, 1991. 8. Oberoi, H., Allewell, N.M. Multigrid solution of the nonlinear Poisson-Boltzmann equation and calculation of titration curves. Biophys. J., 65:48-55, 1993. 9. Holst, M., Kozack, R.E., Saied, F., Subramaniam S. Treatment of electrostatic effects in proteins: Multigrid-based Newton iterative method for solution of the full nonlinear Poisson-Boltzmann equation. Proteins 18:231-45, 1994. 10. Bashford, D., Karplus, M. Multiple-site titration curves of proteins: An analysis of a n exact and approximate methods for their calculation. J. Phys. Chem. 95:9556-9561, 1991. 11. Allewell, N.M. Escherichia coli aspartate transcarbamoylase: Structure, energetics, and catalytic and regulatory mechanisms. Annu. Rev. Biophy. Biophys. Chem. 18:7192, 1989. 12. Herve, G. Aspartate transcarbamylase from Escherichiu coli. In: “Allosteric Enzymes.” Herve, G. (ed.). Boca Raton, F L CRC Press, 1989:61-79. 13. Lipscomb, W.N. Activity and regulation in aspartate transcarbamylase. In: “Regulation of Proteins by Ligand.” The Robert A. Welch Foundation Conferences on Chemical Research XXXVI, Houston, Texas, 1992:103-143. 14. Lipscomb, W.N. Aspartate transcarbamylase from Escherichiu coli: Activity and regulation. Adv. Enzymol. Relat. Areas Mol. Biol. 6857-151, 1994. 313 15. Yuan, X., LiCata, V.J., Allewell, N.M. Effects of assembly and mutation of residues outside the active site on the functional pH dependence of E. coli. aspartate transcarbamylase. J . Biol. Chem. 271:1285-1294, 1996. 16. Pastra-Landis, S.C., Evans, D.R., Lipscomb, W.N. The effect of pH on the cooperative behavior of aspartate transcarbamylase from E. coli. J . Biol. Chem. 253:4624-4630, 1978. 17. Gerhart, J.C., Pardee, A.B. Aspartate transcarbamylase, a n enzyme designed for feedback inhibition. Fed. Proc. 23: 727-735, 1964. 18. Leger, D., Herve, G. Allostery and pK, changes in aspartate transcarbamoylase from Escherichia coli. Analysis of the pH dependence in the isolated catalytic subunits. Biochemistry 27:4293-4298, 1988. 19. Turnbull. J.L.. Walsrou, G.L.. Schachman, H.K. Ionisation of amino acid residues involved in the catalytic mechanism of aspartate transcarbamylase Biochemistry 31:65626569,1992. 20. Allewell, N.M., Friedland, J., Niekamp, K. Calorimetric analysis of aspartate transcabamylase from Escherichia coli. Binding of cytosine 5‘4riphosphate and adenosine 5‘triphosphate. Biochemistry 14:230-234, 1975. 21. Knier, B.L., Allewell, N.M. Calorimetric analysis of aspartate transcarbamylase from Escherichiu coli; binding of substrate and substrate analogs to the native enzyme and catalytic subunit. Biochemistry 17:784-790, 1978. 22. Allewell, N.M., Hofmann, G.E., Zaug, A., Lennick, M. Bohr effect in Escherichiu coli aspartate transcarbamylase. Linkages between substrate binding, proton binding, and conformational transitions. Biochemistry 18:30083015, 1979. 23. Burz, D.S., Allewell, N.M. Interactions of ionizable groups in Escherichia coli aspartate transcarbamylase with adenosine and cytidine 5’-triphosphates. Biochemistry 21: 6647-6655, 1982. 24. McCarthy, M.P., Allewell, N.M. Thermodynamics of assembly of E. coli aspartate transcarbamylase. Proc. Natl. Acad. Sci. USA 80:6824-6828, 1983. 25. Glackin, M.P., McCarthy, M.P., Mallikarachchi, D., Matthew, J.B., Allewell, N.M. Electrostatic interactions in the assembly of Escherichiu coli aspartate transcarbamylase. Proteins 5:66-77, 1989. 26. Glackin, M.P., Matthew, J.B., Allewell, N.M. Electrostatic effects and allosteric regulation in E. coli aspartate transcarbamylase. In: “Conformations and Forces in Protein Folding.” Nall, B.T., Dill, K.A. (eds.). AAAS, 1991:1827. 27. Kosman, R.P., Gouaux, J.E., Lipscomb, W.N. Crystal structure of CTP-ligated T state aspartate transcarbamoylase at 2.5 A resolution: Implications for ATCase mutants and the mechanism of negative cooperativity. Proteins 15: 147-176, 1993. 28. Stevens, R.C., Gouaux, J.E., Lipscomb, W.N. Structural consequences of effector binding to the T-state aspartate carbamoyltransferase: Crystal structure of unligated and ATP- and CTP-complexed enzymes at 2.6 A resolution. Biochemistry 297691-7701, 1990. 29. Briinger, A.T. “X-PLOR, Version 3.1. A System for X-Ray Crystallography and NMR.” New Haven, C T Yale University Press 1992. 30. Powell, M.J.D. Restart procedures for the conjugate gradient method. Math. Progr. 12:241-254, 1977. 31. 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. 4:187-217, 1983. 32. Debye, P., Hiickel, E. On the theory of electrolytes: Freezing point depression and related phenomena. In: “The Collected Works of Peter.” Debye J.W. (ed.). New York: Interscience, 1923:217-263. 33. McQuarrie, D.A. “Statistical Mechanics.” New York: Harper & Row, 1976. 34. Warwicker, J., Watson, H.C. Calculation of electric potential in the active site cleft due to a-helix dipoles. J . Mol. Biol. 155:53-62, 1982. 35. Klapper, I., Hagstrom, R., Fine, R., Sharp, K.A., Honig, B.H. Focusing of electric fields in the active site of Cu-Zn superoxide dismutase: Effects of ionic strength and aminoacid modification. Proteins 1:47-59, 1986. 314 H. OBEROI ET AL. 36. Jayaram, B., Sharp, K.A., Honig, B.H. The electrostatic potential of B-DNA. Biopolymers 28975-993, 1989. 37. Briggs, W.L. “A Multigrid Tutorial.” Philadelphia: SIAM, 1988. 38. Hyman, J.H. Mesh refinement and local inversion of elliptical differential equations. J . Comput. Phys. 23:124-134, 1977. 39. Douglas, C.C. Multi-grid algorithms for elliptic boundaryvalue problems. Ph.D. thesis, Yale University, New Haven, 1982. 40. Douglas, C.C., Douglas, J.A. Unified convergence theory for abstract multigrid or multilevel algorithms, serial and parallel. SIAM J. Numer. Anal. 30:136-158, 1993. 41. Lee, B.K., Richards, F.M. The interpretation of protein structures: Estimation of static accessibility. J . Mol. Biol. 55~379-400, 1971. 42. Bashford, D., Genvert, K. Electrostatic calculations of the pKa values of ionizable groups in bacteriorhodopsin. J . Mol. Biol. 224473-486, 1992. 43. Jorgenson, W.L. Free energy calculations: A breakthrough for modelling organic chemistry in solution. Acc. Chem. Res. 22:184-189, 1989. 44. Tanford, C., Fbxby, R. Interpretation of protein titration curves. Application to lysozyme. Biochemistry 11:21922198, 1972. 45. Yang, A.-S., Gunner M.R., Sampogna, R., Sharp, K., Honig, B.H. On the calculation of pK,s in proteins. Proteins 15252-265, 1993. 46. Beroza, P., Fredkin, D.R., Okamura, M.Y ., Feher, G. Protonation of interacting residues in a protein by a Monte Carlo method: Application to lysozyme and photosynthetic reaction center of Rhodobacter sphaeroides. Proc. Natl. Acad. Sci. USA 88:5804-5808, 1991. 47. Nowlan, S.R., Kantrowitz, E.R. Superproduction and rapid purification of Escherichia coli aspartate transcarbamylase and its catalytic subunit under extreme derepression of the pyrimidine pathway. J . Biol. Chem. 260:1471214716, 1985. 48. Davis, B.J. Disc electrophoresis. 11. Method and application to human serum proteins. Ann. N.Y. Acad. Sci. 121: 404-427,1964. 49. Gerhart, J.C., Holoubek, H. The purification of aspartate transcarbamylase of Escherichiu coli and separation of its protein subunits. J . Biol. Chem. 242:2886-2892, 1967. 50. Tanford C. In: “Electrochemistry in Biology and Medicine.” Shedlovsky,T. (ed.). Hydrogen Ion Titration Curves of Proteins. New York: John Wiley & Sons, Inc., 1955248-265. 51. Weitzman, P.D.J., Wilson, I.B. Studies on aspartate transcarbamylase and its allosteric interaction. J. Biol. Chem. 2415481-5488, 1966. 52. Kleanthouse, C., Wemmer, D.E., Schachman, H.K. The role of a n active site histidine in the catalytic mechanism of aspartate transcarbamylase. J. Biol. Chem. 263:1306213067, 1988. 53. Villoutreix B.O., Spassov V.Z., Atanasov B.P., Hewe G., Ladjimi M.M. Structural modeling and electrostatic properties of aspartate-transcarbamylasefrom Sacchuromyces cereuisiae. Proteins 19:230-243, 1994. 54. Jayaram, B. Free energy of solvation, interaction and binding of arbitrary charge distributions imbedded in a dielectric continuum. J. Phys. Chem. 985773-3777, 1994. 55. Antosiewicz, J., McCammon, J.A., Gilson, M.K. Prediction of pH-dependent properties of proteins. J . Mol. Biol. 238: 415-436, 1994. 56. Gilson, M.K., Sharp, K.A., Honig, B.H. Calculating the electrostatic potential of molecules in solution: Method and error assessment. J. Comp. Chem. 9:327-335, 1988. 57. Tondre, C., Hammes, G.G. Interaction of aspartate transcarbamylase with 5-bromocytidine 5‘-tri-, di- and monophosphates. Biochemistry 13:3131-3136, 1974. 58. Rosenbusch, J.P., Griffin, J.H. Carbamyl phosphate binding to aspartate transcarbamylase. J . Biol. Chem. 248: 5063-5066, 1973. 59. Baker, D.P., Kantrowitz, E.R. The conserved residues glutamate-37, aspartate-100 and arginine-269 are important for the structural stabilization of Escherichia coli aspartate transcarbamylase. Biochemistry 32:10150-10158, 1993. 60. Gilson, M.K., Honig, B.H. The dielectric constant of a folded protein. Biopolymers 252097-2119, 1986. 61. Matthew, J.B., Hanania, G.I.H., Gurd. F.R.N. Electrostatic effects in hemoglobin: hydrogen ion equilibria in human deoxy- and oxyhemoglobin A. Biochemistry 18: 1919-1928,1979. 62. Gilson, M.K., Honig, B.H. Energetics of charge charge interactions in proteins. Proteins 3:32-52, 1988. 63. You, T.J., Bashford, D. Conformation and hydrogen ion titration of proteins: A continuum electrostatic model with conformational flexibility. Biophys. J. 69:1721-1733, 1995. 64. Wendoloski, J.J., Matthew, J.B. Molecular dynamics effects on protein electrostatics. Proteins 5313-321, 1989.