INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

To combat the potentially lethal effects that occur with the freezing of water, organisms have evolved various protective strategies, one of which is to produce macromolecular antifreezes, including antifreeze proteins (AFPs) and antifreeze glycoproteins (AFGPs) that have long been known to inhibit ice formation in physiological conditions. The best characterized of these AFPs come from polar fish and from several different species of insects. Although structurally diverse, these proteins have in common the ability to bind to ice and inhibit its growth (Yeh and Feeney, 1996; Davies and Sykes, 1997). AFPs lower the freezing point below the melting point in a noncolligative manner, a property termed thermal hysteresis, and are believed to act by adsorbing to ice surfaces, causing ice crystals to grow with a submicroscopically curved surface (Raymond and DeVries, 1977; Knight et al., 1991). The addition of water to the curved ice lattice is energetically unfavorable, and additional ice growth is retarded.

Various structurally distinct AFPs have evolved independently (Davies and Sykes, 1997). AFGPs and type I AFPs are rod-like structures with simple repeating sequences (Yang et al., 1988; Sicheri and Yang, 1995). Type II AFPs are lectin-like globular AFPs divided into two main classes, Ca²⁺ dependent and Ca²⁺ independent. Type III AFPs are also globular proteins (Jia et al., 1996; Sönnichsen et al., 1996; Yang et al., 1998), whereas type IV AFPs from longhorn sculpin are predicted to have an -helix bundle structure (Deng et al., 1997). Recent insect AFP structures for spruce budworm AFP (Graether et al., 2000) and Tenebrio molitor AFP (Liou et al., 2000) have been reported, both of which are -helices, though of opposite handedness. Overall, AFPs share very few features in common. Some AFPs have repeating sequences whereas others do not. Some are highly hydrophobic whereas others are more hydrophilic. Indeed, the only few obvious common characteristics among all AFPs are their relatively low molecular mass (which is usually 25 kD) and, for those AFPs with structure determined thus far, a relatively flat and proportionally large surface area for ice binding (Jia and Davies, 2002).

Historically, our understanding of AFP activity has been derived almost exclusively from extensive studies of type I AFP. Initial structural analysis and mutagenesis experiments had suggested that the more hydrophilic face of the AFP is involved in ice binding and, consequently, that hydrogen bonding between ice and these hydrophilic residues, particularly the regularly spaced Thr residues, was the dominant force anchoring AFPs to ice (Yang et al., 1988). Additional experiments involving the mutagenesis of the key Thr residues to Ser and Val, however, have cast doubt on the importance of hydrogen bonding for AFP-ice interactions (Chao et al., 1997; Haymet et al., 1998). Recently, a report of mutagenesis experiments done on the highly conserved Ala-rich face has suggested that in fact it is this hydrophobic face of type I AFP that interacts with ice (Baardsnes et al., 1999). Similarly, the structure of the globular type III AFP has generated much speculation as to its binding mechanism. Structural characterization and mutagenesis results of this AFP have identified one particularly flat, amphipathic surface of this protein as the likely ice-binding face (Chao et al., 1994; Jia et al., 1996). Several novel computational studies have been done in an effort to characterize the possible binding mechanism of this face. These include one by Chen and Jia (1999) that found the hydrophilic forces were not the dominant AFP-ice interactive forces and others by Yang et al. (1998) that found the flatness of this face primarily responsible for ice binding and Graether et al. (1999) that used a neural network to conclude that this plane's hydrophobicity was ultimately responsible for AFP activity. It should be pointed out that the results by Chen and Jia (1999) demonstrated the usefulness of random ice in the computational studies involving AFP. In the case of the insect AFPs, the highly ordered Thr residues along one face of the -helical T. molitor AFP have been found to match ice-oxygen orientations on the prism plane of an ice crystal, making this surface the prime candidate for ice binding (Liou et al., 2000). Furthermore, mutagenesis experiments performed on spruce budworm AFP have clearly demonstrated that the key Thr residues in the T-X-T motif (X stands for any residue) play a critical role in thermal hysteresis activity (Graether et al., 2000).

