INTRODUCTION

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

For many years, one of the main applications of the hydrodynamic properties of macromolecules in solution has been the determination of the conformation or shape of proteins. Before the wide availability of x-ray crystallography, solution properties such as diffusion or sedimentation coefficients and intrinsic viscosities were customarily employed to determine the size and shape of proteins, in terms of a simple hydrodynamic model, the revolution ellipsoid. From one or more properties, it is possible to determine the anisometry (axial ratio) and the degree of hydration of the protein. This classical approach is well described in textbooks (Tanford, 1961; Cantor and Schimmel, 1980).

In the 1980s that situation changed, owing essentially to two developments. On the one hand, x-ray crystallography emerged as a powerful and widespread tool for the determination of protein structure. This technique currently yields the size and shape of the protein, but not just as a low-resolution structure with some overall dimensions in three spatial directions; instead, the position in space (Cartesian coordinates) of nearly every atom in the protein can be determined. There has always been doubt about whether the structure determined in the crystal is the same as that in solution, i.e., under physiological conditions. More recently, newer techniques based on multidimensional NMR have been developed as alternative tools for structural determination in solution.

On the other hand, the last two decades have seen the development of new theories and computational methodologies for the calculation of hydrodynamic properties of macromolecules (see, for instance, García de la Torre and Bloomfield, 1977a-c, 1978, 1981; García de la Torre, 1989; García de la Torre et al., 1994b). These achievements enabled a more detailed modeling of bioparticles and a more rigorous treatment of the hydrodynamics. Another remarkable circumstance is the improvement, during that period, of experimental techniques such as those that monitor rotational motion by fluorescence spectroscopy (Dale and Dale, 1985) and dynamic NMR (Palmer et al., 1996), or the modern developments of analytical ultracentrifugation (Harding et al., 1992; Schuster and Laue, 1994) and viscometry (Harding, 1997).

It was evident that the advances in the determination of protein structure and those in biomolecular hydrodynamics could be linked in an attempt to relate high-resolution structure and hydrodynamic behavior. Specifically, there was the possibility of using atomic coordinates to build hydrodynamic bead models. This was the motivation for various studies that appeared during that time, particularly over the last decade (vide infra). From a global consideration of all those works, two impressions are apparent to us. First, there are some disparities among the various studies, in the modeling procedure as well as in the physical and computational treatment. This can be understood, because the same diversity indeed exists in macromolecular hydrodynamics, but one should be cautious about which model and treatment are most suitable for each case, particularly for globular proteins. The second impression is that some of those works overlooked previous publications on this problem, thus precluding critical comparisons with other approaches, as one can appreciate looking at some recent publications.

Thus we undertook, as an important aspect of the present work, the task of compiling and summarizing the various studies on predictions of hydrodynamic properties of globular proteins of which we are aware. The aim is not a critical, comprehensive review; rather, we placed our emphasis on the analysis of the methodologies employed. In recent years some advances have been made in the understanding of difficult aspects of bead model hydrodynamics, such as volume corrections for rotation (Carrasco and García de la Torre, 1999b) and viscosity (García de la Torre and Carrasco, 1998) or bead overlapping (Carrasco et al., 1999). Thus it is pertinent to look at previous studies of protein modeling to learn what can be employed and refined and what should be avoided.

In a recently published work, of which this paper is an immediate application, we have discussed several bead modeling strategies (Carrasco and García de la Torre, 1999a). We noted that shell models, in the spirit of the pioneering studies of Bloomfield and co-workers (Bloomfield et al., 1967a; Filson and Bloomfield, 1967a; Bloomfield and Filson, 1968), avoid the defects of other strategies and provide a basically correct and general method of predicting the hydrodynamic properties of arbitrarily shaped particles. In the present work we apply the shell modeling method to globular proteins with atomic-level structural details. The methodology has been computationally implemented up to the point of developing a public-domain computer program, HYDROPRO, which, starting from just a Protein Data Bank (Abola et al., 1987) or similar file containing the atomic coordinates, gives a full set of solution properties. These include the primary solution properties calculated by our HYDRO program (García de la Torre et al., 1994b): translational diffusion coefficient, D_t; sedimentation coefficient, s; rotational diffusion coefficient, D_r; relaxation times, ; intrinsic viscosity, []; and radius of gyration, R_g. This set of calculated quantities can even be greatly expanded by our ancillary program SOLPRO (García de la Torre et al., 1997, 1999).

