INTRODUCTION

Top Abstract Introduction Methods Results Discussion Conclusion References

It is now widely accepted that the ability to specifically bind a particular ligand is determined by the three-dimensional shape of a protein (e.g., Creighton, 1994). The fatty acid-binding proteins (FABPs) have essentially similar backbone structures (Sacchettini and Gordon, 1993; Banaszak et al., 1994; LaLonde et al., 1994a), yet ligand specificity varies from one family member to another (Richieri et al., 1994, 1995, 1996). The detailed mechanism of this discrimination is currently not understood. The availability of high-resolution x-ray structures (Banaszak et al., 1994) and excellent thermodynamic data (Richieri et al., 1994, 1995, 1996) for several members of this family has begun to provide a basis for understanding selectivity. The use of molecular dynamics simulations may extend this initial understanding to the atomic level.

A molecular knowledge of hydrophobic ligand binding will illuminate several outstanding research problems beyond the present system. For example, the ability to identify the effect of a particular set of residues on binding affinity would make possible the rational design of future drug compounds (see Ajay and Murcko, 1995, for a recent review). Moreover, because very little is known regarding the importance of the membrane environment for protein function (e.g., Merz and Roux, 1996; Gennis, 1989), a better molecular understanding of the nuances of protein-lipid interactions would aid structure prediction and structure-function connections for membrane proteins.

The relationship between tertiary structure and thermodynamic measurements is established by statistical mechanics. Computer simulations have been used to generate thermodynamic ensembles in simple liquid systems for many years (e.g., see Allen and Tildesley, 1987), but until recently the computational cost of adequately sampling protein systems was prohibitive (e.g., Brooks et al., 1988). Even with modern computational resources, making the full connection between structure and function remains difficult. Despite these limitations, computer models can elucidate details of molecular motion and provide insights into interesting sites for possible mutagenesis.

The application of computational methods to FABPs has been limited. One study used high-temperature dynamics to suggest possible routes for ligand release (Zanotti et al., 1994). Another assessed the conformational energy of the ligand in the x-ray structure of human muscle FABP (M-FABP) (Young et al., 1994). A third group varied the charge state of the system to examine electrostatic effects (Rich and Evans, 1996).

Initial FABP binding studies suggested micromolar binding affinities (e.g., Matarese and Bernlohr, 1988; Miller and Cistola, 1993; Maatman et al., 1994), which seemed inconsistent with the nanomolar free fatty acid concentration believed to be present in the cell (Richieri et al., 1992, 1993). Furthermore, little discrimination with respect to chain length or saturation state was observed. In contrast, recent binding measurements, utilizing acrylodan derivatized intestinal FABP (ADIFAB) (Richieri et al., 1994, 1995, 1996), suggest affinities in the nanomolar range, with strong chain length and saturation state selectivity. For example, they show a 20-40-fold greater binding affinity of stearic acid for M-FABP relative to adipocyte FABP (A-FABP) (Richieri et al., 1994).

The A-FABP structure has been determined to 1.6-Å resolution (Xu et al., 1992, 1993), and the M-FABP structure has been solved to 1.4-Å resolution (Young et al., 1994). The two structures have a backbone root mean square (RMS) difference of 0.7 Å. The conformation of the bound fatty acid in both structures is very similar for the headgroup region and along the alkane chain up to the C12 methylene group. In addition, the R106-R126-Y128 triad interacting with the headgroup is in nearly the same conformation in both structures. The pattern of ordered waters in the cavity is also similar. The 10 -strands characteristic of all FABPs are referred to with letters from A to J, and the two -helices are referred to as -I and -II (e.g., Banaszak et al., 1994).

Previous theoretical work explicitly assumed that the ligand carboxylate group is charged (e.g., Young et al., 1994; Zanotti et al., 1994; Rich and Evans, 1996). This assumption is also implicit in arguments for the importance of electrostatic effects in specific binding by I-FABP and CRBP-II (Jakoby et al., 1993). This reasoning is plausible in light of experimental evidence from intestinal FABP (I-FABP) suggesting that the pK of the carboxylate is similar to that in aqueous solution (Cistola et al., 1989). However, the headgroup structure of A-FABP and M-FABP is different from that of I-FABP and CRBP-II (Banaszak et al., 1994), and x-ray structures do not determine the positions of hydrogens. The current simulations considered both protonation states of the ligand.