There have been numerous molecular modeling studies on the interactions between AFPs and ice (for a review of this work, see Madura et al., 2000). What emerges from this brief portrait of AFP functional studies is that the specific forces involved in AFP-ice interactions are not well understood but nonetheless that superior ice-binding ability is generally held to be the essential characteristic that both distinguishes AFPs from non-AFPs and is responsible for differences observed in AFP activity levels. Indeed, several theories have been proposed that link ice binding to AFP activity, including the adsorption/inhibition model (Raymond and DeVries, 1977; Knight and DeVries, 1994), the reversible adsorption model (Feeney et al., 1986), and the anisotropic interfacial energies model (Wilson, 1994).

Type II AFPs are 14-24-kD homologs of the carbohydrate-recognition domain of Ca²⁺-dependent lectins (Ewart et al., 1992). In addition, type II AFPs also have similarities with pancreatic stone inhibitor proteins (Bertrand et al., 1996; Patard et al., 1996; Gronwald et al., 1998). These proteins have been categorized into two classes based on their Ca²⁺ requirement. The first class includes type II AFPs from Atlantic herring, which have one Ca²⁺-binding site; both the activity and conformation of these AFPs are Ca²⁺ dependent (Ewart et al., 1996). Site-directed mutagenesis has revealed that the ice-binding site of this AFP corresponds to the carbohydrate-binding site of C-type lectin (Ewart et al., 1998). The second class of type II AFPs includes the Ca²⁺-independent AFP from sea raven. This AFP is the largest globular AFP whose three-dimensional structure is known. Although its NMR structure is of low to medium resolution (Gronwald et al., 1998), the structure has revealed that the global fold of this AFP is similar to C-type lectins and pancreatic stone proteins, despite having only ~20% sequence identity (Sönnichsen et al., 1995). The overall structure is characterized by five conserved disulfide bridges, two -helices, one -strand structure, and extensive loop regions. Mutagenesis experiments have indicated that the ice-binding site of sea raven type II AFP is distinct from the carbohydrate-binding site of the homologous C-type lectin (Loewen et al., 1998). In contrast to the Ca²⁺-dependent herring type II AFPs, where the Ca²⁺-binding site is known to be related to the ice-binding site, the ice-binding site of the Ca²⁺-independent type II AFP from sea raven has yet to be established, despite its NMR structure (Gronwald et al., 1998) and extensive mutagenesis (Loewen et al., 1998). And unlike the extensive AFP-ice interaction studies for type I AFP, there has been only one study on type II AFP-ice interactions reported (Wierzbicki et al., 1997). This modeling study has shown that the fold of type II AFP could facilitate a stereospecific interaction with ice, a result that is in general agreement with ice-etching experiments that have suggested that the specific ice-binding plane of type II AFP is {111} (Cheng and DeVries, 1991).

Traditionally, because of size limitations, most computational investigations of macromolecules use molecular mechanical methods only, typically including energy minimization and molecular dynamics simulation. Though molecular mechanical methods aim to predict various molecular properties as deduced in analytical expressions from empirically designed molecular force fields, they ignore the electronic motions and calculate the energy of a system as a function of the nuclear positions only. Because no attention is paid explicitly to the electronic motions in an assembly of molecules, molecular mechanics cannot provide properties that depend upon the electronic distribution in a molecule. In contrast, quantum mechanics explicitly represents the electrons in a calculation, and so the use of this technique provides the possibility of deriving properties that depend upon the electronic distribution within assemblies of molecules, and in particular it allows for the investigation of essential chemical reactions and interactions between molecules. For example, very recently a mixed quantum/molecular mechanics approach was used to study an enzyme reaction mechanism (Dinner et al., 2001). Unfortunately, many of the problems that we would like to tackle in molecular modeling have been too large to be considered by quantum mechanical methods. However, with the advancement of computer technology, some large systems, which have been intractable even on the most sophisticated computers of the past, are gradually becoming feasible in semi-empirical molecular orbital calculations. The inaccuracy of the approximations inherent in this technique is offset to a degree by recourse to experimental data in defining the parameters of the method. Indeed, semi-empirical methods can sometimes be more accurate than some ab initio methods, which require much longer computation times.