We have applied our methodology to a large group of 13 globular proteins, so that a sound knowledge of protein hydrodynamics can be gathered. In this field, an essential (and elusive) aspect is hydration, which exerts a noticeable influence on hydrodynamic behavior and properties of biopolymers, particularly proteins. The extensive scope of our work, with regard to both the number of proteins studied and the variety of properties that are calculated, along with the consideration of macromolecular structure at the atomic level, permits us to reach some conclusions about this difficult aspect and makes it possible to provide some guidelines for handling hydration, along with atomic-level structures, in the prediction of solution properties of globular proteins.

	REVIEW OF MODELS

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

In this section we mention and comment on various procedures and applications of hydrodynamic modeling of quasirigid macromolecules with structural detail that have been described in the literature. They are grouped by type of model. A chronological list with further information is given in Table 1.

Bead per residue models

An immediate possibility for modeling macromolecules consists of replacing each repeating unit (each amino acid residue in a protein) with a frictional element in the hydrodynamic model. Thus the model is a chain of beads that follows the path of the rigid backbone of the macromolecule.

This approach was attempted many years ago in our group, in an unpublished work (Jiménez, 1980; Tirado, 1982). Amino acid residues were replaced with identical beads centered at the -carbons of some proteins with bead radius 1.8 Å, such that neighbor beads are tangent and do not overlap. An interesting detail was the representation of amino acid side chains by one or two more beads of the same size. This model failed, giving results that indicated that the size of the model was insufficient. It is now clear that the reason was the ignorance of the hydration layer. In further refinement of this model, Tirado et al. (1990) proposed to represent adjacent water molecules with smaller beads, with some improvement in the results.

The most successful applications of bead per residue (BPR) models have been described by Antosiewicz and Pörschke (1989, 1993, 1995). In their models of proteins (as in similar models of tRNA), they replace amino acid residues with beads with a radius of ~5-6 Å. As rotational properties were included in the calculations, the volume correction (García de la Torre and Rodes, 1983) was made in a special form, specifically adapted for models with overlapping beads (Antosiewicz and Pörschke, 1989). A model employed to describe DNA, in which nucleotide units are replaced by tangent, nonoverlapping beads (García de la Torre and Horta, 1976; García de la Torre et al., 1994a), also belongs to this category.

A clear advantage of BPR models is that the number of elements in the model is much smaller than in bead-per-atom (BPA) models (see below), and therefore they can be applied to large molecules, such as some high-molecular-weight proteins (Hellweg et al., 1997). As a disadvantage, we mention the fact of bead overlapping in most versions of this method. The use of the Rotne-Prager hydrodynamic interaction tensor for equal overlapping beads (Rotne and Prager, 1969) may lead to correct results for translational coefficients, and the Antosiewicz-Pörschke volume correction (Antosiewicz and Pörschke, 1989) for rotation works well in the cases in which it has been tested. However, there is no extensive experience with this kind of modeling, and its performance in the calculation of intrinsic viscosities has not been reported. Indeed, volume corrections for intrinsic viscosity are somehow problematic (García de la Torre and Carrasco, 1998), and, in general, the theoretical basis for the treatment of overlapping in bead models is not firmly based (Carrasco et al., 1999).

Bead per atom models

The ultimate level of detail in the hydrodynamic model is reached if each atom in the molecule is represented by a frictional element in the model. Atoms can be extended, in the sense of including their bonded hydrogens. Actually, the BPA strategy has been successfully used to predict translational diffusion coefficients of small molecules in organic solvents, using BPA models that are built directly from their atomic structure (Espinosa and García de la Torre, 1987; Pastor and Karplus, 1988). The friction of the atomic beads must be adequately parameterized, which can be done in two ways: 1) by giving a constant hydrodynamic radius to every extended atom, or 2) by assigning a hydrodynamic radius to each atom that depends on its accessible surface area. Venable and Pastor (1988) have made a detailed study of several variations of the BPA method applied to some small or medium-sized proteins. While the differences in the results arising from different parameterizations of atomic friction were found to be less relevant (with effective atomic stick-hydrodynamic radii of, typically, 1 Å, or up to 1.5 Å), the need to include the effect of hydration was quite evident. This was done by adding further 1.6-Å (stick) beads representing water molecules on the surface of the protein, in a variable number, with a maximum corresponding to a full monolayer. Thus it was possible to adjust well the results to experimental values of translational and rotational diffusion of bovine pancreatic trypsin inhibitor (BPTI), lysozyme, and ribonuclease. When the water beads are included, the Venable-Pastor model is somehow a superposition of the BPA model and the shell model (see below). It is also noteworthy that their effective hydrodynamic stick radius is ~0.7 Å for protein atoms, or 1 Å for waters; therefore bead overlapping is negligible in their models.