The present calculations provide insight into the motional behavior of two fatty acid-binding proteins with similar ligands. These molecular details are not obvious from the x-ray structures and could not have been inferred from visual inspection alone. The results suggest similarities and key differences between the two systems. These may be related to the differences in binding affinity of stearic acid between the two proteins.

	METHODS

Top Abstract Introduction Methods Results Discussion Conclusion References

Molecular dynamics computer calculations require an initial conformation and a protocol for relaxing it to an equilibrium state, after which an ensemble of structures can be generated. Independent simulations were performed for two explicitly solvated FABP:stearic acid complexes. In each case the starting point was a high-resolution crystal structure. For M-FABP, this was the 1.4-Å resolution structure (Young et al., 1994). The 1.6-Å-resolution structure was used for A-FABP (Xu et al., 1992, 1993). The protocol was the same for both structures and so is summarized only once. Version C23f3 of the CHARMm program was used with a recently developed parameter set for protein:lipid systems (Schlenkrich et al., 1996).

Determination of the pK at the carboxylate group of the fatty acid was an essential part of these calculations. Initially, the charged (deprotonated) form was used. This was consistent with the assumption in many discussions of FABPs (e.g., Young et al., 1994; Zanotti et al., 1994; Rich and Evans, 1996). However, these trajectories diverged more than 5 Å C- RMS from the crystal structure within 50 ps after equilibration (described below). This prompted Poisson-Boltzmann calculations of the headgroup pK using DELPHI (Gilson et al., 1988). The calculation used CHARMm partial charges and x-ray structures. The calculated result was an estimated shift of pK at the headgroup by 12 pK units. This strongly suggested a neutral headgroup. When a neutral headgroup was used in the CHARMm calculations, the trajectory remained stable through 1 ns of dynamics. At no point was there any evidence of the instability first seen with the charged headgroup. Similar DELPHI calculations and long molecular dynamics trajectories on the intestinal FABP (I-FABP) (Tychko and Woolf, 1997) suggest that the I-FABP system should have a charged ligand. This is consistent with the only available experimental information regarding the titration of the headgroup (Cistola et al., 1989). Taken together, this evidence supports the conclusion that the region near the headgroup of ligands within A-FABP and M-FABP differs from that of I-FABP.

To date, the DELPHI and CHARMm calculations have explored only the charge state of the ligand and not all possible proton distributions in this region. For example, the carboxylate group, the bridging water, and R106 could dynamically share a proton that may have originated from R106. A full analysis of the possible combinations throughout the set of FABPs is planned. In particular, calculations can provide an estimate of the probability of finding the proton on R106, the bridging water, or the ligand.

A "water-droplet" model was used, solvating each system with a large spherical shell of waters. A restraining potential was applied to water oxygens to maintain the overall structure of the droplet. This potential is zero for oxygens more than 0.7 Å from the surface. A small well of 0.25 kcal/mol is present 27.0 Å from the center of the system. This approach is related to the stochastic boundary potentials first introduced by Brooks and Karplus (1983) and allows effective solvation with fewer waters than are required for full periodic boundary conditions. It is expected that the method will be especially effective for ligand interactions within a protein cavity far from the surface. The motions of peripheral residues might be influenced by the effective surface term, so their behavior was treated with less confidence in the analysis. Future work could use the minimum solvation approaches being developed by Beglov and Roux (e.g., Beglov and Roux, 1994).

The systems were constructed by first adding hydrogens to the heavy atoms of the x-ray structure (pdb files 1lif (A-FABP) and 1hmt (M-FABP)). The structures were then gently relaxed in the CHARMm potential through a series of steepest descent minimizations with decreasing harmonic restraints on all heavy atoms. At the end, the RMS deviation from the x-ray structure was 0.3 Å for C, and the gradient was 0.6 kcal/mol-Å. This was the starting point for solvation and equilibration. The FABP was surrounded by a large cube of preequilibrated water. Water with oxygens outside the 27.0-Å spherical boundary was deleted, along with any waters with oxygen atoms within 2.6 Å of any heavy FABP atom. This resulted in a system size of nearly 7700 atoms.

The system was then minimized and equilibrated. A series of steepest descent minimizations was performed, first with the FABP, x-ray water, and stearic acid atoms fixed, and then with decreasing harmonic restraints on those atoms. No restraints were placed on the bulk waters. This was followed by 2 ps of dynamics with harmonic restraints on the heavy atoms of the system. A further equilibration period of 50 ps was used before conformations were saved. The first part of this period used Langevin dynamics with a frictional coupling constant of 25 ps¹ and a temperature of 300 K. During the second part, velocities were rescaled whenever the temperature deviated from 300 K by more than 5 K. The update time for checking the windows was every 2.5 ps.

