Conformational changes of globular proteins upon adsorption on a hydrophobic surface

This paper presents a study of protein adsorption and denaturation using coarse-grained Monte Carlo simulations with simulated annealing. Intermolecular interactions are modeled using the Miyazawa– Jernigan (MJ) knowledge-based potential for an implicit solvent. Three different hydrophobicity scales are tested for adsorption of fibronectin on a hydrophobic surface. The hydrophobic scale BULDG was chosen for further analysis due to its greater stability during heating and its partial regenerative ability upon slow cooling. Differences between helical and sheet structures are observed upon denaturation – a -helices undergo spreading of their native helical order to an elliptical perturbed shape, while b -sheets transform into random coils and other more structured conformations. Electronic calculations carried out on rebuilt all-atom coordinates of adsorbed lysozymes revealed consistent destabilization of helices, while beta sheets show a greater variety of trends.


Introduction
When a protein is placed outside its natural medium, its structure and function are usually affected. 1 Near a solid surface, proteins tend to adsorb, and depending on the protein and surface properties, they may undergo partial or significant conformational changes. Many modern medical devices rely on the direct contact between bodily fluids and artificial surfaces, and their use is often complicated by immunological reactions induced by adsorption of proteins on the artificial interface. 2 In addition, many biomimetic devices 3 are based on catalytic and reception functions of immobilized proteins, such as biosensors, bioreactors, and artificial organs, which also suffer from the problem of adsorption and consequent denaturation of bioactive materials. Proteins tend to adhere more strongly to nonpolar than to polar substrates since the nonpolar environment destabilizes the protein and thereby facilitates a conformational rearrangement that leads to strong protein-surface hydrophobic interactions. 4 Protein adsorption and denaturation is a common but very complicated phenomenon, and to date there has been inadequate understanding of the adsorption process. 5 It is driven largely by an entropy gain arising from the release of surface adsorbed water molecules in addition to the structural rearrangement of the protein. It is possible to define three basic levels of events: (1) microscale that analyzes events at the atomic level, (2) motions of macromolecular clusters like alpha helices and beta sheets which can be of significance in phase transition processes, and (3) macroscale events that take into account forces acting on the whole protein, for example, hydrodynamics of a viscous fluid. 1,6 Some structural information can be obtained from experimental techniques, 1,3,7,8 such as circular dichroism (CD), Fourier-transform infrared (FTIR) spectroscopy, scanning angle reflectometry (SAR), and atomic force microscopy (AFM) imaging. These techniques provide evidence for structural transitions of proteins in adsorbed states 5 such as b-sheets to a-helices transitions. The dynamics of singlemolecule adsorption and denaturation can now be measured using Förster resonance energy transfer. 9 Adsorption of globular proteins follows 2-step kinetics that lead to partial denaturation of the protein, exposing the hydrophobic residues to the surface. 10 For ambient conditions relevant to biological solutions, proteins most generally approach the surface in their native state, and the initial adsorption is followed by a slow structural rearrangement, corresponding to a relaxation process with macroscopically observable effects including an altered secondary structure. 6 AFM experiments point to conformational changes of the adsorbed proteins that occur over periods extending to minutes or hours. 3,11 This kind of 2-step adsorption kinetics is a manifestation of frustrated macromolecules whose conformational space is restricted to a small number of most stable conformations. 12,13 Notably, cooperative adsorption and protein aggregation in solution and on the surface may influence the adsorption kinetics as well as the resulting protein layer structure, 4,14,15 though, it has recently been shown that proteinprotein interactions have little influence on the structure of the lyz on highly hydrophobic surfaces. 16 During the last decade, the role of computational modeling in structural biology has increased substantially. 17 Since protein adsorption is associated with complex molecular events beyond the nanosecond timescale, fully atomistic models require biasing techniques to overcome free energy minima to allow for denaturation. 18 One of the most promising tools is the coarse-grained (CG) approaches, as they demand much less computational effort than all-atom models, and therefore allow for longer time scales of molecular events. A common CG model of proteins represents each amino acid as having two interacting sites (a-carbon and side chains). CG models of proteins provide good representation of the conformational space of proteins near the folded state region, 19,20 as well as those adsorbed onto a hydrophobic surface. 21 A recent model developed by Wei et al. 22 revealed the power of such parametrized coarse-grained models in predicting protein adsorption behavior.
The problem of protein adsorption deals with a highly complex energetic landscape that is characterized by a large number of local minima, for which standard simulation algorithms result in trapping at local minima for long periods. Following the work in ref. 23, we use simulated annealing to predict the structure of the denatured protein upon adsorption. In this method, we search for equilibrium states of the system through heating followed by slow cooling, in order to surmount local minima that kinetically prevent the system from reaching energetically favorable states at short simulation times. We have previously determined the thermal limits of stability of proteins in the bulk and on the surface, 24 which are used in the design of annealing schedules used in this work. The MC algorithm is based on the coarse-grained model developed by Haliloglu et al. 25 using the Miyazawa-Jernigan knowledgebased potential for residue-residue interactions. 26 The specific surface-amino acid interactions are modeled using different hydrophobicity (hb) scales for fibronectin (fn), lysozyme (lyz), and a small helical peptide. Lysozyme CG coordinates are converted to all-atom presentation and evaluation of intramolecular interactions at secondary structure elements at the B3LYP level of theory is performed.