Filling models

If the rigid macromolecule is regarded as a compact, globular particle with a definite contour, then there is an immediate possibility of building a bead model that consists of filling the volume enclosed by that contour. In the AtoB procedure proposed by Byron (1997, 1999), the contour is discretized by superimposing a cubic grid on the three-dimensional structure of the protein. The sides of the cubelets determine the resolution of the procedure. In those cubelets where the center of mass of at least one amino acid residue falls a bead is placed, the radius of the bead being such that the bead volume equals the anhydrous volume of the residues assigned to the cubelet. Thus the interior of the protein model is filled with beads of various sizes that show an appreciable degree of overlapping. This model still has to be modified to account for hydration, which is done by uniform expansion, increasing bead coordinates and radii by a given factor that is regarded as an adjustable parameter. Then our HYDRO computer program (García de la Torre et al., 1994b) was used for the hydrodynamic calculations. Zipper and Durchschlag (1997, 1998) have employed a similar procedure. This modeling technique has yielded reasonably good results for translational coefficients but fails remarkably in the prediction of rotational coefficients and intrinsic viscosity.

In our recent work we discussed the inadequacy of filling models as hydrodynamic models. Even if the internal beads are equal and nonoverlapping, the application of the standard volume correction (García de la Torre and Rodes, 1983), which was built within HYDRO, yields clearly erroneous results for rotation and viscosity (Carrasco and García de la Torre, 1999a). Bead overlap, particularly between beads of unequal size, adds further problems to the description of hydrodynamic interactions, because the interaction tensor devised for this circumstances (Carrasco et al., 1999) is just an ad hoc patch that should be used when overlap within a model takes place occasionally but not in a generalized manner. The filling strategy may be harmless for translational properties but adds unnecessarily computing effort, because it takes into account internal beads that are highly shielded and do not contribute to the hydrodynamic behavior.

Filling models are, on the other hand, correct and useful in practice, for the calculation of solution properties related to radiation scattering, such as the radius of gyration. This is because the whole volume of the particle contributes to the scattering properties. Thus filling models with equal scattering elements, arranged in a cubic lattice (the so-called cube methods) have been proposed to calculate the angular dependence of scattering intensities and the radius of gyration of globular proteins (Müller, 1983; Pavlov and Fedorov, 1983; Pavlov et al., 1986). In an extension of this methodology to hydrodynamic particles, Müller proposed that all the internal beads should be removed, thus transforming the filling model into a shell model (Müller et al., 1983a,b). Models classifiable as the filling kind have been used by other authors to predict scattering intensities by means of the Debye formula; for a recent example see Chacón et al. (1998).

Shell models

Hydrodynamic friction takes place at the surface of the macromolecule. It is therefore the molecular surface that has to be modeled. This physical idea was proposed in the pioneering bead-modeling studies of Bloomfield et al. and tested for models of simple shapes (Bloomfield et al., 1967a; Filson and Bloomfield, 1967, 1968; Bloomfield and Filson, 1968). Small beads are arranged in a shell that describes as closely as possible the molecular surface. Eventually, results for various bead sizes can be extrapolated to zero bead radii, thus making the model converge to a smooth surface. In our preceding work (Carrasco and García de la Torre, 1999a) we have compared the various bead modeling strategies, showing that a shell model gives the best possible description of the hydrodynamic behavior of the particle.

Indeed, shell modeling was the choice in an earlier prediction by Teller et al. (1979) of hydrodynamic properties of globular proteins from atomic structures. These authors calculated translational properties with shell models constructed with equal, nonoverlapping spheres of radius 1.4 Å that coat the protein surface (this is a nonextrapolated shell model). As the model that we are proposing in this work belongs to the shell (SHE) type, we postpone further descriptions and comments on this methodology.

Finite elements models