Trajectories were calculated for 1 ns. The energy and temperature remained stable throughout the simulation. The SHAKE algorithm was used, allowing a 2-fs step size during dynamics. Conformations were saved every 25 steps (0.050 ps) for later analysis. The nonbonded list was generated to 13 Å. Van der Waals (vdW) interactions were switched from 10 Å to 12 Å, and atom-based shifting was used for the electrostatic interactions. A combination of CHARMm routines and custom scripts was used in the analysis. The simulations typically required 2.5 h/ps on a dedicated R4400 SGI Indigo2 workstation. Thus each 1-ns trajectory represents more than 100 days of CPU time.

An importance scale, based on positional correlation and interaction energy with the ligand, was developed to help identify amino acids most involved in ligand binding. The interaction energies between FABP residues and the heavy atoms of the ligand were normalized, as were the zero-time covariances. A high interaction energy score indicated a strong, favorable interaction energy. A high covariance score meant either strongly correlated or anticorrelated motion. The importance value for each residue was calculated by multiplying the interaction energy score by the covariance score for each heavy atom of the ligand, summing over all atoms of the ligand. These combined scores were renormalized, for each protein, to allow comparison between the two FABP systems.

	RESULTS

Top Abstract Introduction Methods Results Discussion Conclusion References

Molecular dynamics computer simulations of two fatty acid-binding proteins with the same fatty acid ligand are described. The two proteins, adipocyte fatty acid-binding protein (A-FABP) and human muscle fatty acid-binding protein (M-FABP), differ by 0.7 Å RMS for backbone atoms in the x-ray structures. The saturated C-18 fatty acid, stearic acid, has different conformations for the terminal regions of the alkane chain in the two structures. This is shown with ribbon diagrams in Fig. 1, A and B. Both simulations were performed with the CHARMm molecular dynamics computer program and consist of calculated trajectories totaling 1 ns, after equilibration, for both proteins.

View larger version (59K): [in this window] [in a new window]   FIGURE 1   Ribbon diagram of the FABP structure, with a CPK representation for the fatty acid ligand. Both A-FABP and M-FABP are shown. The fatty acid ligand adopts a different structure in the binding cavity of the two proteins, with the adipocyte ligand structure being more extended than the muscle stearic acid structure.

An initial concern was to assess the appropriate charge state of the ligand. It has been assumed by several groups (e.g., Young et al., 1994; Zanotti et al., 1994; Rich and Evans, 1996) that the carboxylate group is charged in these two proteins. When the same assumption was made for the CHARMm calculations, the trajectory rapidly diverged from the crystal structure. A DELPHI calculation of the pK at the carboxylate group (Gilson et al., 1988) suggested that the headgroup should be neutral. With this choice of charge state, stable molecular dynamics runs with 1 Å C- RMS from the crystal structure were obtained. The ligand formed a stable hydrogen bonding network in the binding pocket, as seen in Fig. 2.

View larger version (29K): [in this window] [in a new window]   FIGURE 2   H-bonding pattern observed in the crystal structure and the placement of hydrogens predicted from DELPHI and CHARMm calculations. In particular, the carboxylate group of the fatty acid is predicted to be neutral rather than charged. Similar calculations on the I-FABP predict a charged carboxylate group.

Similar electrostatic calculations for a third protein, intestinal fatty acid-binding protein (I-FABP), suggest that the ligand should be ionized (Tychko and Woolf, 1997). Multinanosecond dynamics calculations on the I-FABP system have remained stable with this choice of charge state (Tychko and Woolf, 1997). This finding is consistent with the only currently available experimental evidence, which suggests a charged carboxylate group (Cistola et al., 1989).

Several properties of the trajectories were analyzed. The initial focus was on trajectory-averaged structural properties, dynamic features, interaction energies, water behavior, and motional covariance. These results were then used to develop a scale that ranks the relative importance of various residues in binding.

Structural properties

Average structural properties of the trajectory can be used to assess the reliability of the results and reveal the most probable conformations that result from dynamic fluctuations. The trajectories were stable and well defined relative to the x-ray structures. Fig. 3 A shows the C- and heavy atom RMS deviations from the trajectory-averaged structure. Similar RMS deviations from the x-ray crystal structure were observed. These results suggest that the trajectories may be representative of the actual nanosecond time scale motions of FABPs. There were no indications of problems with the selected charge state, potential function, or boundary conditions.