In this paper, we have used two semi-empirical quantum mechanical methods, AM1 (Dewar et al., 1985) and PM3 (Stewart, 1989a,b), to systematically study a series of large systems involving type II AFP and ice lattices based on the initial configurations obtained by molecular mechanics.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Ice planes and type II AFP structure

Both an {111}ice plane and a random ice plane were built according to Chen et al. (1999) using Sybyl 6.4 (Tripos, St. Louis, MO). The ice planes have an approximate size of 55.9 × 56.9 × 7.4 Å³, larger than any cross-sectional diameter of type II AFP. Although water molecules are discretionarily located in the random ice slab, the volume was set so as to maintain the same density in the random ice slab as is found in regular ice (P6₃/mmc).

The NMR structure (Gronwald et al., 1998) of type II AFP from sea raven, the only experimental structure of type II AFP, was used in this study.

Choice of AFP surface sites and docking to ice

To systematically study surface sites of type II AFP, we selected 17 surface patches to cover the entire surface of the AFP regardless of whether they had been implicated in ice binding or not. Using the protein surface analysis program GRASP (Nicholls et al., 1993), it was clear that these surface patches overlap and together provided complete surface coverage. Because of the overlap, many residues were contained in multiple patches (Table 1).

Using each surface patch in turn, the protein was docked to both the {111}ice plane and our random ice plane using Sybyl 6.4. To ensure maximal AFP-ice contact, a central residue from each patch was roughly aligned with the center of the ice plane for each docking. After this, the protein was positioned to maximize the number of surface-exposed atoms on a given patch that could interact with the ice face. In the docking procedure, the shortest distances between O atoms of an ice face and surface atoms of the AFP were 2-4 Å, a distance appropriate for interactions such as hydrogen bonding and van der Waals force.

Optimization of the AFP-ice complex systems