Instead of the Kirkwood-Riseman-Bloomfield framework (Riseman and Kirkwood, 1950; Kirkwood, 1954; Bloomfield et al., 1967a,b; Bloomfield, 1968) used in bead modeling, the hydrodynamic properties of rigid particles can be predicted in terms of the Youngren-Acrivos treatment (1975a,b), in which the particle's surface is paved with platelets, the role of which is similar to those of beads in the other treatment. The hydrodynamic description is then made using finite element (FEL) methods, which are common to other physical problems. Thus an interesting feature of FEL models is that they also serve to describe electrostatic properties along with hydrodynamics, which makes them applicable, for instance, to electrophoresis (Chae and Lenhoff, 1995; Allison and Tran, 1995; Allison et al., 1999). Although the two frameworks are apparently quite different, a careful comparison of bead models and FEL models may reveal some equivalences of the two techniques in specific, limiting cases, as illustrated by Allison (1998) with respect to the intrinsic viscosity. An in-depth comparison of bead models and FEL models is a suggestive project, but it is beyond the scope of the present paper, in which we concentrate on the various bead methodologies.

Nonetheless, FEL methods should be cited, in the context of this revision, because they have been applied to predict hydrodynamic properties of globular proteins with atomic-level detail (Brune and Kim, 1993; Chae and Lenhoff, 1995; Allison and Tran, 1995; Zhou, 1995). These methods are somehow analogous to the shell methods in bead modeling, because it is the surface of the particle that is described in the model. Indeed, an extrapolation to zero platelet size has been suggested for FEL models (Allison, personal communication), in clear analogy to extrapolated SHE models. An important observation from the FEL calculations of globular proteins, related to those obtained in the present study (vide infra), is that to reach agreement with experimental data, the hydrodynamic particle is delimited, not by the molecular surface, but by the external surface of a surrounding layer, the thickness of which may be from 1 Å (Zhou, 1995) to 3 Å (Allison and Tran, 1995).

Molecular dynamics simulation

The mobility of a protein molecule, modeled with atomic detail, in a bath of explicit water molecules can be studied, although at a very high computational cost, by molecular dynamics (MD). This approach is far from the hydrodynamic approach that we are considering here. However, it is worth mentioning that the MD approach has been explored and applied to lysozyme and BPTI by Smith and van Gunsteren (1994). The conclusions from that study may be a useful complement to those from hydrodynamic approaches.

Ellipsoidal models

The description of the hydrodynamics of rigid, globular macromolecules as ellipsoids is a classical approach that is well described in textbooks (Tanford, 1961; van Holde, 1971), in which these particles are treated as prolate or oblate revolution ellipsoids. The development of quasianalytical treatments for triaxial ellipsoids with three unequal axes (Harding, 1989, 1995) allows a more precise specification of size and shape. The lengths of the three axes can be related to the atomic structure of the macromolecule (Harding et al., 1999). An obvious advantage of this approach is that it requires less computational effort than bead methods. Another important aspect is that the hydrodynamic properties of general ellipsoids can be exactly described, although this is not a great advantage, because currently the errors resulting from the bead approximation are negligible in most cases. However, the details of the molecular structure are blurred in the smooth, entirely convex ellipsoidal shape. For completeness, we have mentioned this kind of model, but we do not go into further detail because this methodology is beyond the scope of the present study.

	MODELS AND METHODS

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

The information needed to construct the hydrodynamic model that represents the atomic-level, three-dimensional, (supposedly) rigid structure of the macromolecule is a list of the atomic coordinates obtained, for instance, from a Protein Data Bank (PDB) file (Abola et al., 1987). The steps that are followed for model building are schematically illustrated in Fig. 1. First all of the nonhydrogen atoms are represented by identical spheres of radius a. The differences in size between carbon, nitrogen, and oxygen (and the differences produced by the different number of implicit hydrogens) are averaged out for simplicity. Neighbor spheres representing covalently bonded atoms would be nearly tangent if the radii of the spheres were taken as the average covalent radii of the atoms, ~0.7 Å (Fig. 1 A). In this way the resulting model would be practically free of overlap, which is an undesirable aspect in bead modeling. However, it is well known that covalent radii underestimate the effective molecular volume, which is accounted for more adequately by van der Waals radii of the atoms. The typical values of van der Waals radii (including implicit hydrogens) (Bondi, 1964, 1968) are in the range of 1.5-2.0 Å; a typical (average) value would be 1.8 Å, which could be taken for the radius of the atomic spherical elements in the model. While this atomic size gives an adequate description of the molecular size in vacuo, it is expected that this size may be insufficient, particularly for biopolymers in solution, because of hydration effects. Thus larger values of a are to be expected, and their difference from what we consider the minimum value would be attributed to hydration. Indeed, the excess over this minimum value would be equivalent of the thickness of the hydration layer (this aspect is discussed later on). A large value of a also has the benefit of filling small gaps or pockets in the interior of the protein that do not have any influence in the hydrodynamic behavior (see Fig. 1 B). At this stage, we have what we call the primary hydrodynamic particle (PHP), which is a compact cluster of overlapping spheres (Fig. 1 C). This gives a realistic representation of the size and detailed shape of the hydrated protein as a hydrodynamic particle. The radius of the atomic elements (AER), a, will be regarded in principle as a floating, adjustable parameter, for which the proper value will be set a posteriori.