Methods
We base our model on the work of Jernigan and coworkers 19,20 developed for the investigation of conformational dynamics of several globular proteins. In this model, each amino acid residue i is represented by two principal interacting sites, an a-carbon (C a ) atom and a sidechain center S i , as is depicted in Fig. 1.
In this model, the energy of a given conformation F is a sum over all non-bonded interactions (E L ) between non-neighboring residue pairs and bonded interactions (E S ), These interactions are given by a potential of mean force discretized by intervals of 0.4 Å in the range 2.0 r r ij r 12.4 Å, where r ij is the distance between interacting sites, and 101 for the angle-dependent bonded potentials. The simulation Fig. 1 Schematic representation of the model. The conformation of the backbone for a protein of N residues is defined by 3N À 6 variables: N À 1 bonds of length l i connecting C a iÀ1 and C a i , N À 2 virtual bond angles y i of the ith a-carbon, and N À 3 dihedral angles f i due to torsion about bond l i . The conformation of the sidechains with respect to the backbone is defined in terms of analogous variables l s i , y s i and f s i , which refer to the respective virtual bond length, bond angle, and dihedral angle of the sidechain. The simulation initiates with the protein at a random orientation with respect to the surface vector n 0 . Simulation proceeds with conformational moves as well as translation and rotation of the protein relative to the surface.
proceeds through perturbations of all dihedral angles and spatial perturbations of backbone monomers.
The surface potential is taken to be an integrated form of the Lennard-Jones 6-12 potential in cylindrical coordinates, where e w is a hydrophobicity-dependent constant that is specific for each residue and is taken to be equal to the potential interaction of the most hydrophobic residue and the interacting residue for the given hydrophobicity scale. Integration yields a one-dimensional potential that depends on the distance d of the particle from the surface, The effective range of E surf is set equal to the range of longrange interactions, and its value is normalized according to the most hydrophobic residue. Hydrophobic surfaces possess a significant water depletion zone and hence the hydrophilic reference near the surface is not optimal. Future work will aim at improving the chosen force field by introducing a distance-dependent MJ potential with water as a reference. Fig. 2 depicts the proteins that were chosen for our analysis: lysozyme (lyz) for its prevalence of a-helices, fibronectin (fn) for its prevalence of b-sheets, and a short single helix-containing peptide. Since realistic timescales of protein adsorption 27 are beyond computational capabilities of molecular models, we use simulated annealing to overcome local energy minima. We carried out simulations in the reduced temperature (T* = kT/e w ) range of 1.0-5.0, which corresponds to 300-400 K, and covers the most characteristic thermal transitions that result in total denaturation of tertiary and secondary structural elements of these proteins. 28 T* = 1 is assigned to the crystallographic temperature of the protein and T* = 0 corresponds to 273 K, with the following conversion for other temperatures: T (Kelvin) = T (reduced) Â [T (crystallographic) À 273)] + 273. Each simulation begins with equilibration of the protein in the absence of the surface under standard conditions for 10 5 MC steps. Isothermal adsorption and equilibration on the surface follows for additional 10 5 MC steps. Non-periodic boundaries with unlimited box dimensions are used to model infinitely dilute solution conditions. We have previously shown 24 that the force field and MC method cannot predict denaturation under isothermal conditions, and have analyzed the thermal limits of denaturation of three proteins under annealing. Accordingly, simulated annealing simulations with maximal reduced temperatures of 2.5, 3.5 and 5.05 were carried out in this work for lyz, fn, and the peptide, respectively. A combination of rotation, translation and angle perturbations was used to achieve equilibrium conformations. Several independent runs were carried out for each of the annealing simulations. For each run, 105 MC steps were performed at each temperature, of which 7 Â 104 MC equilibration steps were followed by 3 Â 104 MC steps of sampling every 102 steps.