View larger version (39K): [in this window] [in a new window]   FIGURE 3   (A) Root mean squared deviation from the trajectory-averaged structure for C atoms and all heavy atoms. The bars above zero are C, and the negative bars are all heavy atoms. The two simulations remained near the x-ray structures for the 1-ns trajectory production. (B) Root mean squared deviation from the trajectory-averaged structure for the fatty acid ligand. The deviations from the average structure are different between the two fatty acids. In particular, the A-FABP has an increasing set of deviations toward the chain end.

The general nature of the RMS deviations from the average structure was similar in the two simulations. Not surprisingly, the turns and terminal regions were more mobile than the strands and helices. The RMS deviations were generally a bit larger for A-FABP than for M-FABP. In both cases, there was a strong periodicity to the pattern of deviations, with low RMS -strands alternating regularly with more mobile turns.

Considerable differences were apparent between the RMS deviations for stearic acid in the two systems, as seen in Fig. 3 B. In M-FABP, the ligand had the greatest RMS deviations at the headgroup region and a relatively smooth set of lower RMS deviations along the alkane chain. In A-FABP, its RMS deviations increased significantly along the alkane chain, to nearly 3.5 Å at the terminus. These distinctions indicated that a different set of fatty acid motions was present in the two simulations.

Fig. 4 presents the average and RMS deviations for backbone dihedral angles in both proteins. The calculated RMS dihedral deviations were generally less than 10°, consistent with the small -carbon and heavy atom RMS deviations. Both similarities and differences were observed in the two simulations. The greatest flexibility was apparent in the turn regions. A turn with large deviations is seen in A-FABP, D87, and G88, with fluctuations of 56° and 63° in  and 71° and 57° in . These residues are part of the turn connecting the F and G -strands. Similarly, large fluctuations were seen in the M-FABP for the region. Differences were seen in the size of fluctuations in two other turn regions. For example, in A-FABP, larger fluctuations were seen in the turn connecting the G and H -strands, whereas in M-FABP, the H-I turn had larger fluctuations.

View larger version (86K): [in this window] [in a new window]   FIGURE 4   Phi-psi values for all residues in both simulations. The relative deviations from the x-ray values were small. There was no evidence of large conformational changes during the simulations. The deviations from the average values were greater in turn regions than in regions of defined secondary structure.

Dynamic properties

Examining the temporal properties of the trajectories can reveal detailed molecular motions of the system. Fig. 5 shows the time series for the ligand dihedral angles. It is interesting to note that the lipid was not fixed in one conformation, despite strong interactions with the binding pocket. In both simulations dihedral transitions occurred, conserving the overall shape of the ligand. Tables 1 and 2 present the total number of transitions for each dihedral and a breakdown of the i:i + 1, i:i + 2, i:i + 3, and i:i + 4 pairs of transitions observed within 3.0 ps of one another. A meaningful test of concerted transitions (e.g., Brown et al., 1995) requires a larger number of transitions than observed with the current trajectory. In particular, a test for concertedness would require statistics to separate the number of randomly occurring independent pairs from the concerted transitions (Brown et al., 1995). Despite the low numbers of transitions, it is interesting to note that the number of i:i + 2 transitions is often higher than the numbers of other types of observed pairs. This suggests that similar motional behavior conserving the overall shape of the ligand is found for both alkane chains in bilayers and in these holo protein simulations (e.g., Venable et al., 1993).

View larger version (38K): [in this window] [in a new window]   FIGURE 5   Dihedral changes along the fatty acid alkane chains. Transitions conserving the overall shape of the fatty acid were observed in both simulations. The results suggest that the behavior of a fatty acid within the FABP cavity is in some ways similar to the behavior of individual alkane chains within a lipid bilayer.

The pattern of ligand dihedral transitions was different in the two simulations, consistent with the observed differences in ligand RMS deviations. However, the two analyses revealed different trends. In M-FABP, the RMS deviations and the dihedral time series implied similar mobilities along the alkane chain. In A-FABP, the RMS deviations rise smoothly along the alkane chain, whereas the frequency of dihedral transitions was high near the headgroup and the terminus, but not in the middle of the chain.

Interaction energies