View larger version (35K): [in this window] [in a new window]   FIGURE 1   Schematic representation of the construction of the hydrodynamic model, in a two-dimensional analog. (A) The model from atomic coordinates, with spheres representing atoms with their covalent radii. (B) The spheres are expanded up to some value, a, of their radii. (C) The resulting primary hydrodynamic particle (PHP). (D and E) The filling model, superimposed on the PHP, and alone. (F) The rough shell model.

In the second stage of the modeling procedure, we seek an adequate description for the calculation of hydrodynamic properties of the PHP. As the PHP is indeed a bead model (of the BPA type), it could be submitted directly to HYDRO (García de la Torre et al., 1994b), but this would create the above-mentioned difficulties related to bead overlap and volume corrections. From our previous work (Carrasco and García de la Torre, 1999a) we know that a shell model is an adequate description of the hydrodynamics of an arbitrarily shaped particle. Then we apply the computer tools that we have developed for shell-modeling a general particle to the PHP. We first construct the filling model (Fig. 1, D-E), in which the particle is filled with beads arranged in the most closely packed, hexagonal lattice. From this model we shall evaluate the volume and the radius of gyration of the particle, R_g, of the macromolecule. Then all of the beads that are internal, in the sense of being completely surrounded by a number of beads equal to the coordination number of the lattice (12 in our case), are removed. Thus we are left with what we call the rough shell model (Carrasco and García de la Torre, 1999a) (Fig. 1 F), which represents the surface of the particle (where hydrodynamic forces act) as a shell of small beads of radius . For this model, we use HYDRO (García de la Torre et al., 1994b) to calculate hydrodynamic properties, such as the translational diffusion coefficient D_t, the rotational diffusion coefficient D_r, and the intrinsic viscosity []. The discontinuities in the model, because of the discrete size of the beads and the roughness arising from their geometrical, lattice arrangement, are eliminated by extrapolation of values calculated for various 's to the shell model limit, corresponding to  = 0.

