Application of a Spherical Harmonics Expansion Approach for

Article pubs.acs.org/JPCB

Application of a Spherical Harmonics Expansion Approach for Calculating Ligand Density Distributions around Proteins Siddharth Parimal, Steven M. Cramer, and Shekhar Garde* Howard P. Isermann Department of Chemical and Biological Engineering and Center for Biotechnology and Interdisciplinary Studies, Rensselaer Polytechnic Institute, 110 Eighth Street, Troy, New York 12180, United States S Supporting Information *

ABSTRACT: Protein−ligand interactions are central to many biological applications, including molecular recognition, protein formulations, and bioseparations. Complex, multisite ligands can have affinities for different locations on a protein’s surface, depending on the chemical and topographical complementarity. We employ an approach based on the spherical harmonic expansion to calculate spatially resolved three-dimensional atomic density profiles of water and ligands in the vicinity of macromolecules. To illustrate the approach, we first study the hydration of model C180 buckyball solutes, with nonspherical patterns of hydrophobicity/-philicity on their surface. We extend the approach to calculate density profiles of increasingly complex ligands and their constituent groups around a protein (ubiquitin) in aqueous solution. Analysis of density profiles provides information about the binding face of the protein and the preferred orientations of ligands on the binding surface. Our results highlight that the spherical harmonic expansion based approach is easy to implement and efficient for calculation and visualization of threedimensional density profiles around spherically nonsymmetric and topographically and chemically complex solutes.

1. INTRODUCTION Quantifying the interactions of proteins with small molecule ligands is important to many applications, including drug design and discovery,1−6 chromatography-based bioseparations,7−9 and design of biosensors.10,11 These interactions are inherently molecular in nature, with at least three components playing a roleprotein, ligands, and water. Although state-of-the-art biophysical techniques (e.g., X-ray crystallography,12 nuclear magnetic resonance,13−17 or atomic force microscopy18) can provide molecular level information, application of these techniques is time-consuming, and interpretations are sometimes difficult. Further, fundamental physicochemical understanding of these interactions requires the understanding of the role of water structure and dynamics in mediating protein− ligand interactions. To this end, molecular simulations provide a useful complement to experiments.19 Further, the exponential growth in computing power and concomitant advances in experimental methodology have reduced the gap between the time and length-scale that can be simulated on computers and probed by state-of-the-art experiments. Binding of proteins to small molecules is a complex process, and its characteristics depend on various direct and indirect interactions between a protein, ligands, and water molecules. The nature of these interactions can vary from being strong and highly specific to weak and nonspecific. The entire spectrum of interactions from strong, intermediate, to weak is of interest from different applied and fundamental perspectives. For example, to be effective as a drug molecule, one expects interactions between a drug molecule and its target protein to be strong (much more than thermal energy, kT), resulting from © 2014 American Chemical Society

structural and chemical complementarity of the interacting regions near the binding site.20 In contrast, interactions of proteins with solution additives such as urea, guanidine, trimethylamineoxide solutes, or salt ions are examples of weak, nonspecific interactions, which broadly affect the folding−unfolding equilibrium.21,22 Between the strong and weak interactions is the intermediate range, characterizing the interactions between proteins and complex ligands, such as those in multimodal chromatography, which contain multiple chemical moieties and may bind to several sites/regions on the protein surface with intermediate binding strengths. The multiple interactions that occur in multimodal chromatography have been shown to produce unique and important selectivities in protein separation systems.23 Several specialized molecular simulation methods (e.g., steered molecular dynamics and various umbrella sampling methods) are available to study the strong binding scenario. When the binding is weak and nonspecific, the interactions can be studied efficiently by performing the so-called free solution simulations, in which ligand molecules are present at a sufficiently low concentration in aqueous solution surrounding the protein. If the protein−ligand interaction is not many times stronger than thermal energy, one expects efficient sampling of the protein surface by free ligands in solution. The local concentration of the ligand in the vicinity of the protein can then be analyzed by performing various statistical mechanical Received: July 9, 2014 Revised: August 29, 2014 Published: September 8, 2014 13066

dx.doi.org/10.1021/jp506849k | J. Phys. Chem. B 2014, 118, 13066−13076

The Journal of Physical Chemistry B

Article

analyses (e.g., by calculating radial distribution functions,24 preferential interaction coefficients,24−26 and so on). For complex multisite ligands, which have multiple binding sites around the protein, similar free solution simulations can be carried out and calculations of spatial distribution of ligands can be performed based on the probability of finding a ligand atom (or center of mass) within a certain distance from the atoms/ residues on the protein.27,28 These calculations are, however, limited to a local region around the protein, and the directional dependence of ligand distribution in the three-dimensional (3D) space around the proteins becomes difficult to quantify. Grid-based calculations can also be performed by keeping track of the location of ligand atoms during the simulation. However, keeping counts in a spatially resolved manner, especially given that protein itself is a nonspherical soft object, and visualizing them become a complex task. Spherical harmonics, which are sometimes called the Swiss Army knife of mathematical physics, provide a convenient formalism to model three-dimensional distributions of data and have been employed in fields such as spectroscopy,29−32 geophysics,33,34 computer graphics,35,36 and medical imaging.37 For example, the heat flow field of the earth has been represented in terms of spherical harmonics to analyze its global distribution.33 Here we use an approach based on spherical harmonic expansion of atom densities, where the coefficients of the expansion are calculated using simulations of protein in an aqueous solution of free ligands. The main idea is similar to that used by Garde et al.38 to describe water density profiles around complex solutes. We demonstrate applications of such an approach to describe water density around roughly spherical buckyball solutes with nonspherical patterns of hydrophobicity/-philicity on their surface. We also extend the method to calculate density profiles of increasingly complex ligands around the protein ubiquitin in aqueous solutions. Our results show that a spherical harmonic expansion based description of density is efficient and easy to implement. It can be extended to different spatial resolutions and provides information that may be used to gain meaningful physical insights into protein−ligand interactions.

respectively, as shown in Figure 1b. In B1/2, half of the solute is hydrophobic and the other half hydrophilic, the B1/4 pattern

Figure 1. (a) Snapshot of a C180 solute (gray) solvated in water (blue and red wireframe). (b) Hydrophilicity patterns for different C180 solutes (cyan = hydrophilic; brown = hydrophobic). (c) Radial distribution function g(r) of water around different versions of the C180 solute. For comparison, g(r) for completely hydrophobic (dotted brown line) and a completely hydrophilic (dotted cyan line) C180 solute are shown.