Before energy evaluation by quantum mechanical methods, it was necessary to optimize all AFP-ice docking complexes to avoid unnecessary and computing-intensive complex energy computations. Each optimization contained three steps, similar to those previously reported (Chen et al., 1999). The AFP-ice complex was first subjected to energy minimization using the conjugate gradient method with the AMBER force field (Weiner et al., 1984) parameterized in 1991 (see http://www.amber.ucsf.edu/amber/ff91.html). During the energy-minimization process, all AFP atoms were allowed to move freely, with the exception of the oxygen atoms in the ice. After energy minimization, molecular dynamics simulations were carried out at the constant temperature of 273 K without constraints. A final round of energy minimization completed the procedure. For more details, see our previous report (Chen et al., 1999). These procedures were applied to all 17 surface patches, each one docked to both ice planes, and the final interaction energy for each docking was calculated. The program Discover 2.9 in the package Insight II 95.0 (Biosym/MSI Co., San Diego, CA) was used for the calculations.

Semi-empirical quantum mechanical calculations

Because each system contains a large number of atoms (a total of 3632 atoms with 1559 nonhydrogen atoms), we used the AM1 and PM3 semi-empirical molecular orbital methods to evaluate the energy and other properties of the AFP on its own, both ice planes on their own, and the 34 AFP-ice complex systems. The calculations were applied to the systems that resulted after molecular mechanics optimization as described above. Through comparisons of the energies for the AFP and the ice planes in each of the 34 AFP-ice complexes, the interaction energies between AFP and ice were obtained.

	RESULTS AND DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

AFP model from optimization

Molecular dynamics simulations can result in large changes in protein structure, particularly when such simulations are undertaken without constraints, as was the case in this study. To estimate structural change, all AFP structures obtained from the AFP-ice complexes after optimization were compared with the original AFP NMR structure. All atoms of each optimization-derived model were superimposed on the NMR structure by means of least-squares fitting. The results are tabulated in Table 2, showing that root mean square deviations range from 2.5 to 2.9 Å, which are reasonable. The N-terminal region has the largest difference, which is not unexpected because it is very likely to be intrinsically flexible. The optimized patch 1 found in the AFP-{111}complex is depicted in Fig. 1.

View larger version (30K): [in this window] [in a new window]   FIGURE 1   The optimized patch 1 system in the AFP-{111} ice series. Type II AFP is shown as a green ribbon. The O atoms and H atoms in the ice are shown in red and blue, respectively.

Molecular mechanical calculations

Though the main purpose of using molecular mechanical methods in this study was to provide optimized systems for later evaluation using quantum mechanics, we also calculated interaction energy by molecular mechanical methods after optimization. Table 3 lists the interaction energy of the random ice series. Patch 1 was found to have the most favorable interaction energy. More importantly, there is an obvious energy discrimination between patch 1 and the next best patch whose interaction energy is only 67.5% of that of patch 1, corresponding to a decrease of 45.65 kcal/mol. The interaction energies of other patches decrease very quickly, with the last patch having only 7.4% of the interaction energy of patch 1. Furthermore, we have divided the interaction energy into van der Waals interaction energy (I-vdW), hydrogen bond interaction energy (I-HB), and coulomb interaction energy (I-C) and calculated them respectively. Consistent with many recent findings, the results listed in Table 3 show that I-HB is not the dominant AFP-ice interaction energy (Chao et al., 1997; Chen et al., 1999; Haymet et al., 1998). Using the regular {111} ice plane as an interacting target (see Table 4), the results again indicate that, despite some reordering among the patches with respect to their interaction energies, patch 1 remains the best ice-binding patch, and again by a significant margin. Note, however, that because the distributions of water molecules in the random and regular ice slabs are different, it is impossible to compare the interaction energies directly between the two series.

View this table: [in this window] [in a new window]   TABLE 4   Interaction energy between AFP-II and {111} ice plane using molecular mechanics method

Quantum mechanical calculations

Molecular mechanical calculations deal with interactions and other properties from a purely classical mechanical viewpoint. Therefore, they cannot be used to investigate bond formation and breaking, or a system in which electronic delocalization or molecular orbital interactions plays a major role in determining geometry and/or properties of the given system. In contrast, quantum mechanical methods directly investigate electronic movements and deal with the detailed features of interactions between molecules, thereby giving rise to insights into the nature of molecular properties. The two calculation methods we used (AM1 and PM3) are suitable for studying systems that contain hydrogen-bonding interactions. For example, in our previous study (Chen et al., 2002) we demonstrated that the AM1 method might actually be as good as the higher-level B3LYP method (Becke, 1993; Lee et al., 1988). Another excellent recent example was the computational study of the reaction mechanism of uracil-DNA glycosylase by the AM1 method (Dinner et al., 2001).

The interaction energies in the{111} and the random ice series were calculated using the AM1 and PM3 semi-empirical molecular orbital methods; the results are shown in Table 5. For both the regular and the random ice series, patch 1 has the most favorable interaction energy among all patches. Although there are some discrepancies between the interaction energies as calculated by the molecular mechanical method and the two quantum mechanical methods, there is certainly a common trend in these energies among all three methods (Fig. 2).

View larger version (19K): [in this window] [in a new window]   FIGURE 2   Interaction energy curve. (A) The random ice series; (B) The regular ice series. It is evident that patch 1 consistently has the best interaction energy. In general, all three different methods generate a similar trend in terms of energy change.

Weak molecular orbital overlap interactions

In linear combination of atomic orbitals-molecular orbital approximation, the molecular-orbital function in a molecule can be expressed as a linear combination of atomic orbitals. Here the localized molecular orbital (LMO) methods (Von Niessen, 1972) were used. If the square of the combination coefficient of the wave-function of an atom is larger than 0.0002 (i.e., the combination coefficient of the wave-function on an atom is larger than 0.01414) then it is taken as the criterion of the essential contribution toward the LMO in the system. Those wave-functions were omitted for some atoms whose squares of combination coefficients are less than 0.0002. In light of these conditions, in the AFP-ice complex systems the LMOs may be divided into three types, that is, AFP alone, ice alone, and AFP-ice complex. The LMOs composed of combinations of wave-functions of atoms in both AFP and water molecules of the ice lattice were considered as interactive LMOs. In fact, the combination coefficients of atomic orbitals of some atoms in interaction area are usually ~0.03-0.10. Therefore, the interactive LMOs represent weak orbital overlap interactions, which are, however, important in helping us understand the mode of AFP action. Tables 6 and 7 list the number of the occupied interactive LMOs (N), the sum of their orbital energies (E), and the proportion of these energies in the total occupied energies of the whole 5024 LMOs in the system. These data can reflect the size and degree of the orbital interaction between AFP and ice. Patch 1 has the most interactive LMO numbers and largest interactive LMO energy in both the random and regular ice slabs regardless of the semi-empirical molecular orbital method used. For example, the results derived from the PM3 calculations show that the number of interactive orbital levels of patch 1 is 121 in the {111}ice series, the largest value among all patches. However, the occupied orbital energy of patch 1 is 4456.0 eV, which is substantially lower than 1173.5 eV and 3594.1 eV of patch 9 and patch 17, respectively. From quantum chemistry theory, the more the orbital overlap, the better the interaction between molecules. Given this, patch 1 again appears to provide the strongest docking surface for ice.

View this table: [in this window] [in a new window]   TABLE 6   Number of occupied interactive LMOs (N), sum of their orbital energies (E), and its proportion in the {111} ice series

Bond order

Semi-empirical molecular orbital methods can offer a measure of the strength of the bonds in a molecule, which is known as bond order. Bond order change actually reflects a change of (strong and weak) bond charge density. Therefore, it will impact on molecular conformation. Although it is very difficult to investigate the bond order of every atom in such a large system, it is nevertheless possible to evaluate the net difference of the whole system before and after complex formation. In the analysis of bond orders, we found an important observation in all systems investigated: the bond order in ice lattices decreased significantly upon AFP binding. For example, in the patch 1 with the {111}ice system, the bond order of the ice lattice decreased by 0.296 and by 0.483 using AM1 and PM3, respectively. Usually the bond order of hydrogen bonds within ice lattices is in the range of 0.013-0.015 (AM1) and 0.014-0.024 (PM3). Therefore, this analysis suggests that a minimum of 20 hydrogen bonds in the ice lattice would be broken upon interaction with type II AFP. It is possible that the disruption of the internal hydrogen bonds within the ice lattice could result in a partial melting of the surface of an ice crystal and prevent it from further growth. Perhaps more significantly is the possibility of the partially weakened/distorted ice surface layers that would grow and result in a new ice crystal with a slightly twisted lattice as compared with the original lattice, thereby causing ice surface poisoning that would require extra energy to overcome before further growth can occur.

	CONCLUSIONS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Using 17 surface patches of type II AFP and both random and {111} ice planes as interacting targets, we have determined that patch 1 of fish type II AFP is the most favorable ice-binding site. The results obtained by molecular mechanics and semi-empirical quantum mechanics are consistent. This surface patch contains 19 residues (K37, L38, Y20, E22, Y21, I19, L57, T56, F53, M127, T128, F129, R17, C7, N6, P5, G10, Q1, and W11). For the first time, the detailed information afforded by semi-empirical quantum mechanical calculation provides new insights into the interactions between type II AFP and ice, including the weak orbital overlay phenomenon and the change of bond orders of the ice. The results obtained from two semi-empirical approaches are consistent. They both imply that there are weak orbital interactions between type II AFP and ice. Importantly, weakened bond order and the disruption of otherwise regular hydrogen bonds within the ice lattice would result. We suggest that this finding may be closely related to how AFPs perform their antifreeze function.

In this paper we have provided a semiquantitative picture on the antifreeze-ice interaction. For such a large system, current computing power would not allow the use of better and more sophisticated quantum mechanical methods. Using semi-empirical quantum mechanical methods we have been able to make a number of novel observations, particularly that of the weakening of the internal hydrogen bonding interaction with the ice lattice upon AFP binding, which would be otherwise impossible to obtain with the classical molecular mechanical methods.