For this hydrodynamic calculation we employ the latest version of HYDRO (hydrod11.for), in which the matrix inversion subroutine employed in older versions has been replaced with subroutines from LAPACK and BLAS public-domain libraries (http://www.netlib.org/lapack), where we have found highly efficient subroutines that are applicable when the 3N × 3N supermatrix is not only symmetrical but positive definite. This is only warranted for nonoverlapping beads, either equal or unequal, the hydrodynamic interaction of which is described using the Rotne-Prager-Yamakawa or the Garcia de la Torre-Bloomfield interaction tensor, respectively, and which is described, for equal overlapping beads, with the Rotne-Prager tensor (Rotne and Prager, 1969; Yamakawa, 1970; García de la Torre and Bloomfield, 1977b). The ad hoc procedure that has been proposed (Zipper and Durchschlag, 1997; Carrasco et al., 1999) for unequal overlapping beads does not fulfill the condition of positive definiteness. Another noteworthy aspect of the calculations is that, as demonstrated in our previous work (Carrasco and García de la Torre, 1999a), the volume correction does not have to be included in shell modeling computations; its inclusion is irrelevant because its contribution vanishes in the shell model limit, but it adds some inclination to the linear extrapolations.

As the number of frictional elements is very high for small  (up to N  3000 beads in some cases), we considered the need to do the calculations with double precision. However, we found that the results with simple and double precision coincide practically for translational and rotational coefficients; calculation of the intrinsic viscosity results requires double precision only for more than 1000 beads. As the extrapolations can be estimated from results with a moderately small N (see below), the computationally expensive double precision is not really necessary.

For the construction of the filling model and the shell model we employ the general subroutines described elsewhere. This subroutines, along with HYDRO, and another home-written program that extracts the atomic coordinates from the PDB file, have been collected in a single piece of software, HYDROPRO. This program accepts the PDB file and the simple supplementary data needed by HYDRO (temperature, solvent viscosity, etc). For a given value of the AER, a, and for a set of user-supplied values of the radius of the shell beads, , the calculations are carried out. The shell model extrapolations are included within HYDROPRO, so that the user obtains directly the final values of the solution properties.

	RESULTS

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

The above methods have been applied to a set of 13 proteins that are listed in Table 2. The range of molecular masses is 6-230 kDa.

In the following, lysozyme is chosen as an example to display the results of the calculations for each protein. Indeed, lysozyme is the protein most frequently considered in the previous works that we reviewed above. Images of the PHP and the SHE model of lysozyme are given in Fig. 2. The shell extrapolations for the three hydrodynamic properties are illustrated in Fig. 3 (the radius of gyration is practically independent of ). The values of the AER for which we repeat the calculation are a = 2, 3, 4 Å, and, for larger proteins, up to 5 Å. Computed results for lysozyme are plotted versus a in Fig. 4. Experimental values in Table 2 were taken mainly from data compilations. As literature citations, we usually mention the source from which we obtained the value, which may not the original publication. When two or more references are given, they may correspond to different values, of which we will take the mean.

View larger version (90K): [in this window] [in a new window]   FIGURE 2   (A) A bead-per-atom (BPA) model of lysozyme, which we take as the primary hydrodynamic particle (PHP) that represents this protein. The atomic element radius (AER) is a = 3 Å. (B) A shell model (SHE), derived from the PHP, used for hydrodynamic calculations. The radius of the small beads in this case is  = 0.8 Å.

View larger version (22K): [in this window] [in a new window]   FIGURE 3   Shell-model extrapolations to  = 0 of the hydrodynamic properties of lysozyme with a = 3 Å. The shell-model limit is obtained from the intercept of a least-squares-fitting straight line.

View larger version (21K): [in this window] [in a new window]   FIGURE 4   Determination of the experimental values ag, at, ar, and a of lysozyme by extrapolations in property versus a plots and estimation of the computed properties for the final value of a.

For each property, the value of the AER, a, for which the calculated result would match the experimental result is found by Newtonian interpolation, as illustrated in Fig. 4. The resulting values, a_t, a_r, and a, obtained, respectively, from D_t, D_r, and [] (3.0, 3.1, and 1.9 Å for lysozyme), are condensed into a single value, which is their average, a_h (2.7 Å for lysozyme), and can be regarded as a best fit for the set of hydrodynamic properties. Similarly, a value a_g is obtained separately from the radius of gyration (2.5 Å for lysozyme). The intention is to compare its result with that of a_h, to analyze possible differences in the effects of hydration on scattering and hydrodynamics (this will be discussed later); for lysozyme we note that a_g and a_h are nearly coincident. All of the available values of the AER, including a_g, can be combined into a single value a, which would be the average of all of them (2.7 Å for lysozyme). The results for each protein are listed in Table 3. With the values of a, we can again calculate by interpolation the values of the solution properties of the proteins, whose percentage deviations from the experimental values are listed in Table 3. On average, the deviations are ~2% for D_t and R_g and ~5% for D_r and []. For the set of 13 proteins, the average of the values of the AER is a = 3.2 Å, although the individual values show a noticeable variability, which, as discussed below, can be meaningful. Despite this variability, it is interesting to test the possibility of describing all of the proteins with a single AER of 3.2 Å. We again calculate the properties with this value of a and obtain the percentage deviation from the experimental results (data not shown). The average over all of the proteins of the absolute values of the percentage deviations is ~4% for D_t and R_g, 8% for D_r, and 16% for [].

As discussed in our previous work (Carrasco and García de la Torre, 1999a), the drawback of the shell-modeling strategy is computing time. The number, N, of beads with radius  needed to cover a given surface is proportional to ². The cpu time needed for the HYDRO calculation is proportional to N³ and, therefore, to ⁶. Thus the computational cost of, say, the two latest points in the extrapolations (Figs. 3) may be as much as that of all the others. Table 4 gives some cpu times for typical cases.

All of the results reported so far have been obtained from a full shell model calculation; it seemed necessary to have the maximum accuracy in the calculations for a relevant comparison with experimental data. Afterward we investigated computational shortcuts that would reduce the cpu time without an important loss of accuracy. Fortunately, the slope of the extrapolations is quite small, and the numerical values for not too small 's are quite close to the extrapolated ones. Then, if the slope b_p in the linear extrapolation of any property, p() = p(0) + b_p, could in some way be estimated, then the shell limit value could be obtained from a single datum, corresponding to a given , as p(0) = p()  b_p. Even if b_p is not very accurately approximated, it is still useful, because the b_p term is small. In our preceding study of shell models of geometric bodies, we found that when all of the quantities were expressed in a conveniently reduced (nondimensional) form, the reduced slopes b_p* in the linear variation p*() = p*(0) + b_p** take a numerical value that depends on the shape of the object but not on its size. As the overall shapes of the various globular proteins are not much different, we can expect that a single value of b_p* can be used for any of them. Reasoning by way of the reduced quantities, we arrive at the conclusion that the slopes for R_g, D_t, D_r, and [] could be put in the forms b_g = q_g, b_t = q_tD_t², b_r = q_rD_r^4/3, and b = qMD_t², where q_g, q_t, q_r, and q are numerical constants. From the full extrapolations that we have carried out for each of the 13 proteins, we notice that despite the wide range of molecular weights, the numerical values of the constants are of the same order of magnitude and do not deviate too much from each other. Taking the mean over the 13 proteins and the various values of a, we arrive at the following values for the constants: q_g  0 (as mentioned above, the extrapolation of R_g is nearly horizontal), q_t = 9.7 × 10³, q_r = 1.11 × 10⁴, and q = 2.44 × 10⁶, where D_t is expressed in cm²/s, D_r is expressed in s¹, [] is given in cm³/g, and  is in Å.

With these values, the slopes for the estimated extrapolations can be obtained. We have checked the performance of the estimated extrapolations, comparing their results with those from the full extrapolations. We start from a single datum, obtained for a value of  (different for each protein) for which the number of beads in the shell is close to 1000, and therefore the cpu time is a few minutes (see Table 4). Then the extrapolated values of the properties are estimated using the above values of q_g, q_t, q_r, and q. Table 5 shows how the estimated values differ very little from the true ones. The differences are usually smaller than the errors that we can expect in the experimental value of the properties.

Thus the estimated extrapolation is a computational shortcut that avoids the computationally expensive computations for a very small bead radius and places the strategy that we propose for the calculation of solution properties of globular proteins within the reach of a personal computer.

	DISCUSSION

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

The comparison of calculated and experimental results requires some comments about the latter. From our literature review, we have gained the impression that, in many instances, the experimental data may be, say, 40 years (or more) old. The information on rotational diffusion coefficients is scarce, although the inclusion of this property in an analysis like the present one is most valuable, owing to its sensitivity to macromolecular structure.

There are also difficulties in the determination of intrinsic viscosities. For globular proteins with a shape that does not deviate much from spherical, [] is lower than for any other macromolecules. The theoretical lower limit for the intrinsic viscosity of a spherical particle obeying stick hydrodynamic boundary conditions is the Einstein [] = 2.5 cm³/g value. (As kindly pointed out by a referee, the limit is different for slip boundary conditions, but there is a general consensus in that proteins in water are better described by stick rather than slip behavior.) The experimental values for the various proteins are slightly above the Einstein [] = 2.5 cm³/g limit for an anhydrous sphere, and the difference is due in part to the (uncertain) effect of hydration. Being so small, the experimental determination of [] is prone to errors, which may particularly affect most of the values used here, which are certainly old (the modern refinements in the viscometric techniques, as described by Harding (1997), may improve this situation in the future). Despite these difficulties, the inclusion in our analysis of the intrinsic viscosity, [], is a must, because it is a classical property that represents the important aspect of hydrodynamics in shear flows. A reliable calculation of the intrinsic viscosity is also useful for discarding anomalous, too small values of a, for which [] may fall below the Einstein limit. In this context, it is remarkable that, even if the fit for [] is not entirely satisfactory, we have greatly improved the predictions of some previous works (Byron, 1997; Zipper and Durchschlag, 1997), which employed hydrodynamic models and procedures that were inadequate for this property.

These comments suggest the potential interest of new experimental measurements of the solution properties of the most common globular proteins, which would be proper references for testing theories and for a better knowledge of protein hydrodynamics. Meanwhile, we have to proceed with the available data.

A clear conclusion from our work is that up to four solution properties, R_g, D_t, D_r, [], of a given globular protein can be predicted with expected errors of ~2%, ~2%, ~5%, and ~5%, respectively, with a single, adjusted value of the radius of the atomic elements, a, in the primary hydrodynamic model. The predictive capability of our methodology is similar, and in some cases appreciably better, than that of the procedures proposed by other workers. The use in our method of an adjustable parameter is in common with the other studies, where quantities such as the thickness of the hydration layer (Allison and Tran, 1995; Zhou, 1995) or the number of water molecules within it (Venable and Pastor, 1988) or the hydration expansion factor (Byron, 1997; Zipper and Durchschlag, 1997, 1998) have been adjusted to fit the experimental properties.

An obvious question is whether the values of the adjustable parameter (whatever it is) should be constant or variable for various globular proteins. Unlike most previous studies, we can address this point because we have considered a large set of proteins. As commented above, the values of our AER vary in the range of a = 2-5 Å, with an average of 3.3 Å. An evident source of a good part of the observed scatter in a is the quality of the experimental data, which has already been discussed. It should also be remarked that this variability must be judged in terms of the changes that it may introduce in the overall particle size; thus changing a from 2 to 4 Å changes the (equivalent) hydrodynamic radius of lysozyme from 19 to 21 Å, which is an increase of only 10%. The importance of a similar change for a larger protein is even smaller. However, such changes produce remarkable changes in the calculated values of, mostly, the rotational coefficient and the viscosity; in the case of lysozyme, such a change modifies these two properties in about 30%.

Apart from the uncertainty in the experimental data and the different sensitivities of properties to modeling details, we think that the variability in the AER has a physical, molecular origin in protein hydration. Apart from the ambiguity in defining hydration, it is well known that, for a given definition, the amount of hydration varies from protein to protein. For instance, in the classical description of proteins as ellipsoidal particles (or generally, in models where hydration is treated by uniform expansion), the amount of hydration expressed as , grams of water per gram of protein, presents values in a wide range (say,  = 0.2-0.5 g/g). In our model, as described above, the effect of hydration is a contribution to the AER, which is expected to be larger than the van der Waals (average) radius a_vdW  1.8. The difference, l = a  a_vdW, could be regarded as the thickness of the hydration layer that coats the protein surface. Then, according to our results, a typical value of l would be 1.5 Å, corresponding to the average AER of 3.3 Å, although individual values may vary appreciably, reaching 3 Å in some instances. Workers who employ the FEL strategy (which has in common with our procedure the focus on the macromolecular surface) also reach conflicting conclusions about the hydration layer; while Zhou (1995) reports a thickness of scarcely 1 Å, Allison and Tran (1995) obtain larger estimates of up to 3 Å.

In our results we notice some tendency of l to increase with the size of the proteins, as seen for the proteins with molecular mass, M, above 100 kDa (see Table 3). Thus although the evidence is not strong enough, it is possible to speculate about a possible increase of l with M. Indeed, a constancy of the thickness would not be fully compatible with a constancy of the  parameter.

Another aspect related to hydration is its effect on the particle size determined by scattering, measured by R_g. Whether the scattering effective hydration is the same or smaller than hydrodynamically effective hydration is an unsolved question. Looking at the results in Table 3, we see that the AER obtained from R_g is very similar to that obtained from hydrodynamic properties, and this is true not only for the average values, but also for the individual values for most of the 13 proteins. It seems, therefore, that the hydration effects are similar, and the same model parameter serves for hydrodynamics as well as for scattering properties.

The question of which of the measures of hydration (i.e., l, , etc.) is most adequate (less variable) for globular proteins cannot be answered yet. It should be pointed out that the definition or at least the evaluation of any hydration-measuring parameter relies on the specificities of a hydrodynamic model. We hope that the modeling strategy proposed in this work will be useful in the future, in studies including many proteins with reliable data, for gaining a better understanding of the hydration and hydrodynamics of globular proteins.

	COMPUTER PROGRAMS

TOP ABSTRACT INTRODUCTION REVIEW OF MODELS MODELS AND METHODS RESULTS DISCUSSION COMPUTER PROGRAMS REFERENCES

The bead and shell models displayed in this work were visualized using the public-domain POLYRAY raytracing software (http://ftp.tu-clausthal.de/pub/TEXT/mirror/povray/polyray). The computer program HYDROPRO, which has been used in all of the calculations in this project, has been built from the modules described in our previous work (Carrasco and García de la Torre, 1999a) for shell models of arbitrary particles and includes a subroutine to extract the atomic coordinates from a PDB file. HYDROPRO is in the public domain and will be available for downloading, both as FORTRAN source code as well as in executable forms for several plataforms, at our web site (http://leonardo. fcu.um.es/macromol).