has a quarter of the atoms in a contiguous patch hydrophilic and the rest hydrophobic, whereas the Bpatches pattern contains 67 (i.e., 37%) hydrophilic atoms, distributed in five separate patches (one of which is a quadrant in B1/4), distributed randomly on the surface. Intramolecular bonds, bond angles, and dihedral interactions for the C180 solutes were taken from AMBER parm-94 force field for carbon type CA.39 As stated earlier, our model C180 solute contains two types of atoms, hydrophobic and hydrophilic. The hydrophobic atoms are purely repulsive WCA (Weeks−Chandler−Andersen) solutes, whereas the hydrophilic ones have twice the attractive interactions relative to the CA atom type in AMBER. The WCA40,41 scheme separates Lennard-Jones (LJ) potential into purely repulsive and purely attractive contributions as follows:

2. METHODS We report analysis from two sets of simulationsone for spherically symmetric buckyballs in water and another for the protein ubiquitin in aqueous solutions of different ligands. The first set validates the spherical harmonics algorithm to calculate water density around buckyballs with increasingly complex distributions of hydrophobicity/-philicity on their surface, whereas the second set applies the algorithm to calculate ligand distributions around ubiquitin. 2.1. MD Simulations of a C180 Solute in Water. We chose a spherically symmetric buckyball, C180, containing 180 carbon atoms as our reference solute. By manipulating the strength of interaction between water and a specific subset of carbon atoms, we create hydrophobic and hydrophilic patterns on the buckyball surface. These patterns are expected to affect the local density of water. We capture the three-dimensional response reflected in the vicinal water density by calculating the coefficients of spherical harmonics expansion of water density distribution and relating the features of water density to the underlying pattern of hydrophobicity/-philicity on the buckyball surface. We constructed three different hydrophobicity patterns on the C180 solute, referred to as B1/2, B1/4, and Bpatches,

wca att usw (r , λ) = usw (r ) + usw (r )

(1a)

where ⎧ u (r ) + ϵ , r ≤ 21/6σ ⎪ lj sw sw wca usw (r ) = ⎨ ⎪ r > 21/6σsw ⎩ 0,

(1b)

⎧−λ ϵ , r ≤ 21/6σ ⎪ sw sw =⎨ 1/6 ⎪ λulj(r ), r > 2 σsw ⎩

(1c)

att usw (r )

⎡⎛ σ ⎞12 ⎛ σ ⎞6 ⎤ ulj = 4ϵsw ⎢⎜ sw ⎟ − ⎜ sw ⎟ ⎥ ⎝ r ⎠⎦ ⎣⎝ r ⎠

(1d)

.Here λ = 1 with σsw = 0.3284 nm and εsw = 0.4837 kJ/mol corresponds to the typical buckyball carbon−water parameters (atom type CA in AMBER39) used in previous simulations.42 We used λ = 0 to describe interactions of hydrophobic atoms, and λ = 2 to describe interactions of hydrophilic atoms with water oxygens. Each solute (B1/2, B1/4, and Bpatches) was solvated in approximately 2000 SPC/E43 water molecules in a 3D cubic 13067



Article

from MD simulations. We view density as a function of r, θ, and φ in a spherical coordinate system placed at the center of the reference solute. In this coordinate frame, density in a given thin spherical shell located at a distance, r, is a function of θ and φ and can be expressed in terms of orthogonal basis functions in θ and φ space:

periodic box (∼4 nm × 4 nm × 4 nm). Molecular dynamics (MD) simulations were performed in the isothermal−isobaric (N, P, T) ensemble using the GROMACS 4.5.3.44−46 Temperature (300 K) and pressure (1 bar) were maintained using the Nose−Hoover thermostat47,48 and Parrinello− Rahman barostat,49 respectively. A cutoff of 1 nm was used for all LJ interactions. The particle mesh Ewald (PME) method50 was used to calculate electrostatic interactions. A time step of 2 fs was used for all MD simulations. Each simulation run was for 1 ns. For comparison, we performed similar simulations for two chemically homogeneous solutes: purely hydrophobic and purely hydrophilic C180 buckyballs. 2.2. MD Simulations of Ubiquitin in Aqueous Solutions of Ligands. We performed all-atom MD simulations of the protein ubiquitin solvated in dilute aqueous solutions of three different ligands: (a) benzene, (b) guanidinium ion, and (c) Capto MMC (a multimodal ligand). The structure of each of these ligand molecules is shown in Figure 5. The initial structure for ubiquitin was taken from the RCSB Protein Data Bank (PDB ID: 1D3Z). The protonation states of ionizable residues were assigned at a pH of 5that is, all ionizable carboxyl (−1e), ammonium (+1e), and histidine (+1e) residues are charged. At this pH, ubiquitin has a charge of +1e. The force field parameters for ubiquitin were taken from AMBER parm-94.39 For guanidinium ions, OPLS model was used,51 which has parameters similar to the guanidinium group of arginine in AMBER and has been used in the literature for GdmCl simulations.52−54 For benzene and the MMC ligand, the LEAP program was used to generate their topology files,55 using parameters taken from similar groups in the AMBER parm-94 force field. For each simulation, one copy of ubiquitin was dissolved with multiple copies of the ligand and explicit water molecules in a ∼7 nm × 7 nm × 7 nm 3D periodic box. Sodium and chloride counterions were added to maintain electroneutrality, wherever necessary. The details on numbers of solutes and water molecules are given in Table 1. MD simulations were performed with the same runtime parameters as mentioned previously for the C180 solute systems.

N

f (θ , φ ) =

n=0 m=0

Capto MMC benzene guanidinium ion

no. of ligands 24 70 150

no. of waters 10636 10225 10058

no. of ions 23 1 151

(2a)

where Cnm and Snm are coefficients of the expansion. Two sets of orthogonal functions are combined in this expansion to capture the dependence on the polar (θ) and the azimuthal (φ) angles−the associated Legendre polynomials are used to describe the θ dependence, and sine and cosine functions capture the dependence on φ. Thus, e Y nm (θ , φ) = Pnm(cos θ) cos(mφ)

(2b)

o Y nm (θ , φ) = Pnm(cos θ) sin(mφ)

(2c)

where the superscripts e and o indicate the even and odd harmonics related to the cosine and sine functions, respectively, and are really the real and complex components of the general spherical harmonic function: Ynm(θ , φ) = ( −1)m

2n + 1 (n − m)! Pnm(cos θ )eimφ 4π (n + m)! (2d)

.Spherical harmonic expansion has been used in a number of applications in the past, for example, to describe the heat flow field around the earth,33 in spectroscopy,29,31 in computer graphics,35,36 and in medical imaging.37 To describe the r dependence of density, we make the coefficients of the expansion to vary with r. Hence, we write the three-dimensional density profile, ρ(r), of the molecules of interest as an expansion in spherical harmonics truncated at the Nth order,38,56 ρ (r ) = ρ0