Throughout the course of a dynamics trajectory, a given subset of atoms will experience a range of interactions. The instantaneous interaction energies between groups can be calculated and binned over the trajectory to produce a probability histogram. These histograms contain information that is not available through inspection of a single crystal structure. A converged histogram qualitatively describes the interaction enthalpies that are involved with a particular set of conformations. For this reason, they give indications of the average energetic connections between components of the system. Another way to view the histogram information is that it is related to a linear response model for the thermodynamics of transfer (e.g., Aqvist et al., 1994; Aqvist and Hansson, 1996). In this approach, the changes in average interaction energies, separated into vdW and electrostatic terms, are used to computationally estimate the thermodynamics.

The energetic interactions of the ligand with the environment were separated into headgroup and acyl tail components. Fig. 6 A shows the interaction energy distributions for these two groups. The interactions were roughly twice as strong for the acyl tail as compared to the headgroup: 40 kcal/mol versus 20 kcal/mol. Whereas the headgroup interactions were largely electrostatic, the largest contribution to the total interaction energy arises from the many small favorable vdW interactions along the acyl chain.

View larger version (29K): [in this window] [in a new window]   FIGURE 6   (A) Interaction energies between the carboxylate headgroup and the tail with the rest of the system. The energies were calculated for all conformations saved during the trajectory production and then binned, histogrammed, and normalized for this figure. The result shows a strong interaction energy for the tail relative to the head. (B) Interaction energies for the individual methylene groups along the fatty acid for both FABPs. Notice the differences in distributions for the two proteins. In the M-FABP, the C-18 group is the most strongly interacting, whereas in the A-FABP, the C-18 is, on average, much more weakly interacting.

A further analysis of the ligand interaction energies calculated the contribution of the individual methylene groups to the total. Fig. 6 B shows that these interaction energies were not uniform along the length of the fatty acid chain. The strength of interaction varied from 5 to 0 kcal/mol, with an average interaction energy of 2.5 kcal/mol. There was a range of distribution widths from 2 to 5 kcal/mol. This diversity of distribution widths and means is related to the differences in the binding cavity between the two proteins. An average strong interaction energy between a methylene group and its surroundings indicates a favorable set of interactions mediated by vdW effects. For example, the M-FABP has a C-18 group, on average, with 4.5 kcal/mol of interaction energy. In contrast, the A-FABP simulation suggests that this same group has an average of 2.0 kcal/mol with a much broader distribution width of nearly 5 kcal/mol. The difference is related to the relatively extended form of the ligand in the A-FABP structure, compared to the curved form of the ligand in the M-FABP structure.

A consideration of the interaction energies for individual amino acids provides further insight into the possible functional role of particular residues. Fig. 7 illustrates the individual amino acid interactions with the ligand. A strong amino acid interaction with the fatty acid may indicate a role of the side chain in fatty acid recognition and binding. Breaking the interactions into vdW and electrostatic terms further elucidates the type of enthalpic connection between amino acid and ligand. It is not immediately obvious from the crystal structure which residues have strong interactions. Furthermore, the set of conformations from the trajectory can provide insight into the range of energetic fluctuations, whereas a single conformation could incorrectly suggest a much stronger or weaker interaction. For example, the crystal structure would suggest equally strong roles for the triad near the ligand headgroup of R106-R126-Y128. But the amino acid interaction energies suggest that R126 is much stronger in the M-FABP than the A-FABP binding interaction, whereas the R106 is stronger in the A-FABP than in the M-FABP.

View larger version (29K): [in this window] [in a new window]   FIGURE 7   Individual interaction energies for the amino acids with the ligand. The grey bars show the average, with the black tips added to the average to give an idea of the RMS deviations. This assumes a roughly Gaussian distribution for the interaction energies. The interaction energies are further analyzed in terms of the electrostatic and vdW contributions in the two lower panels of the figures. For example, the interaction of R126 in M-FABP is largely electrostatic, as can be seen from the breakdown.

Water behavior

The water found in the FABP-binding cavity is intrinsically involved in the binding process. For example, a hydrogen-bonding network involving water is inferred from the crystal structure (Fig. 2). The motional behavior of this water has not been measured experimentally. However, this sort of information is readily available from a molecular dynamics trajectory.