Results and discussion
Hydrophobicity scales define the relative hydrophobicity of amino acid residues, where the more hydrophobic residues have higher positive values. Out of 37 hb scales examined by Cornette et al., 29 three (BULDG, 19,20,30 ABODR 21,31 and MIJER 32 ) were deemed to be most relevant for residue surface interactions, and were tested in adsorption simulations of fibronectin. The BULDG scale is based on the effect of varying concentrations of each amino acid on the surface tension of water. The ABODR scale measures the mobility of the amino acids on chromatography paper and might therefore be appropriate for adsorption processes. The MIJER scale is based on contact energies between amino acids within proteins taken from crystallographic data, and is the basis of the knowledge-based potential set used in this simulation. Most of the remaining hb scales are based on  Leu) appear to be in much better agreement among the three scales than the hydrophilic residues. But generally the statistical MIJER hb scale assigns more positive indices than the experimental hb scales.
The impact of the hb scale on the average fraction of native contacts, f nc , and the radius of gyration (R g ) is shown in Fig. 4 for fibronectin. Both f nc and R g reveal significant differences among the hb scales during the course of the simulation. For each scale, fluctuations in these measures seem to stabilize for T* o 3.0 during the cooling cycle. This temperature corresponds to the sharp decrease in the perpendicular component of R g (Fig. 5B), beyond which the local structural rearrangement occurs. R g under the ABDOR hp scale seems to be least affected by the variation in temperature and MIJER is the most affected, especially during heating. The MIJER scale also reveals the least stable native contacts. The BULDG scale leads to the greatest reduction in R g on the one hand, while retaining the most number of native contacts upon heating. Moreover, the BULDG scale shows the greatest ability for regeneration of native contacts for the adsorbed protein upon cooling. The BULDG scale has a high amphipathic index, which encompasses both hydrophobic and hydrophilic elements of native amphipathic a-helices. 29 It is likely that there is a close similarity between intermediate folding stages when secondary structural elements are formed in solution and advanced stages of protein adsorption when secondary structures begin to denature. Hence, we concentrate on the BULDG hb scale below.
The fraction of native contacts and the components of the radius of gyration as a function of annealing temperature for BULDG scale are plotted in Fig. 5 for the peptide, lyz, and fn, and the changes in R g are summarized in Table 1. Both the fraction of native contacts and the radius of gyration change significantly during the course of annealing. The adsorbed proteins settle at a conformational landscape characterized by 30-50% reduction in native contacts around T* = 3.0, following the collapse of the perpendicular component of R g and expansion in the parallel direction. These values are retained in the cooling cycle, suggesting that we are sampling a characteristic conformational landscape of the adsorbed protein. AFM data reveal high statistical variability in the dimensions of the adsorbed proteins, 33 but the detected shape that has been most commonly reported is ellipsoidal with an aspect ratio a r (the ratio of lateral to transverse radii components with respect to    View Article Online the surface) ranging from 0.6-0.8, which is somewhat lower that what we observe in our simulation. Contact maps showing the average contact probability between different residues before and after annealing, are shown in Fig. 6. For the peptide, the helix persists throughout adsorption, but the internal helical morphology is altered at residues 12-15, where most of the hydrophobic elements are located ( Fig. 6A and B). Alteration of the native helical structure causes the appearance of some non-native tertiary contacts as well (between monomers 7-9 and [17][18][19]. For fn, we observe significant denaturation and appearance of several secondary structure elements (Fig. 6C and D). Beta structures 1, 5 and 6 denature completely while 3 new secondary order clusters appear from monomers 65 to 90, one with a characteristic signature of the b sheet. Beta structure 4 undergoes significant perturbations as well, but retains its initial position. In contrast, the contact maps of lyz reveal a more robust protein where most secondary alpha structures and the beta 2 sheet sustain the annealing process, whereas beta 1 completely denatures ( Fig. 6E and F). Still, the shape and the area of tertiary contacts on the final map are altered and a new tertiary contact a 6-7,8 appears.
Ramachandran plots allow further examination of the morphological alterations of the internal secondary structures, which are not fully traceable by contact maps. In these plots, the distributions of torsional angles for secondary structure elements are mapped. 34 In the coarse-grained model, we cannot use such plots since each residue is represented by two sites, alpha carbon and side chain. Alternatively, we consider the bending-torsion distribution upon annealing, shown in Fig. 7-9. The spread of most secondary structures shifts and widens either in the bending or torsion dimensions, or both. b2 is the most stable element, maintaining nearly the same contour before and after annealing, as was also apparent from the contact maps (Fig. 7). For lysozymes, on the other hand, we observe that the helical elements are still confined to their characteristic regions on the torsion-bending plot, but undergo a significant stretching of the bending angles. b1 completely shifts while b2 remains relatively unchanged. The conformational transformations during annealing of all a-helices in lysozymes as well as the single long peptide a-helix demonstrate consistent deformation upon reduction of the average torsion angle and an increase in the average bending angle. The same pattern of helical deformations that all a-helices share suggests opening and radial stretching. As the hydrophobic residues approach the surface, the helices tend to be more ellipsoidal with different curvatures. The b-sheets of fibronectin demonstrate higher structural variation during annealing than its helices. However, the degree of deformation does not depend on the average hydrophobicity of the lyz helix or change in the surface potential during annealing. In addition to its characteristic length, the decisive factor in defining helical spreading is its stabilization within the protein structure by long-range interactions with nearest neighbor elements. Similar to the helical elements of the lysozyme, the peptide helix undergoes partial denaturation upon adsorption, characterized by significant torsional untwisting.  The changes in the average distance from the surface hDi of the residues for the proteins and peptide are shown in Fig. 10-12. After adsorption, lyz monomers are on average 43% closer to the surface relative to the initial distances, while fn monomers are on average 57% closer. However, the final average distance of the monomers is similar for the two proteins (hDi lyz = 17.5 Å and hDi fn = 16.7 Å), indicating greater structural rearrangement of fn. But a correlation between the hb index of the secondary structure elements and the decrease in the average distance during annealing could not be determined. Fig. 11 reveals that lyz residues in the approximate range of 50-100 are closer to the surface than others. This range corresponds to the coiled region and the nearby helices and sheets. There is theoretical evidence that this portion of the protein corresponds to the weaker portion of the protein 35 and denatures first upon adsorption. 22 There appears to be a relationship between the total decrease in the distance (%) and the hb indices of the monomers for the peptide (Fig. 12). It is clearly seen that the hydrophobic amino acids 17-19 (VAL and VAL respectively) approach the surface. Since no such relationship is noticeable for the proteins, this strengthens the notion that the conformational complexity of the protein, which is absent in peptides, determines adsorption and denaturation. For proteins, the long-range monomer-monomer interactions compensate for the surface hydrophobic potential.
Structural changes in the protein upon absorption may crucially affect the activity of the protein. Enzymatic activity is highly affected by the reorganization of the electronic structure of those atoms involved in the reaction. We obtained a rough indication of the potential effects on enzymatic activity due to structural changes upon adsorption by reconstructing the atomistic coordinates of the protein followed by analysis of the electronic properties of the helices and sheets. We used the program PULCHRA 36 for rebuilding the atomic coordinates of the adsorbed lysozyme, averaged over equilibrated annealed conformations, as described by Rotkiewicz and Skolnick. 37 Lysozyme coordinates 1JSF listed in the RCSB Protein Data Bank PDB 38,39 were used as reference. The a-carbons of the Fig. 9 Torsion-bending plots of the peptide before (darker oval) before and after annealing (lighter oval). Fig. 10 Histogram of the average distance (Å) of fibronectin residues from the surface for initial adsorption prior to annealing (light gray) and after annealing (dark gray). The average HB indices of the beta sheets are assigned according to BULDG hydrophobicity scale. Fig. 11 Histogram of the average distance (Å) of lysozyme monomers from the surface is shown using initial adsorption data prior to annealing (light gray) and after annealing (dark gray). The average HB indices of the alpha helices are assigned according to BULDG hydrophobicity scale. final annealing step were submitted as PULCHRA input for generating atomic coordinates of the adsorbed protein, while the a-carbons of the original PDB conformation were used for rebuilding the protein in solution. The 3-dimensional presentations of optimized adsorbed and native conformation, shown in Fig. 13 are plotted with the Structural Identification algorithm (STRIDE 40 ) incorporated in the VMD package. 41,42 Ramachandran plots of the absorbed protein and the rebuilt protein in solution versus original PDB conformation are shown in Fig. 14. The plot of the protein obtained from PDB coordinates show a smaller angle distribution than the rebuilt conformations. Rebuilding combined with MM optimization leads to a more relaxed structure. This is also due in part to the reduced water content in the protein structures obtained from crystallographic data. 37 The observed distortion of helical secondary elements and dissolution of two beta-sheets during the adsorption process are consistent with torsion-bending angle distributions (Fig. 8).
We further analyzed the conformational changes upon adsorption using first principle studies. Single Point Energy (SPE) calculations have been performed using the Jaguar 7.9 tool of the Schrödinger package 43,44 in the gas phase at the B3LYP/6-31G* hybrid DFT level of theory using 6D functions at fully analytic accuracy for the regenerated atomistic coordinates. The dimensionless energy data listed in Table 2 have been calculated according to e* = f * (E i , surf À E i , bulk )/|E pdb À E bulk | where f * is   the normalized average fraction of the energy term, and the quantity |E pdb À E bulk | represents the noise inherent to the rebuilding procedure. The surface potential was neglected during these calculations. From Table 2, we observe that most of the helical structures show a clear decrease in repulsive contributions to the Hamiltonian upon adsorption and an overall increase in the internal energy for these secondary elements, suggesting destabilization of the adsorbed structures. Beta sheet 1, in contrast, shows the opposite trend. Indeed, adsorption of lysozyme on hydrophobic surfaces reveals the appearance of beta sheets as well as a decrease in helical content. 7 The number of hydrogen bonds, determined using the Maestro tool of the Schrödinger package, is listed in Table 3 for the secondary structures and for the entire protein.
Hydrogen bonding within the secondary structures follows a similar trend to that observed in the angle distribution of the adsorbed protein (Fig. 7-9). A significant decrease in the total number of hydrogen bonds in the adsorbed protein corresponds to the observed distortion of the secondary structures of the proteins upon adsorption. The difference in amount of hydrogen bonds is 6% between native PDB conformation and the rebuilt one, giving an indication of the rebuilding efficiency.

Conclusions
In our previous work we analyzed adsorption of globular proteins at constant temperature using a coarse-grained model. Our analysis revealed minor conformational changes for single proteins adsorbed on a hydrophobic surface under the MJ force field. In this work, we performed simulated annealing to obtain the equilibrium adsorbed conformation that would otherwise require unmanageable computational timescales under realistic conditions. We carried out a detailed analysis of protein adsorption on a hydrophobic surface using three different hb scales, and chose to concentrate on the BULDG due to the greater stability of native contacts of this hb scale. A correlation between the hb indices and the distance of the monomers from the surface upon adsorption could not be found. The observed decrease in energy following adsorption was in the 35-75% range and depends on conformational stiffness of the native protein. Differences in the behavior of secondary structures were apparent, in particular, a-helices undergo spreading as the hydrophobic residues approach the hydrophobic surface while the b-sheets lose their structure due to their higher hydrophobicity compared with a helices.
Contributions to the intramolecular interactions were analyzed using all-atom rebuilding of the adsorbed CG model. Dissimilarities in the behavior of helices and sheets were characterized at the DFT level of theory pointing to the appearance of more dense and stabilized beta sheets with a higher hb index. These calculations confirm that secondary structures are destabilized by adsorption. Such a multiscale approach applied to the active site of, e.g., immobilized proteins may potentially be used to predict their functional properties as well as be applied for modeling, e.g., biosensors and bioreactors.