Table 1. System Details for Protein Simulations ligand

n

e o (θ , φ) + SnmY nm (θ , φ)] ∑ ∑ [CnmY nm

N

n

e o (θ , φ) + Snm(r ) Y nm (θ , φ)] ∑ ∑ [Cnm(r) Y nm n=0 m=0

(2e)

time (ns)

.Using

200 50 50

Natoms

ρ(r) =

∑

δ(r − ri)

i=1

To obtain sufficient sampling of the space by ligands around the protein, we performed 50 ns long simulations each for ubiquitin in solutions of benzene and guanidinium ions, and 200 ns long simulations for ubiquitin in Capto MMC ligand solution. The Capto MMC ligand is larger, has lower diffusivity, and therefore takes longer to equilibrate. In simulations of benzene and guanidinium, we used first 5 ns as equilibration time, whereas a much longer window of 50 ns was used for equilibrating the Capto MMC solution. The remaining trajectories (with configurations saved every 1 ps) were used to perform analysis of spherical harmonic coefficients as described in later text. 2.3. Calculation of Density Distributions of Water and Ligands Using Spherical Harmonic Expansion Approach. We calculate the three-dimensional density profiles of the molecules of interest (e.g., water, ligands) around a reference solute (e.g., a C180 solute or a protein) using data

(2f)

and the orthonormality relation 2π

∫0 ∫0

π

e or o [Y nm (θ , φ)]2 sin θ dθ dφ =

(n + m)! 4π (1 ± δm0) 2(2n + 1) (n − m)!

(2g)

we obtain the coefficients Cnm(r) and Snm(r) for atoms of interest from simulation trajectories. Specifically, each simulation configuration provides vector locations ri for atoms of interest relative to the center of the reference solute. Angles θ and φ, and |ri| (distance from reference solute center) are used to obtain the contribution of that configuration to the ensemble average of coefficients. We perform block averaging to obtain error estimates for coefficients. Once the coefficients are available, the desired density at any location around the reference solute can be calculated using the expansion in eq 2e. As discussed below, the dominant features of the density profile are captured in the lowest order coefficients, and increasing 13068



Article

Figure 2. (a) C180 solutes, (b) water density profile at a radial distance of 9 Å from the center of the C180 solutes (red = high; blue = low), (c) spherical harmonic coefficients for hydrophilic atom distribution of C180 solutes, and (d) spherical harmonic coefficients for water around the C180 solutes. We selected overlapping bins of 1 Å width, centered every 0.1 Å radially.

them are fully captured by the spherically symmetric measures, such as g(r) or C00(r). However, rdfs cannot account for threedimensional density distributions around solutes with nonspherical hydrophobicity patterns on their surface. Indeed, rdf profiles in Figure 1c are unable to distinguish between the mixed solutes (especially B1/2 and Bpatches) which not only have different amounts of hydrophobicity/hydrophilicity but also have different chemical patterns on their surface. To quantify the differences in water density profiles around these solutes, we employ the spherical harmonic expansion approach. Three-Dimensional Water Density Profile around C180 Solutes Using a Spherical Harmonic Expansion. Figure 2 summarizes characteristics of the three-dimensional density profile of water oxygens surrounding three C180 solutes with nonspherical chemical patterns, B1/2, B1/4, and Bpatches. As discussed in Methods, we used MD simulation trajectories to calculate the distance-dependent coefficients of the spherical harmonic expansion, eq 2e. We note that the values of coefficients do not depend on the order of truncation. However, higher order coefficients capture increasingly complex angular dependence and contain larger errors. We calculated the coefficients up to fifth-order, similar to that used by Garde et al.38 Note that once coefficients up to Nth order are available, calculations of density can be done using any truncation up to the Nth order. E.g., truncating at N = 0 gives the spherically averaged C00(r) as the density profile. Figure 2a shows snapshots of three C180 solutes that highlight the hydrophilicity patterns on their surface. Figure 2d shows the corresponding selected spherical harmonic coefficients of water density up to the third order. A more complete picture of all coefficients is available in the Supporting Information. These coefficients show oscillations (about zero) with a period consistent with layering of water surrounding the solutes, and approach zero beyond ∼1 nm from the solute

complexity, especially the angular dependence, is incorporated by including the higher order terms in the expansion.

3. RESULTS AND DISCUSSION 3.1. Describing the Hydration of Various C180 Solutes. Comparing Radial Distribution Functions. We studied the hydration of five different versions of C180 buckyball solutes. Although these solutes are topographically spherically symmetric, only the fully hydrophobic and hydrophilic versions have spherically symmetric chemistry. Thus, we expect the hydration of partially hydrophilic versions, B1/2, B1/4, and Bpatches (Figure 1b) to display nonsymmetric hydration patterns. Figure 1c shows C180−water−oxygen radial distribution functions (rdfs), which are the same as the coefficient C00(r) in the spherical harmonic expansion. Water layers near all five solutes, as reflected in the well-defined oscillations in the rdf profiles. The local density of water near the fully hydrophobic C180 solute is the lowest, consistent with the expected dewetting of such a solute.57−62 Note that the size of C180 solute is in the “crossover region”.62,63 Further increasing the solutes size will reduce the vicinal density of water,57 eventually leading to a sigmoidal, vapor−liquid interface-like profile near much larger solutes. With increasing hydrophilicity, the local density of water increases, and the first hydration shell also becomes more welldefined. Specifically, the first peak of rdf increases from hydrophobic to B1/4, to B1/2 (∼Bpatches), and reaches ∼2.5 for the fully hydrophilic version of C180. A quantitative study of the dependence of hydration on the strength of solute−water attractions has been presented previously by Athawale et al.42 Our results for purely hydrophobic and hydrophilic versions in Figure 1c are in excellent agreement with that study. The purely hydrophobic and purely hydrophilic buckyballs are spherically symmetric, and water density profiles around 13069



Article

surface. That these coefficients oscillate about zero is consistent with the fact that C00(r) captures the spherically averaged density in a shell, whereas the other coefficients represent perturbations from that average representing the nonspherical nature of the true density profile. The hydrophobicity/-philicity pattern on the C180 surface is expected to influence the water density profile surrounding it. To explore whether the symmetries of that underlying pattern are also reflected in the water density profile, we also calculated the spherical harmonic coefficients of the hydrophilic chemical pattern on the solute surface, as shown in Figure 2c. The B1/2 solute exhibits symmetrical hydrophilic atom distribution about its equator. This distribution can be solely captured by harmonics which are circularly symmetric, or the zonal harmonics with m = 0. In such cases, the even harmonics, Yen0(θ, φ), reduce to associated Legendre polynomials without any dependence on the azimuthal angle, φ: Y ne0(θ , φ) = Pn0(cos θ )

(3a)

and the odd harmonics, Yon0(θ, φ), disappear. Indeed, we find that all the even coefficients with m = 0 are nonzero. All other coefficients, which describe nonsymmetric distributions about the equator, do not contribute to the hydrophilic atom distribution in B1/2. Since the water density profile around these solutes is expected to mimic the underlying hydrophilic atom distribution, the coefficients for water density distribution around B1/2 also show significant values only if they correspond to m = 0. The patterns for B1/4 and Bpatches solutes are more complex relative to the B1/2 solute. As such, it is difficult to deconvolute the contributions of individual coefficients in these cases. Still, we observe that for each solute, similar coefficients exhibit dominance in the hydrophilic atom distribution and the corresponding water density distribution. We used the spherical harmonic coefficients up to N = 5 to calculate the spatially resolved water density around the three solutes. This water density distribution in the first hydration shell of the solute is shown in Figure 2b using a blue (low) to red (high) color map. For B1/2 and B1/4 solutes, the underlying hydrophilicity patterns are clearly visible in the water density map. Near the Bpatches solute, an additional interesting feature the asymmetry of hydrophobic and hydrophilic solvationfirst pointed out by Acharya et al.,64 also becomes apparent in the water density distribution. Specifically, introduction of even a small hydrophilic patch in a hydrophobic background is able to enhance the local water density significantly. The ability to efficiently calculate three-dimensional density profiles of atomic sites of molecules can also provide some insights into orientations of those molecules in the vicinity of the reference solute. This can be especially useful if the density profile is not spherically symmetric. Figure 3 demonstrates this idea for water as a ligand by resolving its oxygen (OW) and hydrogen (HW) density profiles around three versions of C180 solutes. The features of spherically averaged density profiles of OW and HW profiles in Figure 3a are similar near the three solutes; namely, the first peak of HW and OW densities is approximately at the same distance, suggesting that water molecules adopt orientations in which the HOH plane of water is parallel (or tangential) to the solute surface. The tail of the hydrogen density profile at low solute−hydrogen separations suggests a small tendency of water molecules to point their hydrogen atom(s) toward C180 solutes, which is a well-known feature of hydrophobic hydration.65−67

Figure 3. (a) Solute−oxygen and solute−hydrogen rdfs (i.e., C00 coefficients) around the three C180 solutes. (b) ρ/ρ0 values for water oxygen calculated using the spherical harmonics expansion in the plane that bisects the solute as shown (red = high density; blue = low density). (c) Local charge density (qOW(ρ/ρ0,OW) + 2qHW(ρ/ρ0,HW)) calculated using the spherical harmonics expansion in the same plane as in b . Red and blue colors indicate negative and positive charge densities, respectively.

Figure 3b shows the oxygen density distribution in a plane that bisects the solute vertically (as shown), and Figure 3c shows the partial charge density, (qOW(ρ/ρ0,OW) + 2qHW(ρ/ ρ0,HW)), in the same plane. The enhanced layering of water molecules near hydrophilic regions is clear from these distributions. Interestingly, even near hydrophilic patches (e.g., in the Janus B1/2 solute or in other B1/4 and Bpatches solutes), the underlying orientational preferences of water remain largely unaffected. This observation suggests that understanding the orientational preferences of water near hydrophilic solutes that have high van der Waals interactions (high λ in our scheme) but no directional electrostatic interactions (e.g., hydrogen bonding) with water is an interesting problem worth pursuing in an independent study in the future. 3.2. Solvation of Ubiquitin in Aqueous Solutions of Ligands. We applied the spherical harmonics expansion approach to calculate ligand density profiles around the protein ubiquitin using simulations of ubiquitin in aqueous solutions of three different ligands (see Table 1). Ubiquitin is a small 76 amino acid protein that is ubiquitously expressed in eukaryotes.68 High-resolution structures of wild-type ubiquitin are available.69 One face of ubiquitin has a cluster of several hydrophobic residues: L8, I44, V70, and L73, flanked by a number of positively charged residues: K6, R42, K48, H68, R72, and R74 (Figure 4). Thus, multimodal ligand molecules, which contain both hydrophobic and negatively charged groups, would be expected to interact favorably with this face of the protein. To understand the molecular details of such interactions, we selected three ligands with varying degrees of chemical complexity, (Figure 5). Figure 5a shows the active ligand of the chromatographic resin “Capto MMC”, which is known to exhibit high selectivity in multimodal chromato13070



Article

Figure 4. (a) Cartoon representation of ubiquitin (showing the secondary structural elements). (b, c) Charged (blue = positive; red = negative) and hydrophobic (green) residues on the front and back face of ubiquitin.

We also calculated density distributions for two other molecules, benzene and guanidinium ions (Figure 5b,c) around ubiquitin. Benzene is a planar, neutral, hydrophobic molecule whereas guanidinium is a planar, bulky, positively charged ion. As such, benzene is expected to interact primarily with hydrophobic residues of ubiquitin, probably near the hydrophobic cluster mentioned earlier. While guanidine has been shown to interact strongly with hydrophobic polymers, driven by van der Waals interactions,54 the hydrophobic cluster on the face of ubiquitin is surrounded by positively charged residues, which may repel guanidinium ions. For ubiquitin, one might expect guanidinium ions to interact primarily with the negatively charged residues on the back face of the protein. Our application of the spherical harmonics expansion to calculate the spatially resolved density of these ligands (and their constituents in the case of Capto MMC) provides insights into the nature of protein−ligand interactions. Like that for C180 solutes, we used 1 Å wide bins placed radially from the protein center of geometry. The bin centers were placed every 0.1 Å, leading to overlapping bins, and a fine-grained coverage of the space surrounding the protein surface. Interactions of Capto MMC Ligand Molecules with Ubiquitin. Once the coefficients of the spherical harmonic expansion are obtained as a function of the distance from the protein center, the density of a ligand can be calculated at any point in space surrounding the protein. For physically meaningful results, it is best to define a surface that envelopes

Figure 5. Small molecule ligands as probes for interactions with ubiquitin: (a) Capto MMC ligand, (b) benzene, and (c) guanidinium ion (green = carbon; red = oxygen; blue = nitrogen; yellow = sulfur; white = hydrogen).

graphic systems. This Capto MMC ligand has multiple chemical moieties, each with a different mode of interaction: benzene (hydrophobic), amide (hydrogen bond), carboxylate (electrostatic), and aliphatic (hydrophobic). Previous studies using chromatographic and NMR experiments and MD simulations have highlighted the importance of interactions of this ligand with the front face of ubiquitin.8,28

Figure 6. Effect of increasing order of spherical harmonic coefficients on the density of the Capto MMC ligand (center of geometry) near the protein surface. With increasing order from N = 1 to N = 5, higher spatial resolution is achieved. Red indicates regions of high ligand density (ρ/ρ0 ∼ 200) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2. 13071



Article

the protein and represents the ligand density profile on this surface using a convenient color map. Such a surface could be defined using a variety of methods, e.g., by selecting a locus of points whose proximal distance to the protein is the same. Here we use the method of Willard and Chandler70 to define the protein surface envelope as a locus of equidensity points having a small fixed value of coarse-grained density of protein atoms. The coarse-graining procedure of Willard and Chandler replaces each heavy atom with a truncated Gaussian density distribution. This coarse-grained density field of protein atoms tails off away from the protein surface. Selecting points in space surrounding the protein which have the given same (small) value of coarse-grained density provides such an envelope. The smaller this value, the farther is the envelope from the protein surface. A coarse-grained value equal to half the bulk density in the protein interior, ρ′ = ρprotein/ρprotein,core = 0.5, would correspond to the “instantaneous surface” as defined by Willard and Chandler. This enables us to readily examine the behavior of ligands at different distances from the protein surface. Figure 6 shows the density for the center-of-geometry of the Capto MMC ligand near ubiquitin using the spherical harmonic expansion truncated at different levels from N = 1 to 5. In the lowest order description (N = 1), the ligand density appears quite homogeneous around the protein. As we include a higher orders harmonics, the spatial resolution of the density increases. Density distribution is very similar for N = 4 and N = 5, where we observe well-defined high-density regions around the protein for N = 5. To maintain such spatial resolution, we chose N = 5 for further analysis. Figure 7 shows Capto MMC density calculated for ρ′ = 0.2, 0.1, and 0.02 interfaces around ubiquitin, capturing the ligand

Figure 8. ρ/ρ0 profile for the MMC ligand in the local domain of ubiquitin: (a) front, and (b) back. Red indicates regions of high ligand density (ρ/ρ0 ∼ 200) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2.

We also calculated density profiles for different subgroups present within the MMC ligand molecule to characterize the different modes of the multimodal ligand−protein interactions, as shown in Figure 9. The profiles for the benzene, carboxylate,

Figure 9. ρ/ρ0 values plotted for different subgroups of the MMC ligand in the local domain of ubiquitin (front face): (a) benzene group, (b) amide group, (c) carboxylate group, and (d) aliphatic group. Red indicates regions of high ligand density (ρ/ρ0 ∼ 200) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2. (e) Grid points with ρ/ρ0 > 100 for specific modes. The grid points for each mode are colored differently: green = benzene; blue = amide; red = carboxylate; yellow = aliphatic. Inset: Possible binding conformation of the ligand in the high-density site.

Figure 7. Ligand density values at different protein density interfaces.

density in regions very close to the protein and moving out toward the bulk region. As expected, ligand density is most affected in the local region, leading to well-defined density patterns that become coarse-grained and more homogeneous as the surface moves away from the protein. Figure 8 shows the density profile of Capto MMC near the front and back face of ubiquitin. The front face of the protein (Figure 8a) is in the same orientation as shown in Figure 4. The Capto MMC ligand displays a high-density hot-spot on the front face, surrounded by a larger region of weak binding, which coincides with the hydrophobic cluster surrounded by positively charged residues (Figure 4b). There are two more regions on the back face of ubiquitin, where the ligand binds weakly. This is consistent with the idea that complex, multimodal ligands, such as Capto MMC, can have multiple interaction sites of varying strengths around a protein. Our spherical harmonic expansion based approach is able to characterize these sites based on their location and strengths.

and amide moieties contain well-defined high-density subregions in the main binding region, whereas the high-density region for the aliphatic group is relatively less localized. The carboxylate group mainly interacts with the positively charged residues K6 and H68, and benzene and aliphatic groups interact with the hydrophobic patch consisting of residues L8, V70, and I44 (see Figure 4b). Comparison of the centers of high density for the different modes suggest a preferred orientation of the Capto MMC ligand in the binding region, with the benzene and aliphatic group oriented toward the center of the hydrophobic cluster and the charged group extending toward the charged region made up of residues K6 and H68 (Figure 9e). Density values plotted at another protein interface 13072



Article

(Supporting Information, Figure S4) indicate that the ligand can also adopt a second orientation in which the carboxylate group interacts with the residue R42. The primary high-density site of the Capto MMC ligand on the front face of ubiquitin results from the favorable interactions of all of the subgroups of the ligand with protein, as indicated by the high-density regions for those groups near the binding site (Figure 9). This may suggest synergistic interactions of the multimodal Capto MMC ligand with the chemically heterogeneous surface of ubiquitin. In contrast, on the back face of ubiquitin (Figure 10), interactions of almost

Figure 11. ρ/ρ0 values plotted in the local domain of ubiquitin: (a) benzene subgroup of Capto MMC ligand, and (b) benzene molecule. Red indicates regions of high ligand density (Capto MMC, ρ/ρ0 ∼ 200; benzene molecule, ρ/ρ0 ∼ 16) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2.

effect on protein selectivity in multimodal cation exchange chromatographic systems.23 Accordingly, it was of interest to examine protein−guanidinium chloride simulations and the results are presented in Figure 12. There are three high-density

Figure 10. ρ/ρ0 values plotted for different subgroups of the MMC ligand in the local domain of ubiquitin (back face): (a) benzene group, (b) amide group, (c) carboxylate group, and (d) aliphatic group. Red indicates regions of high ligand density (ρ/ρ0 ∼ 200) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2.

every group of the ligand are weaker, with the carboxylate group being the least involved, consistent with the weaker overall binding of the ligand on that face.8 Interactions of Benzene Molecules with Ubiquitin. We are interested in understanding the differences between the binding of free benzene versus the binding of the benzene moiety in the Capto MMC ligand to ubiquitin. Figure 11b shows the highdensity regions for interactions of benzene with ubiquitin, highlighting that free benzene molecules interact primarily with I44 and H68. In addition, benzene also interacts, albeit weakly, with F45. In contrast, the binding of the benzene moiety in the MMC ligand is limited to a narrow region on the protein surface consisting primarily of L8 and V70, likely due to the additional constraint of the other interactions by the MM ligand. These results clearly indicate the effect of chemical and topographical context of the ligand for determining its interaction with a protein surface. Interactions of Guanidinium Ions with Ubiquitin. Guanidine has been previously shown to have a significant

Figure 12. ρ/ρ0 values plotted for guanidinium ions in the local domain of ubiquitin: (a) front and (b) back. Inset shows the positive (blue) and negative (red) residues on ubiquitin. Red indicates regions of high ligand density (ρ/ρ0 ∼ 10.5) while blue indicates regions of low ligand density (ρ/ρ0 ∼ 1). Values are plotted at the protein density interface with ρ′ = 0.2.

regions for guanidinium ions on the back face of ubiquitin. A look at the location of negatively charged residues on ubiquitin (see inset) reveals that while guanidinium interactions with ubiquitin are strongly affected by electrostatics, the locations of high-density sites for guanidinium do not always precisely match with the locations of negatively charged residues. Sites 2 and 3 are centered at neutral regions flanked by negatively charged residues (E21/E24 and D32/E34, respectively). On the other hand, at the bottom of the back face of ubiquitin, 13073



Article

there are multiple negative and positive charges (E16, E18, and K29) which seem to result in minimal interaction of guanidinium ions with this particular region. Interestingly, while previous studies have shown association of guanidinium ions to hydrophobic polymers through significant van der Waals interactions,54 purely hydrophobic association of guanidinium ions to the protein surface is not observed in these simulations. This is likely due to the significant electrostatic interactions with negatively charged residues coupled with the lack of contiguous hydrophobic patches on ubiquitin’s surface (without the presence of positively charged residues nearby to repel the ions) and the relatively low concentration of guanidinium ions used in these simulations.

becoming available from MD simulations on faster computing platforms.

■

ASSOCIATED CONTENT

S Supporting Information *

Spherical harmonic coefficients for the hydrophilic atom distributions in B1/2, B1/4, and Bpatches in Figures S1, S2, and S3, respectively, and density profile for the carboxylate group of the Capto MMC ligand at different protein density interfaces around ubiquitin in Figure S4. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel.: (518) 276-6298.

4. CONCLUSION Performing the analysis of atomic density profiles in a spatially resolved manner in the vicinity of complex solutes that contain topographical and chemical heterogeneity can present challenges in MD simulations. We employed a method to calculate the spatial density distribution of ligands using an expansion in terms of spherical harmonics. We applied the method to study the hydration of three different versions of a buckyball C180 solute containing hydrophobic/hydrophilic patterns on the surface. Atomic and partial charge density profiles around these buckyball solutes capture the physics of hydration, e.g., the layering of water molecules near the hydrophilic regions of the C180 solutes, water orientations, and response of water density to the underlying pattern of hydrophobicity. We also performed simulations of a protein ubiquitin surrounded by aqueous solutions of three different ligands: Capto MMC, benzene, and guanidine. Application of the spherical harmonic expansion provides the overall distribution of ligand centers around the protein as well as the distributions of various ligand sites relative to each other revealing information about preferred ligand binding orientations and regions of strong and weak binding. The Capto MMC ligand displays a preference for the front face of ubiquitin, which corroborates results from previous studies.8,28 The collective information obtained from binding maps of different ligands to the protein surface obtained here is akin to that obtained by multiple crystal solvent structure method of Ringe et al.71 Understanding of binding of specific ligands to protein surfaces is especially useful in bioseparations, drug design and biomaterials development. Our results on the differential bindings of the front and back faces of ubiquitin suggest a preferred orientation of ubiquitin as it interacts with ligands immobilized on a chromatographic resin. This preference can be manipulated either by mutating the protein or by adding cosolvent molecules which compete with ligand−protein interactions. This may be useful in practical applications, such as the separation of close protein variants and the design of high-affinity ligands. Further, a multimodal ligand (like Capto MMC) can be immobilized to expose its different chemical moieties for interactions with the protein, providing an additional handle on protein−ligand interaction. Calculation of distance-dependent spherical harmonic coefficients represents a contraction of a large amount of MD trajectory data into smaller, yet physically meaningful, representation. This factor and the ability to easily visualize the three-dimensional information will become increasingly important in light of the larger amount of trajectory data

Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We gratefully acknowledge financial support from National Science Foundation Grant CBET-1160039. We also thank Srivathsan Vembanur, Hari Acharya, and Vasudevan Venkateshwaran for numerous useful discussions.

■

REFERENCES

(1) Kitchen, D. B.; Decornez, H.; Furr, J. R.; Bajorath, J. Docking and scoring in virtual screening for drug discovery: Methods and applications. Nat. Rev. Drug Discovery 2004, 3 (11), 935−949. (2) Rees, D. C.; Congreve, M.; Murray, C. W.; Carr, R. Fragmentbased lead discovery. Nat. Rev. Drug Discovery 2004, 3 (8), 660−672. (3) Gilson, M. K.; Zhou, H. X. Calculation of protein-ligand binding affinities. Annu. Rev. Biophys. Biomol. Struct. 2007, 36 (1), 21−42. (4) Lu, Y.; Wang, Y.; Zhu, W. Nonbonding interactions of organic halogens in biological systems: Implications for drug discovery and biomolecular design. Phys. Chem. Chem. Phys. 2010, 12 (18), 4543− 4551. (5) Kuhn, B.; Mohr, P.; Stahl, M. Intramolecular hydrogen bonding in medicinal chemistry. J. Med. Chem. 2010, 53 (6), 2601−2611. (6) Reynolds, C. H.; Holloway, M. K. Thermodynamics of Ligand Binding and Efficiency. ACS Med. Chem. Lett. 2011, 2 (6), 433−437. (7) Melander, W. R.; el Rassi, Z.; Horvath, C. Interplay of hydrophobic and electrostatic interactions in biopolymer chromatography. Effect of salts on the retention of proteins. J. Chromatogr. A 1989, 469, 3−27. (8) Chung, W. K.; Freed, A. S.; Holstein, M. A.; McCallum, S. A.; Cramer, S. M. Evaluation of protein adsorption and preferred binding regions in multimodal chromatography using NMR. Proc. Natl. Acad. Sci. U. S. A. 2010, 107 (39), 16811−16816. (9) Holstein, M. A.; Chung, W. K.; Parimal, S.; Freed, A. S.; Barquera, B.; McCallum, S. A.; Cramer, S. M. Probing multimodal ligand binding regions on ubiquitin using nuclear magnetic resonance, chromatography, and molecular dynamics simulations. J. Chromatogr. A 2012, 1229, 113−120. (10) Geschwindner, S.; Carlsson, J. F.; Knecht, W. Application of optical biosensors in small-molecule screening activities. Sensors 2012, 12 (4), 4311−4323. (11) Famulok, M.; Mayer, G. Aptamer modules as sensors and detectors. Acc. Chem. Res. 2011, 44 (12), 1349−1358. (12) Gray, J. J. The interaction of proteins with solid surfaces. Curr. Opin. Struct. Biol. 2004, 14 (1), 110−115. (13) Harner, M. J.; Frank, A. O.; Fesik, S. W. Fragment-based drug discovery using NMR spectroscopy. J. Biomol. NMR 2013, 56 (2), 65− 75. (14) Vuignier, K.; Schappler, J.; Veuthey, J. L.; Carrupt, P. A.; Martel, S. Drug-protein binding: A critical review of analytical tools. Anal. Bioanal. Chem. 2010, 398 (1), 53−66.

13074



Article

(15) Renaud, J. P.; Delsuc, M. A. Biophysical techniques for ligand screening and drug design. Curr. Opin. Pharmacol. 2009, 9 (5), 622− 628. (16) Homans, S. W. NMR spectroscopy tools for structure-aided drug design. Angew. Chem., Int. Ed. 2004, 43 (3), 290−300. (17) Roberts, G. C. NMR spectroscopy in structure-based drug design. Curr. Opin. Biotechnol. 1999, 10 (1), 42−47. (18) Friddle, R. W.; Battle, K.; Trubetskoy, V.; Tao, J.; Salter, E. A.; Moradian-Oldak, J.; De Yoreo, J. J.; Wierzbicki, A. Single-molecule determination of the face-specific adsorption of Amelogenin’s Cterminus on hydroxyapatite. Angew. Chem., Int. Ed. 2011, 50 (33), 7541−7545. (19) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications (Computational Science), 2nd ed.; Academic Press: San Diego, CA, USA, 2001. (20) Jorgensen, W. L. The Many Roles of Computation in Drug Discovery. Science 2004, 303 (5665), 1813−1818. (21) Canchi, D. R.; Jayasimha, P.; Rau, D. C.; Makhatadze, G. I.; Garcia, A. E. Molecular mechanism for the preferential exclusion of TMAO from protein surfaces. J. Phys. Chem. B 2012, 116 (40), 12095−12104. (22) Canchi, D. R.; Paschek, D.; Garcia, A. E. Equilibrium study of protein denaturation by urea. J. Am. Chem. Soc. 2010, 132 (7), 2338− 2344. (23) Holstein, M. A.; Parimal, S.; McCallum, S. A.; Cramer, S. M. Mobile phase modifier effects in multimodal cation exchange chromatography. Biotechnol. Bioeng. 2012, 109 (1), 176−186. (24) Kang, M.; Smith, P. E. Preferential interaction parameters in biological systems by Kirkwood−Buff theory and computer simulation. Fluid Phase Equilib. 2007, 256 (1−2), 14−19. (25) Baynes, B. M.; Trout, B. L. Proteins in Mixed Solvents: A Molecular-Level Perspective. J. Phys. Chem. B 2003, 107 (50), 14058− 14067. (26) Shukla, D.; Shinde, C.; Trout, B. L. Molecular Computations of Preferential Interaction Coefficients of Proteins. J. Phys. Chem. B 2009, 113 (37), 12546−12554. (27) Morrison, C. J.; Godawat, R.; McCallum, S. A.; Garde, S.; Cramer, S. M. Mechanistic studies of displacer-protein binding in chemically selective displacement systems using NMR and MD simulations. Biotechnol. Bioeng. 2009, 102 (5), 1428−1437. (28) Freed, A. S.; Garde, S.; Cramer, S. M. Molecular simulations of multimodal ligand-protein binding: elucidation of binding sites and correlation with experiments. J. Phys. Chem. B 2011, 115 (45), 13320− 13327. (29) Nomura, S.; Kawai, H.; Kimura, I.; Kagiyama, M. General description of orientation factors in terms of expansion of orientation distribution function in a series of spherical harmonics. J. Polym. Sci., Part A-2: Polym. Phys. 1970, 8 (3), 383−400. (30) Pottel, H.; Herreman, W.; van der Meer, B. W.; Ameloot, M. On the significance of the fourth-rank orientational order parameter of fluorophores in membranes. Chem. Phys. 1986, 102 (1−2), 37−44. (31) Sourisseau, C. Polarization Measurements in Macro- and MicroRaman Spectroscopies: Molecular Orientations in Thin Films and Azo-Dye Containing Polymer Systems. Chem. Rev. 2004, 104 (9), 3851−3892. (32) Rodriguez, V.; Lagugné-Labarthet, F.; Sourisseau, C. Orientation Distribution Functions Based upon Both upon ⟨P1⟩, ⟨P3⟩ Order Parameters and upon the Four ⟨P1⟩ up to ⟨P4⟩ Values: Application to an Electrically Poled Nonlinear Optical Azopolymer Film. Appl. Spectrosc. 2005, 59 (3), 322−328. (33) Pollack, H. N.; Hurter, S. J.; Johnson, J. R. Heat flow from the Earth’s interior: Analysis of the global data set. Rev. Geophys. 1993, 31 (3), 267−280. (34) Novák, P. High Resolution Constituents of the Earth’s Gravitational Field. Surv. Geophys. 2010, 31 (1), 1−21. (35) Westin, S. H.; Arvo, J. R.; Torrance, K. E. Predicting reflectance functions from complex surfaces. SIGGRAPH Comput. Graph. 1992, 26 (2), 255−264.

(36) Ramamoorthi, R.; Hanrahan, P. A signal-processing framework for reflection. ACM Trans. Graph. 2004, 23 (4), 1004−1042. (37) Alexander, D. C.; Barker, G. J.; Arridge, S. R. Detection and modeling of non-Gaussian apparent diffusion coefficient profiles in human brain data. Magn. Reson. Med. 2002, 48 (2), 331−340. (38) Garde, S.; Hummer, G.; Garcia, A. E.; Pratt, L. R.; Paulaitis, M. E. Hydrophobic hydration: Inhomogeneous water structure near nonpolar molecular solutes. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 1996, 53 (5), R4310−R4313. (39) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. A Second Generation Force Field for the Simulation of Proteins, Nucleic Acids, and Organic Molecules. J. Am. Chem. Soc. 1995, 117 (19), 5179−5197. (40) Weeks, J. D. Role of Repulsive Forces in Determining the Equilibrium Structure of Simple Liquids. J. Chem. Phys. 1971, 54 (12), 5237−5246. (41) Chandler, D.; Weeks, J. D.; Andersen, H. C. Van der Waals picture of liquids, solids, and phase transformations. Science 1983, 220 (4599), 787−794. (42) Athawale, M. V.; Jamadagni, S. N.; Garde, S. How hydrophobic hydration responds to solute size and attractions: Theory and simulations. J. Chem. Phys. 2009, 131 (11), No. 115102. (43) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91 (24), 6269− 6271. (44) Berendsen, H. J. C.; van der Spoel, D.; van Drunen, R. GROMACS: A message-passing parallel molecular dynamics implementation. Comput. Phys. Commun. 1995, 91 (1−3), 43−56. (45) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4 (3), 435−447. (46) Pronk, S.; Pall, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; Shirts, M. R.; Smith, J. C.; Kasson, P. M.; van der Spoel, D.; Hess, B.; Lindahl, E. GROMACS 4.5: A high-throughput and highly parallel open source molecular simulation toolkit. Bioinformatics 2013, 29 (7), 845−854. (47) Evans, D. J.; Holian, B. L. The Nose−Hoover thermostat. J. Chem. Phys. 1985, 83 (8), No. 4069. (48) Nosé, S. A unified formulation of the constant temperature molecular dynamics methods. J. Chem. Phys. 1984, 81 (1), 511−519. (49) Parrinello, M.; Rahman, A. Crystal Structure and Pair Potentials: A Molecular-Dynamics Study. Phys. Rev. Lett. 1980, 45 (14), 1196− 1199. (50) Ewald, P. P. Evaluation of optical and electrostatic lattice potentials. Ann. Phys. (Leipzig, Ger.) 1921, 64, 253−287. (51) Jorgensen, W. L.; Tirado-Rives, J. The OPLS [optimized potentials for liquid simulations] potential functions for proteins, energy minimizations for crystals of cyclic peptides and crambin. J. Am. Chem. Soc. 1988, 110 (6), 1657−1666. (52) Mason, P. E.; Neilson, G. W.; Enderby, J. E.; Saboungi, M. L.; Dempsey, C. E.; MacKerell, A. D., Jr.; Brady, J. W. The structure of aqueous guanidinium chloride solutions. J. Am. Chem. Soc. 2004, 126 (37), 11462−11470. (53) O’Brien, E. P.; Dima, R. I.; Brooks, B.; Thirumalai, D. Interactions between hydrophobic and ionic solutes in aqueous guanidinium chloride and urea solutions: Lessons for protein denaturation mechanism. J. Am. Chem. Soc. 2007, 129 (23), 7346− 7353. (54) Godawat, R.; Jamadagni, S. N.; Garde, S. Unfolding of hydrophobic polymers in guanidinium chloride solutions. J. Phys. Chem. B 2010, 114 (6), 2246−2254. (55) Case, D. A.; Darden, T. A.; Cheatham, T. E., III; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Merz, K. M.; Wang, B.; Pearlman, D. A.; Crowley, M.; Brozell, S.; Tsui, V.; Gohlke, H., Mongan, J.; Hornak, V.; Cui, G.; Beroza, P.; Schafmeister, C.; Caldwell, J. W.; Ross, W. S.; Kollman, P. A. AMBER 8; University of California: San Diego, 2004. 13075



Article

(56) Arfken, G. Mathematical Methods for Physicists, 3rd ed.; Academic Press: San Diego, CA, 1985. (57) Stillinger, F. H. Structure in aqueous solutions of nonpolar solutes from the standpoint of scaled-particle theory. J. Solution Chem. 1973, 2 (2−3), 141−158. (58) Lum, K.; Chandler, D.; Weeks, J. D. Hydrophobicity at Small and Large Length Scales. J. Phys. Chem. B 1999, 103 (22), 4570−4577. (59) Hummer, G.; Garde, S. Cavity Expulsion and Weak Dewetting of Hydrophobic Solutes in Water. Phys. Rev. Lett. 1998, 80 (19), 4193−4196. (60) Ashbaugh, H. S.; Pratt, L. R. Colloquium: Scaled particle theory and the length scales of hydrophobicity. Rev. Mod. Phys. 2006, 78 (1), 159−178. (61) Huang, X.; Margulis, C. J.; Berne, B. J. Dewetting-induced collapse of hydrophobic particles. Proc. Natl. Acad. Sci. U. S. A. 2003, 100 (21), 11953−11958. (62) Rajamani, S.; Truskett, T. M.; Garde, S. Hydrophobic hydration from small to large lengthscales: Understanding and manipulating the crossover. Proc. Natl. Acad. Sci. U. S. A. 2005, 102 (27), 9475−9480. (63) Chandler, D. Interfaces and the driving force of hydrophobic assembly. Nature 2005, 437 (7059), 640−647. (64) Acharya, H.; Vembanur, S.; Jamadagni, S. N.; Garde, S. Mapping hydrophobicity at the nanoscale: Applications to heterogeneous surfaces and proteins. Faraday Discuss. 2010, 146, 353−365. (65) Lee, C.; McCammon, A.; Rossky, P. J. The structure of liquid water at an extended hydrophobic surface. J. Chem. Phys. 1984, 80 (9), 4448−4455. (66) Stirnemann, G.; Rossky, P. J.; Hynes, J. T.; Laage, D. Water reorientation, hydrogen-bond dynamics and 2D-IR spectroscopy next to an extended hydrophobic surface. Faraday Discuss. 2010, 146, 263− 281. (67) Trudeau, T. G.; Jena, K. C.; Hore, D. K. Water Structure at Solid Surfaces of Varying Hydrophobicity. J. Phys. Chem. C 2009, 113 (46), 20002−20008. (68) Joazeiro, C. A. P. BIOCHEMISTRY: UbiquitinationMore Than Two to Tango. Science 2000, 289 (5487), 2061−2062. (69) Cornilescu, G.; Marquardt, J. L.; Ottiger, M.; Bax, A. Validation of protein structure from anisotropic carbonyl chemical shifts in a dilute liquid crystalline phase. J. Am. Chem. Soc. 1998, 120 (27), 6836− 6837. (70) Willard, A. P.; Chandler, D. Instantaneous liquid interfaces. J. Phys. Chem. B 2010, 114 (5), 1954−1958. (71) Mattos, C.; Bellamacina, C. R.; Peisach, E.; Pereira, A.; Vitkup, D.; Petsko, G. A.; Ringe, D. Multiple solvent crystal structures: Probing binding sites, plasticity and hydration. J. Mol. Biol. 2006, 357 (5), 1471−82.

13076


Application of a Spherical Harmonics Expansion Approach for

Recommend Documents