The full trajectory was examined to identify which waters were most important in ligand binding. These were determined by counting how many times each water passed within 4.0 Å of any heavy atom in the ligand over the course of the simulation. The effective diffusion constants of the 40 waters most frequently in contact with the ligand were calculated from the time derivative of the mean squared displacement function: Interestingly, these effective diffusion constants ranged over two orders of magnitude. Whereas some waters stayed in the binding pocket throughout the simulation, others migrated in and out. The least constrained waters near the ligand were roughly half as mobile as waters in bulk molecular dynamics simulations (Brooks et al., 1988). Other waters were up to 100 times less mobile.

Tables 3 and 4 list the most frequent ligand neighbors for these waters, in A-FABP and M-FABP, respectively. For example, the first water listed in Table 4 bridged R106 and the headgroup. Its most frequent neighbors were the headgroup atoms O1, O₂, C1, and C2. The total number of contacts, rather than a frequency, is presented, because the count was performed on the basis of number of neighbors in each conformation rather than the number of times the water was near the ligand, regardless of neighbor. The most restricted water diffused ~10 times more slowly than the water most frequently in contact with the fatty acid.

Covariance analysis

It has been commented that the average motional properties of an individual atom in a molecular dynamics trajectory are often reasonable, but that the determination of coupling between atomic motions is more difficult to describe accurately and to connect with experiment (e.g., Clarage et al., 1995). This coupling of atomic motions is an important point for the application of molecular dynamics methods to the analysis of the molecular motions. To address the coupling of molecular motions in the current simulations, the zero-time covariance function was calculated as an average over the full trajectory: (where A and B are the fluctuations from the average of the atom position and the brackets indicate an average over the full trajectory).

This expression gives an estimate of the immediacy of coupled motion between sets of atoms. The two sets chosen for analysis consist of the protein or the subset of 40 waters, each coupled with the fatty acid ligand motion. This is shown in Fig. 8. The coupling suggests that quite different motions are involved within the two systems.

View larger version (89K): [in this window] [in a new window]   FIGURE 8   Covariance analysis of the motions of protein and waters with the stearic acid ligand. (A) The protein sequence along the bottom and the ligand as a column. The color scheme runs from correlated motion in red to anticorrelated motion in blue. No correlation in motion is presented in yellow. A schematic view of the protein structure is shown along the bottom of A. The differences in correlated motional behavior of protein and FA are striking. Notice, for example, the sets of red and blue correlations with a subset of the stearic acid for the M-FABP, in contrast to the weaker correlations for the A-FABP that run the full length of the fatty acid. (B) A similar plot of the correlated motions of the stearic acid ligand, with neighboring waters in the binding cavity. The waters are presented along the x axis, with ligand along the y axis. Again, the scale runs from red correlations to blue anticorrelations.

The M-FABP system is seen to have correlated and anticorrelated motions that alternate between the headgroup and along the length of the alkane chain of the fatty acid. The A-FABP system has correlations along the length of the chain. The stearic acid in A-FABP has much stronger interactions with the water than with the protein. This suggests that a different mechanism of selectivity could be involved between the two systems. The water could play a larger role in the binding for the A-FABP system than for the M-FABP.

	DISCUSSION

Top Abstract Introduction Methods Results Discussion Conclusion References

The present molecular dynamics computations provide initial insights into the molecular motions, average structural and dynamic properties, and interaction energies of the same ligand within two different FABP structures. The results are important for several reasons. First, the finding of a neutral headgroup for the fatty acid in these two proteins suggests that analysis of the selectivity based solely on a charged headgroup is probably not correct. Second, the strong hydrophobic interactions of the tail region with the protein suggest that hydrophobic effects are important in the binding selectivity. Third, the motional differences within the two proteins suggest that there are differences in the mechanisms for selectivity between the two proteins. For example, the pattern of correlated and anticorrelated motions along the length of the fatty acid ligand is strikingly different between M-FABP and A-FABP.

The prediction of the protonated carboxylate group of the fatty acid is testable either by NMR methods similar to that used for I-FABP (Cistola et al., 1989) or by Fourier transform infrared methods (Gericke et al., 1997). Similar molecular dynamics calculations applied to the I-FABP system suggest a charged carboxylate group, consistent with the experimental information (Cistola et al., 1989). It is intriguing to note that the extra proton suggested by the present results could, in principle, be distributed throughout the headgroup region. The extra proton could be found, at times, on the bridging water (H₃0⁺) or return to the R106 side chain. An additional test of this concept would be to calculate the predicted proton locations for the hexadecanesulfonic acid ligand in the A-FABP structure (LaLonde et al., 1994b). The effective free energy barriers for proton transfers throughout the headgroup region may be calculated with an extension of the current molecular dynamics calculations (e.g., Pomes and Roux, 1996).

The differences in zero-time covariance functions between the two systems suggest that water plays a much greater functional role in the A-FABP system than in the M-FABP. This is intriguing in that changes in the protein residues lining the interior cavity could then have both a direct effect on the ligand and an indirect effect mediated through the internal water.

Suggestions for mutagenesis

It is plausible to assume that residues with strong energetic interactions and strong positional covariances with the ligand are important in binding and that such residues are good targets for mutagenesis. The residues selected in this way are not immediately obvious from the x-ray crystal structures. A method combining interaction energy and covariance analysis could thus be of general utility in selecting functionally important residues.

A general scoring procedure for this approach was developed that combined the interaction energy and zero-time covariance on the same scale. The method normalized the interactions and covariances, and then summed over the full set of ligand heavy atoms for a normalized value to describe the relative importance of the individual amino acid to the system. This second normalization was necessary to allow comparison between the two protein systems. The results of this analysis are presented in Fig. 9. The upper panel of the figure shows the individual scores, whereas the middle and lower panels show the sum and differences of the individual scores. Another way of thinking about the approach is that it attempts to determine the largest interaction free energy between amino acid residues and ligand. In this regard, the interaction energies are an interaction enthalpy, and the covariance estimate is an interaction entropy. The importance scale then attempts to weigh the terms equally, but a more elaborate formalism could be imagined that describes a scale based on thermodynamic arguments.

View larger version (39K): [in this window] [in a new window]   FIGURE 9   The relative importance of each residue to the binding function of the protein was qualitativly predicted by using a combination of the interaction energy and the zero-time covariance analysis. The resulting normalized scale is presented with both individual scores and the sum and difference of the scores. The normalization was necessary to compare the two proteins, because the covariance and interactions differed in overall range and magnitude. The sum may provide an indication of residues that are important for both systems, whereas the difference may show residues that are more important for specificity in one protein than another.

Several regions of the two proteins show intriguing differences in their importance scores. The triad of residues near the headgroup region had very different scores in the two systems. A-FABP had very high importance scores at R106 and Y128, and a fairly high score for R126. M-FABP had a very high score at R126, and insignificant scores for R106 and Y128. This may be related to the stronger motional coupling with water observed in the A-FABP relative to the M-FABP and to the shifts in average interaction energies seen between the two FABPs.

A second intriguing region was the first -helix. The combined scale indicated that F16 and M20 (conserved in both structures) were important. These residues are near the portal region and may be important for maintaining the integrity of the internal binding site. Their main contribution was through vdW interactions in both simulations. It is interesting to note that Y19 is the phosphorylation site that is conserved in both proteins (Buelt et al., 1992). It may be that a change in phosphorylation state of this residue could alter the dynamics of F16 and M20, thereby influencing the kinetics or stabilization of binding.

Another interesting region that was highlighted as being important for binding in both systems is the turn between -strands E and F. This turn region has four conserved charged residues in a row in both structures (D76, D77, R78, and K79). This set of charges is preceded by a conserved A75. The A75 and D76 residues show especially high scores on the combined scale and suggest that the properties of this turn region may be important functionally. This region is near the bend in the alkane chain seen in the M-FABP structure. The D76 is the residue nearest to direct contact with the ligand. The D77 and R78 are much stronger on the combined scale for M-FABP than for A-FABP. This suggests that changes in other regions of the protein structure could be altering the affinity via changes in the E-F turn.

The second -helix also interacts with the alkane chain of the ligand. The interaction was much stronger for M-FABP than for A-FABP. In particular, P38 (conserved) and T36 (Met in A-FABP) showed high combined scores. The P38 interaction was largely a vdW interaction with the alkane chain. A correlated motion with the headgroup was seen, along with an anticorrelated motion with the middle and ends of the alkane chain. This implies that the two -helices may play a significant role in helping to create the binding site for the hydrophobic ligand and that their roles could differ between the two proteins.

The L117 in M-FABP played a much bigger role in energetics and combined motion than did the C117 in A-FABP. This site, in the Ith -strand, may suggest an individual amino acid that is interacting differently between the proteins. The main sites of vdW interaction for L117 are the C3 and C4 atoms of the ligand. The C117 is less strongly interacting and significantly less strongly correlated with the fatty acid ligand.

An additional point of comparison between the two proteins was a consideration of the internal cavity side chains that could provide sites for water interaction. This leads to four sites of difference between the two proteins. The sites are K44, N100, L108, and G111 in M-FABP, which change to V44, K100, R108, and E111 in A-FABP. The sites are relatively removed from the ligand itself, but are near presumed entrance/exit sites for water and may thus have an effect on the relative importance of water:protein motions and interactions within the system. Thus the intriguing result that the water motions are more tightly coupled in A-FABP could be related to the changes at these four sites. Mutations at these locations could thus have an indirect effect on the binding affinity by changing the relative behavior of water within the cavity.

It needs to be emphasized that the proposed importance scale is not expected to be accurate relative to detailed lamba-coupled relative free energy calculations (e.g., Kollman, 1993). The scale could give some qualitative insights into candidate sites for mutagenesis and suggest ways that the site might be involved with function. A more detailed analysis and comparison of free energy differences produced by mutagenesis will require quantitative calculation using free energy perturbation methods (e.g., Kollman, 1993). Alternatively, qualitative free energy methods, which may produce estimates for rank-ordering of mutations to confidence greater than the importance scale, can be used (e.g., Aqvist et al., 1994).

Relation to previous simulations

Zanotti et al. attempted to construct a path for ligand entry/exit in I-FABP by a set of short high-temperature molecular dynamics runs (Zanotti et al., 1994). The current simulations used longer simulations for M-FABP and A-FABP with larger solvent surroundings, but were not directed at exploring the mode of ligand entry/exit. Attempts to determine reaction coordinates in fully detailed atomic structures are difficult (e.g., Brooks et al., 1988). An additional difficulty is that the time scales for ligand entry/exit are long compared to molecular dynamics trajectories. For example, recent experimental measurements were made of k_off (for leaving the binding site) of 1.8/s for oleate and palmitate in I-FABP (Richieri et al., 1996). The types of motion observed in the current trajectories do provide some suggestions for sites of the protein with greater motional flexibility. Thus the trajectories show larger motions of the D-E gap and the -I and -II helices relative to the rest of the structure. This is consistent with the results of Zanotti et al. and supports the idea that these regions of the protein could act as portal sites for the entry/exit of ligands (e.g., Hodsdon and Cistola, 1997).

Rich and Evans (1996) examined electrostatic effects in the A-FABP system. An assumption was made that the ligand was normally charged and changes were made in electrostatics of the binding cavity to account for the importance of charge on the binding process. The present calculations examined two possible states for the ligand, but did not vary the charges throughout the ligand-binding site region. The Rich and Evans approach was to assume that either all charged forms should be used for amino acids near the binding site and the ligand, or that all neutral forms should be used. This is equivalent to using the DISCOVER simulations to address two possible states of a much larger grid of possible protonation states. Their simulations used very short equilibration times (10 ps), short total dynamics (100 ps), a small nonbonded cutoff (10 Å), a distance-dependent dielectric, and a small solvation model. The current calculations suggest that the headgroup regions in the M-FABP and A-FABP systems are different from that in the I-FABP system. Thus direct electrostatic effects may be less important in the M-FABP and A-FABP systems than in the I-FABP system.

Young et al. attempted to rationalize the binding affinity of three ligands in M-FABP by calculation of the change in conformational energy of the ligand bound in the protein relative to solution (Young et al., 1994). Their approach was secondary to the x-ray structural results that they report. The calculation used a single conformation for the solution state and the holo protein x-ray structure conformation for the bound state. But the free energy of transfer is related to changes in entropy and enthalpy of the ligand. A single conformation does not sample the entropy. Furthermore, the protein can provide interaction enthalpy that dramatically changes the free energy surface regarding the conformations of the ligand. Hence a Boltzmann-weighted ensemble of structures for the ligand in the two states would be expected to differ.

	CONCLUSION

Top Abstract Introduction Methods Results Discussion Conclusion References

Molecular dynamics calculations are presented on M-FABP and A-FABP complexed with stearic acid. The results provide initial insights into the molecular motions and ligand interactions characteristic of this family of proteins. In particular, the results suggest that vdW interactions along the hydrophobic cavity may be quite important for specificity. This finding could be related to the ability of FABPs to discriminate based on fatty acid chain length and saturation state. In the long term, this work hopes to help determine the relationship between tertiary structure and binding affinity, leading to the eventual rational design of new proteins capable of binding specific fatty acids.