Direct Determination of Site-Specific Noncovalent Interaction

Apr 28, 2017 - A novel approach to accurately determine residue-specific noncovalent interaction strengths (ξ) of proteins from NMR-measured fast sid...
1 downloads 0 Views 2MB Size
Subscriber access provided by University of Colorado Boulder

Article

Direct Determination of Site-specific Noncovalent Interaction Strengths of Proteins from NMR-derived Fast Side Chain Motional Parameters Rajitha Rajeshwar Tatikonda, and Marimuthu Krishnan J. Phys. Chem. B, Just Accepted Manuscript • Publication Date (Web): 28 Apr 2017 Downloaded from http://pubs.acs.org on April 28, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Direct Determination of Site-specific Noncovalent Interaction Strengths of Proteins from NMR-derived Fast Side Chain Motional Parameters Rajitha Rajeshwar T. and Marimuthu Krishnan∗ Center for Computational Natural Sciences and Bioinformatics, International Institute of Information Technology, Gachibowli, Hyderabad 500 032, India E-mail: [email protected] Phone: +91-40-6653-1447. Fax: +91-40-6653-1413

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract A novel approach to accurately determine residue-specific noncovalent interaction strengths (ξ) of proteins from NMR-measured fast side chain motional parameters 2 (Oaxis ) is presented. By probing the environmental sensitivity of side chain conforma-

tional energy surfaces of individual residues of a diverse set of proteins, the microscopic 2 connections between ξ, Oaxis , conformational entropy (Sconf ), conformational barriers

and rotamer stabilities established here are found to be universal among proteins. The results reveal that side chain flexibility and conformational entropy of each residue decrease with increasing ξ and that for each residue type there exists a critical range of ξ, determined primarily by the mean side chain conformational barriers, within which 2 flexibility of any residue can be reversibly tuned from highly flexible (with Oaxis ∼ 0) to 2 highly restricted (with Oaxis ∼ 1) by increasing ξ by ∼3 kcal/mol. Beyond this critical

range of ξ, both side chain flexibility and conformational entropy are insensitive to ξ. The interrelationships between conformational dynamics, conformational entropy, and noncovalent interactions of protein side chains established here open up new avenues to probe perturbation-induced (for example, ligand-binding, temperature, pressure) changes in fast side chain dynamics and thermodynamics of proteins by comparing their conformational energy surfaces in the native and perturbed states.

Introduction The noncovalent interactions play a decisive role in determining the function and stability of proteins. 1–8 Their importance to the thermodynamics and kinetics of protein folding, optimal packing of side chains in the native and non-native states of proteins, molecular recognition, binding and transport of ligands to and from protein active sites has long been recognized. 9–11 Nevertheless, a quantitative experimental characterization of how noncovalent interactions of individual residues influence folding, structure, dynamics, and function of proteins proved difficult. Although it is known that noncovalent interactions of individual residues of a

2

ACS Paragon Plus Environment

Page 2 of 51

Page 3 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

protein differ in strengths and chemical origins, an accurate measure of interaction strengths of individual residues and their contributions to protein stability are not readily amenable to experimental investigations. 9,12,13 Most experimental techniques aiming to probe interaction strengths of individual residues mainly rely on either chemically modified proteins (i.e., proteins subjected to covalent modifications including residue-specific chemical substitution of atoms, mutation of residues, and chemical cross-linking of proximal residues) 14–17 or proteins located in modified environments (i.e., proteins that are either immobilized on solid substrates or ionized to the gas phase or bound to ligands/substrates). 18–21 For instance, site-directed mutagenesis experiments belong to the former class of techniques while the surface plasmon resonance spectroscopy (SPRS), and electrospray ionization mass spectrometry (ESI-MS) fall into the latter class of techniques. 18–21 Through chemical modifications, a specific interaction of interest (for example, the hydrogen bonding between residues, or the salt bridges between the anionic carboxylate and the cationic ammonium of basic amino acids, or the disulphide bonds between cysteine residues or the artifical protein cross-links) can be either removed from or introduced into the protein. 22 The measured difference in the stabilities of the unmodified and modified proteins that differ only by the presence of a specific interaction serves as a metric of the contribution of the chosen specific interaction or residue to protein stability. The structure-based mutational studies compare the structures of a protein in the wildtype and mutant forms to calculate the degree of mutation-induced structural changes in the protein, which is then used as a metric to quantify the contribution of the specific mutated residue to protein stability. 23–26 Such mutational studies have demonstrated that mutation of some residues significantly alter the protein structure while some mutations do not perturb the structure much. The observed variation in the degree of mutation-induced structural changes among different residues indicates the heterogenous contributions of individual residues to the overall stability of the protein. Although the aforementioned experimental techniques that rely on modified proteins provide energetic insight into protein stability, a

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

major shortcoming of these techniques is that the difference in the Gibbs free energies of the unmodified (GW ) and modified (GM ) proteins, GW -GM , measured by these techniques only quantify the contribution of a specific interaction to the relative stability of the protein with respect to the modified state. However, none of these experimental techniques can measure the absolute contributions of individual residues to the overall stability of the unmodified protein. 15,23–26 To circumvent this problem, some experimental techniques use protein intrinsic probes to determine the interaction strengths of individual residues of proteins without chemically modifying them. For instance, the vibrational spectroscopic techniques use the vibrational frequencies of specific bonds (for example, C-D bonds) of proteins as probes to determine the environmental effects on the chosen bonds while the vibrational stark effect (VSE) experiments quantify the local electric fields around the carbonyl (C=O) and nitril (N≡C) bonds of proteins in the presence of a weak external electric field. 27–30 Similarly, the fluorescence spectroscopy exploits the sensitivity of the quantum yield, fluorescence intensity and stokes shift of the intrinsic fluorophores (Tryptophan, Phenylalanine, Tyrosine residues) in a protein to examine the strength of flurophore-environment interactions. 31,32 In H/D exchange experiments, the rate of H/D exchange of individual hydrogen atoms between the protein and D2 O is used as a sensitive probe to characterize the microenvironment around the exchangeable hydrogens of proteins. 33 However, it is non-trivial to establish a quantitative relationship between the interaction strengths of individual residues and the physical quantities measured from these experimental techniques. The solution nuclear magnetic resonance (NMR) spectroscopy has emerged as a powerful experimental technique to characterize the microscopic origin and diversity of site-specific fast sub-nanosecond dynamics of backbones and side chains of proteins and protein complexes of significant size in their native and less-populated excited states. 34–53 Using advanced isotopic enrichment and biosynthetic labeling strategies, the hydrogen and/or carbon atoms of side chain methyl groups of a protein are labeled with deuterium (2 H) and 13 C, respectively.

4

ACS Paragon Plus Environment

Page 4 of 51

Page 5 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The relaxations of these labelled methyl spin probes (13 CH3 , CD3 ,13 CHD2 ,13 CH2 D) are examined to quantify the fast motions of individual methyl-containing side chains of isotopically enriched proteins. Specifically, the Lipari-Szabo “model-free” formalism is used to extract the relaxation times (τ ) and amplitudes (O2 ) of motions of side chains from the measured Nuclear Overhauser effects (NOE), longitudinal (T1 ) and transverse (T2 ) relaxation times of 2 H- and/or

13

C-containing C-H bond vectors. Similarly, the motional parameters

obtained from the NMR relaxation data of the amide 15 N-H bond vectors quantify the backbone dynamics of proteins. 34–42,49–52,54–59 To probe fast side chain dynamics, the principal quantity measured from NMR experiments is the square of the generalized order parameter (O2 ) for each methyl-containing residue of a protein and it quantifies the degree of spatial restriction of a given methylbearing side chain. Given that the relaxations of methyl C-H bond vectors are primarily due to the reorientational dynamics of the methyl symmetry axis and the rotation about 2 2 2 the symmetry axis, O2 can be written as a product of Orot and Oaxis , where Orot is the 2 order parameter for methyl rotation about the symmetry axis and Oaxis is that for motion 2 =1/9 for a methyl group that rotates comof the symmetry axis. It can be shown that Orot 2 is assumed in pletely and possesses ideal tetrahedral geometry and this ideal value of Orot 2 can vary between 0 (maximum all NMR experimental studies. 35,43 The magnitude of Oaxis 2 flexibility) and 1 (maximum rigidity); the higher the value of Oaxis the more restricted is the 2 mobility of the side chain. By measuring Oaxis of individual residues and their distributions 2 )) for various proteins, NMR experiments have demonstrated that protein backbone (P(Oaxis

is dynamically more restricted than the side chains and that the fast side chain dynamics 2 is heterogenous in nature in all proteins. 39,42,43,46,60–62 The broader the distribution of Oaxis ,

the higher is the degree of dynamical heterogenity of side chain motions. 42,43,60,63,64 2 2 Using Oaxis as a dynamical proxy for protein conformational entropy, an Oaxis -based

entropy meter was developed for accurate determination of the conformational entropy of proteins in their native and excited states. 43,60,63–65 Using a combination of calorimetric and

5

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2 NMR experiments, the ligand-induced changes in Oaxis of individual residues of proteins

can be related to changes in protein conformational entropy and the total protein-ligand 2 -based microscopic connection between conformational binding entropy. 43,60,63,64 The Oaxis

dynamics and conformational entropy of proteins, once established, is necessary for a better understanding of the essential role of conformational entropy in protein thermodynamics, function, protein-ligand binding and molecular recognition. 43,60,63,64,66–68 Although it is gen2 erally believed that Oaxis of a side chain is highly sensitive to its local environment, a precise

understanding of how local environmental factors such as packing density, degree of solvent 2 and conformational entropy is still lacking. exposure, and depth of burial influence Oaxis 2 , NMR experiments on various proteins have shown no or weak correlations between Oaxis

Sconf and the aforementioned environmental factors. 43,61 For instance, it is believed that 2 than tightly packed side chains of proteins are dynamically more restricted with higher Oaxis

loosely packed side chains but some studies have demonstrated that fast dynamics of some residues are catalyzed by tight crystal packing of proteins. 69,70 In the present work, we demonstrate that strengths of noncovalent interactions of individual residues of proteins can be used as a suitable collective variable to probe the combined influence of various aforementioned environmental factors on fast side chain dynamics and conformational entropy of proteins. Drawing from the fact the side chain dynamics and thermodynamics are dictated by the underlying conformational energy surface and using the fundamental principles of statistical mechanics, the experimentally measurable metrics of 2 side chain flexibility (Oaxis ), conformational entropy (Sconf ), and strengths of noncovalent

interactions of individual residues of proteins are calculated accurately from their conformational free energy surfaces. The sensitivity of relative energies of side chain conformational states and the side chain conformational barriers to noncovalent interactions and their crit2 ical influences on Oaxis and Sconf are investigated by examining the correlation between 2 Oaxis , Sconf , activation barriers, relative rotamer stabilities, and strengths of noncovalent

interactions of side chains in proteins.

6

ACS Paragon Plus Environment

Page 6 of 51

Page 7 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Theory Estimation of Strength of Noncovalent Interactions of Residues in Protein We consider an interacting system consisting of a protein solvated in water and a few neutralizing ions. Let N be the total number of atoms in the system and their positions and momenta are given by the sets {r} = (r1 , r2 , ..., rN ) and {p} = (p1 , p2 , ..., pN ), respectively. The Hamiltonian of the system is given by H({r}, {p}) = T({p}) + U({r}), where U({r}) and T({p}) represent the potential and kinetic energies of the system, respectively. Given these microscopic details of the system, we seek to develop a theoretical framework to quantitatively examine how the strength of noncovalent interactions of individual residues with their surrounding microenvironment influences their fast side chain dynamics and conformational entropy. In particular, the primary aim of this theoretical excercise is to use statistical mechanical principles to establish a quantitative relationship between the fast side chain motional parameters obtained from site-specific NMR experiments and the strength of noncovalent interactions of methyl-containing residues of proteins. Once established, such a quantitative relationship would enable us to estimate the site-specific noncovalent interaction strengths of individual residues directly from the NMR order parameters. Although 2 Oaxis -based experimental and computational studies continue to enrich our understanding of

protein dynamics, thermodynamics and function, a simple method to estimate site-specific 2 noncovalent interactions from Oaxis has not been reported. To this end, the system of in-

terest is divided into two sub-systems (sub-system-I and sub-system-II) such that n atoms belong to the sub-system-I and the remaining N − n atoms constitute the sub-system-II (see Figure 1). Although the theoretical approach proposed here may be applicable to any sub-systems, in the present context of site-specific fast side chain dynamics of proteins, the sub-system-I corresponds to a methyl-containing residue and all the remaining atoms belong to the sub-system-II. Our aim is to determine the strength of noncovalent interaction of 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 51

sub-system-I with sub-system-II. Given that the conformational free energy landscape of a residue dictates its fast side chain dynamics, we first derive the microscopic relationship between the side chain conformational free energy and noncovalent interaction between the side chain of interest (subsystem-I) and its microenvironment (sub-system-II). Then we proceed to examine the relationship between the motional parameter and free energy profile that can be used to con2 struct free energy profile of a side chain from its Oaxis . 71,72 Let {rI } = (r1 , r2 , ...rn ) and

{pI } = (p1 , p2 , ..., pn ) denote the positions and momenta of atoms in the sub-system-I and {rII } = (rn+1 , rn+2 , ...rN ) and {pII } = (pn+1 , pn+2 , ..., pN ) denote the positions and momenta of the atoms in the sub-system-II. The total potential energy of the system U({rI }, {rII }) can be expressed in terms of the interactions within and between the subsystems as follows:

I−II I II U ({rI }, {rII }) = Uintra ({rI }) + Uintra ({rII }) + Uinter ({rI }, {rII })

(1)

II I ({rII }) correspond to the intra-system potential energies of ({rI }) and Uintra where Uintra

sub-system-I and sub-system-II, respectively. The label ”intra” denotes that these potential I−II energies include only atomic interactions within each subsystem. Uinter ({rI }, {rII }) is the

potential energy arising from atomic interactions between sub-system-I and sub-system-II. We define a conformational reaction coordinate φ ({rI }), which depends only on the coordinates of atoms of the sub-system-I, to describe various conformational states of the sub-system-I. The potential of mean force, F(φ), quantifying the free energy of the system along the chosen reaction coordinate φ is used as the central quantity to relate the fast side chain dynamics and the features of the conformational free energy surface. 71,72 Once F(φ) is determined accurately, the basic statistical mechanical principles are applied to derive the motional parameters, thermodynamic properties, the strength of noncolvalent interactions of sub-system-I directly from F(φ). To determine F(φ), the Cartesian coordinates and

8

ACS Paragon Plus Environment

Page 9 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

momenta of the atoms of sub-system-I are first transformed to the generalized coordinates (φ, Q) and generalized momenta (pφ , pQ ), where Q and pQ denote the other (n-1) generalized coordinates excluding φ and their corresponding generalized momenta, respectively, and pφ denotes the generalized momentum associated with φ. Given the Hamiltonian of the system expressed in terms of these generalized coordinates and generalized momenta H(φ, Q, {rII }, pφ , pQ , {pII }), the φ-dependent partial partition function, Z(φ), and potential of mean force, F(φ), can be written as follows: 1 Z(φ) = 3N h

Z Z

Z Z

Q



pQ

{rII }

Z

e−βH(φ,Q,{rII },pφ ,pQ ,{pII }) JdQdpφ dpQ d{rII }d{pII }

(2)

{pII }

F (φ) = −kB T ln Z(φ)

(3)

where J is the Jacobian associated with the transformation from (r,p) to (φ, Q, pφ , pQ ), h is the Planck’s constant, kB is the Boltzmann constant and T is the temperature. Using Eqs. 1-3, the derivative of F(φ) with respect to φ can be written as dF (φ) = dφ

 

=

∂U (φ, Q, {rII }) ∂ ln |J| − kB T ∂φ ∂φ

 (4) φ

I−II II I ({rII }) ∂Uinter (φ, Q, {rII }) ∂ ln |J| ∂Uintra (φ, Q) ∂Uintra + + − kB T ∂φ ∂φ ∂φ ∂φ

 (5) φ

where the angular brackets h i denote an ensemble average over all accessible configurations of the system and the subscript φ denotes that these averages are determined at a given value II of φ. Since Uintra {rII } is independent of φ,

II ∂Uintra {rII } ∂φ

= 0. Consequently, the above equation

can be written as dF (φ) = dφ



I ∂Uintra (φ, Q) ∂φ



 + φ

I−II ∂Uinter (φ, Q, {rII }) ∂φ



  ∂ ln |J| − kB T ∂φ φ φ

(6)

Let FP (φ) and FG (φ) denote the potentials of mean force calculated in the presence and 9

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 51

absence of the sub-system-II, respectively. In the absence of the sub-system-II, the subI−II system-I is said to be in the gas phase (G) with Uinter =0 and

I−II ∂Uinter ∂φ

= 0. The actual protein

environment (P) is taken into account in the presence of sub-system-II. Using Eq. 6, the derivatives of FP (φ) and FG (φ) with respect to φ can be written as 

I ∂Uintra ∂φ





I−II ∂Uinter ∂φ

  ∂ ln |J| + − kB T ∂φ φ,P φ,P φ,P  I    dFG (φ) ∂Uintra ∂ ln |J| = − kB T dφ ∂φ ∂φ φ,G φ,G dFP (φ) = dφ



(7) (8)

where the subscripts G and P indicate that the statistical averages are calculated in the gas phase and in the protein environment, respectively. D E D E ∂ ln|J| Since kB T ∂ ln|J| = k T , the difference in derivatives of the potentials of B ∂φ ∂φ φ,P

φ,G

mean force can be written as dFP (φ) dFG (φ) − = dφ dφ



I−II ∂Uinter ∂φ

"

 + φ,P

I ∂Uintra ∂φ



 − φ,P

I ∂Uintra ∂φ



# (9)

φ,G

By applying the common approximation used in vibrational stark effect experiments that the residue-specific intramolecular force constants (and intramolecular properties in D I E ∂Uintra 28 general) remain constant as the environment changes, we can assume that = ∂φ φ,P D I E ∂Uintra and show that ∂φ φ,G

dFP (φ) dFG (φ) − = dφ dφ Note that

D

I−II ∂Uinter ∂φ

E



I−II ∂Uinter ∂φ

 (10) φ,P

is a measure of the mean force on sub-system-I due to its noncovalent φ,P

interactions with sub-system-II, which upon integration with respect to φ yields the mean strength of noncovalent interaction between sub-system-I and sub-system-II. If ξ denotes the strength of noncovalent interaction between sub-system-I and sub-system-II, using Eq. 10,

10

ACS Paragon Plus Environment

Page 11 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

it can be shown that 1 ξ= 2π

Z

π

[FP (φ) − FG (φ)]dφ

(11)

−π

The above equation suggests that given FP (φ) and FG (φ) for a given residue, the corresponding ξ can be determined accurately.

Simulation Details Molecular dynamics (MD) simulations on set of eight proteins (eglin C (1EGL), ubiquitin (1UBQ), fibronectin (FNFn10), FYN tyrosine kinase SH3 domain (1SHF), adipocyte lipid-binding protein (1LIB), myoglobin (1J52), gammaF crystallin (1A45), HIV protease (3KF1)) 39,73–79 with different secondary structures and molecular weights were performed using NAMD2.8 80 with CHARMM27 all-atom 81 and TIP3P 82 water force-fields. The choice of six of the eight proteins (FNFn10, 1UBQ, 3KF1, 1EGL, 1SHF, 1LIB) studied here was based on the availability of previous experimental and molecular dynamics (MD) simulation studies. We have extended the study to two more proteins selected from the two extreme ends of the spectrum of the secondary structure composition: myoglobin (primarily an α-helix protein) and Gamma F crystallin (primarily a β-sheet protein) so that the selected set of eight proteins reasonably covers a moderate range of secondary structures and molecular weights while remaining within the limits of computational capabilities. 72 Each model system is solvated in a TIP3P water box of suitable dimensions: (54×47×49 ˚ A

3

3 3 3 (1EGL), 40×42×46 ˚ A (1UBQ), 45×59×61 ˚ A (FNfn10), 44×47×42 ˚ A (1SHF), 61×55×47 3 ˚ ˚3 (1J52), 58×54×80 A ˚3 (1A45) and 54×48×58 A ˚3 (3KF1). All the A (1LIB), 44×57×60 A

systems were subjected to energy minimization using the conjugate gradient method followed by 3 ns NPT equilibration MD runs at 1 atm pressure and at 300 K. Five independent 10ns production MD runs were performed in the NPT ensemble using a Langevin thermostat and barostat with a damping coefficient of 5 ps−1 . The equations of motion were integrated with 11

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

a time step of 1fs. The non-bonded pair interaction potential was truncated at 12 ˚ A and smoothed between 10 and 12 ˚ A using cubic switching function. Periodic boundary conditions were applied. Electrostatic interactions were computed using the particle mesh Ewald (PME) method with a real space cutoff of 13 ˚ A, and the reciprocal space interactions were computed on grids (using sixth-degree B-splines) of suitable dimensions: 58×50×54 (1EGL), 44×46×50 (1UBQ), 50×64×66 (FNfn10), 48×50×46 (1SHF), 64×58×50 (1LIB), 45×60×64 (1J52), 60×56×82 (1A45) and 56×50×60 (3KF1). The gas phase simulations were performed for each individual methyl-containing residue (ILE, VAL, LEU, MET, THR) without using periodic boundary conditions. Using the Collective Variable Module 83–86 and the adaptive biasing force (ABF) method implemented in NAMD2.8, the conformational free energy surfaces as a function of the methyl torsional angle, φ, as the reaction coordinate (Figure 1) of all methyl-containing side chains were determined in the gas phase and in the protein environment. The mean force as a function of the chosen reaction coordinate (φ) was accumulated in bins of width 1◦ from a total simulation time of 65ns for each system. The free energy profiles for all methylcontaining residues were recorded at regular intervals of 7.2 ns along the ABF trajectory. The convergence of free energy surfaces is monitored by estimating the root mean square deviation (RMSD) between the consecutive free energy profiles obtained. The conformational sampling is considered complete when RMSD is less than a set threshold value (0.05 kcal/mol) for all residues. 72 To asses the roboustness of the results, additional MD and ABF simulations of one of the proteins (ubiquitin) under study were performed using the AMBER ff99SBnmr forcefield. 87–89 The system size (total number of atoms and box dimensions), initial configuration, and simulation parameters were the same as in CHARMM simulations of ubiquitin. The AmberTools16 90 was used to generate the initial configuration and interaction parameters needed for the AMBER simulations and NAMD was used for MD and ABF runs.

12

ACS Paragon Plus Environment

Page 12 of 51

Page 13 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2 Estimation of Oaxis and Conformational Entropy

The angular distribution function, P(θ,φ) (θ and φ are defined in Figure 1), is calculated from F(φ) using the following relationship: e−βF (φ) e−

P (θ, φ) = R R π π −π

0

Here, θ0 is the equilibrium value of θ,

e−βF (φ) e− √σ β

β(θ−θ0 )2 2σ 2

β(θ−θ0 )2 2σ 2

(12)

sin θdθdφ

is the width of the Gaussian distribution char-

acterizing the fluctuations of θ around θ0 , β =

1 , kB T

kB is the Boltzmann’s constant and T

is the temperature. The values of θ0 and σ for each residue were determined by fitting the distribution of θ obtained from the MD trajectories with a Gaussian function. The Cartesian coordinates of the unit vector along the symmetry axis of a given methyl group can be expressed in terms of θ and φ as follows: x = sin(θ) cos(φ), y = sin(θ) sin(φ), and z = − cos(θ). 2 using Given P(θ,φ) and x, y, z in terms of θ and φ, it is straightforward to determine Oaxis

the following equation

2 Oaxis =

 1 3 < x2 >2 + < y 2 >2 + < z 2 >2 +2 < xy >2 +2 < xz >2 +2 < yz >2 − (13) 2 2

The side chain conformational entropy can be computed using the following equation: Z

π

Z

π

Sconf = −kB

P (θ, φ) ln P (θ, φ) sin θdθdφ 0

(14)

−π

P(θ,φ) defined in Eq. 12 can be used in the above equation to calculate Sconf for each residue.

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 51

Results and Discussion Side Chain Conformational Energy Profiles The calculated side chain conformational free energy surfaces of a few representative methylcontaining residues in the gas phase (FG (φ)) and in the protein environment (FP (φ)) are shown in Figure 2. Both FG (φ) and FP (φ) of all residues (except alanine residues) consist of three stable rotamer states: gauche+ (denoted as g+ ) with φ ∼ π3 , gauche− (g− ) with φ ∼

−π 3

and trans (t) with φ ∼ ±π. The relative energies of the rotamer states and the activation barriers separating them differ significantly among different residues both in the protein and in the gas phase. The MD-derived θ/φ scattered plots reveal clustering of data at φ values corresponding to minima on the conformational free energy surfaces and the absence of data at φ values corresponding to the free energy barriers indicating consistency between between MD and ABF results. The calculated difference FP (φ) − FG (φ) for a given residue is a measure of the degree of modulation of its side chain conformational energy surface by its noncovalent interactions with the surrounding microenvironment. For instance, a residue with FG (φ) = FP (φ) is likely to be insensitive to noncovalent interactions or to have gas phase-like microenvironment around them (Figure S1(a) - Supporting Information). For some residues, depending upon the strength of their noncovalent interactions, the relative positions and energies of rotamer states, and the activation barriers for rotameric transitions can be significantly different in the protein environment than those in the gas phase (Figure S1(b-f) - Supporting Information). A schematic of a few possible modifications of side chain conformational surfaces induced by noncovalent interactions is shown in Figure S1 (Supporting Information).

14

ACS Paragon Plus Environment

Page 15 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Dependance of Side Chain Conformational Barriers and Rotamer Stabilities on Noncovalent Interactions The limitation of the previously established quantitative relationship (Theory section) between ξ and the side chain conformational free energy surfaces is that both ξ and free energy surfaces are not readily amenable to experimental techniques. It is therefore neccessary to relate ξ to readily measurable experimental physical parameters. One such quantities is the activation barriers for side chain conformational transitions in proteins; NMR and neutron scattering experiments provide accurate measures of the side chain activation barriers on per-residue basis. 91–93 Here, we seek to understand the relationship between ξ and the activation barriers for side chain conformational transitions. The other measurable quantities 2 and Sconf and we have established a parametric relationship that can be related to ξ are Oaxis 2 and between ξ and Sconf in the following sections. In this section, we between ξ and Oaxis

derive an exact relationship between ξ, the activation barriers and the relative energies of side chain rotamer states by considering the following three different cases of side chains. Case-I: In this case, we consider side chains whose rotamer states are isoenergetic and are separated by barriers of same height. Although it is generally expected that both the relative energies and barriers of rotamer states and their sensitivity to ξ are likely to be different for different rotamer states, the simplification of isoenergetic rotamer states with the same barrier height used here enables us to neglect the effect of noncovalent interactions on relative rotamer energies and allows us to investigate only the dependance of barriers on ξ. In Case-II and Case-III, the generalised results involving rotamer states with different relative energies and barriers are derived using two different methods. Case-II uses the results of Case-I and considers a piecewise integration of the free energy profiles of individual residues to estimate ξ while Case-III uses the truncated Fourier expansion of free energy profiles of individual residues to determine ξ. 94 In what follows, Case-II will be referred to as the piecewise integration approach and Case-III will be referred to as the Fourier expansion method. 15

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 51

Case-I: Side Chains with Isoenergetic Rotamer States Separated by Barriers of Same Height The side chain conformational energy surfaces with three isoenergetic rotamer states (g+ at φ ≈ π3 , g− at φ ≈

−π , 3

and t at φ ≈ ±π ) that are separated by barriers of same height are

represented by a three-fold symmetric periodic function of the following form

F (φ) =

∆E (1 − cos(3φ − δ)) 2

(15)

where ∆E is the activation barrier, δ is the phase factor, and the range of φ is −π ≤ φ ≤ π. This periodic conformational energy surface (with a periodicity of 2π) consists of three isoenergetic minima at φ =

δ 3



2π , 3

φ =

δ , 3

and φ =

δ 3

stable rotamer states (for δ = π; g+ at φ ≈ π3 , g− at φ ≈

2π 3

corresponding to the three

−π , 3

and t at φ ≈ ±π) and they

+

are separated by activation barriers of height ∆E. Let us assume that the conformational energy surface of a residue remains three-fold symmetric both in the gas phase (FG (φ)) and in the protein environment (FP (φ)). Let ∆EG and ∆EP denote the activation barriers in the gas phase and in the protein environment, respectively and the corresponding phase factors are denoted by δG and δP , respectively. For this simple case of three-fold symmetric surface, it can be shown that Z π 1 [FP (φ) − FG (φ)] dφ ξ = 2π −π  Z  1 π ∆EP ∆EG = (1 − cos(3φ − δP )) − (1 − cos(3φ − δG )) dφ 2π −π 2 2 ∆EP − ∆EG = 2 where we have used the fact that

Rπ −π

(16)

(17)

cos(3φ − δ)dφ = 0

This simple relationship between ξ and the activation barriers provides interesting insights into how the activation barriers are influenced by ξ. Firstly, it is known that ∆EP and ∆EG are finite for all residues of all proteins. Consequently, based on the Eq. 17, it 16

ACS Paragon Plus Environment

Page 17 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

can be concluded that ξ must also be finite for all residues of all proteins. Moreover, Eq. 17 also dictates the allowed range of ξ and this range of ξ for a given residue type is critically controlled by the lower and upper bounds of the side chain conformational barriers of residues of that residue type. Let ∆EP,min and ∆EP,max be the minimum and maximum activation barriers for a given residue type in the protein environment, respectively, and the corresponding values in the gas phase are denoted by ∆EG,min and ∆EG,max , respectively. From Eq. 17, it is evident that the parameters ∆EP,min , ∆EP,max , ∆EG,min , and ∆EG,max determine the range of ξ, which is ∆EP,min − ∆EG,max ≤ ξ ≤ ∆EP,max − ∆EG,min . Secondly, it is evident that ξ can be negative or positive or zero depending upon the values of ∆EP and ∆EG . ξ of a residue would be negative when ∆EP < ∆EG (i.e., the activation barrier is less in the protein environment than in the gas phase). This implies that the fast side chain conformational dynamics of methyl-containing residues with negative ξ is catalyzed by noncovalent interactions resulting in a reduction of the activation barrier in the protein environment with respect to that in the gas phase. The nonbonded environments of residues with positive ξ (∆EP > ∆EG ) restrict the side chain dynamics by enhancing their activation barriers in the protein environment. The spatially heterogenous fast motions of side chains observed in NMR studies on various proteins suggest that noncovalent environments of some residues enhance their conformational barriers while some other residues may experience barrier reduction in the protein environment to maintain a balance between structural stability and functionally-relevant flexibility.

Case-II: Side Chains with Rotamers of Different Stabilities and Barriers - A Piecewise Integration Approach The three-fold symmetric free energy function used in Case-I may not truly represent the actual conformational energy surfaces of various side chains of proteins. Given the heterogenity of microenvironments and diverse chemical nature of side chains, the actual side chain conformational energy surfaces are expected to be asymmetric with rotamers of dif-

17

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 51

ferent relative energies separated by barriers of differing heights. A schematic of a typical conformational energy surface of a side chain is shown in Figure 3 for illustrative purposes. There are three rotamer states at φ ≈

π 3

(g+ rotamer state), φ ≈

−π 3

(g− rotamer state), and

φ ≈ ±π (t rotamer state) and each rotamer state is surrounded by a pair of barriers giving rise to a total of six barriers per side chain. Let E1rot , E2rot , and E3rot denote the energies of the t, g − , and g + rotamer states, respectively and ∆E1 , ∆E2 , ∆E3 , ∆E4 , ∆E5 , ∆E6 denote the six barriers of a side chain (see Figure 3). It was demonstrated earlier that the distribtion of rotamer barriers of side chains follow a trimodal distribution for many proteins. 72 Given the variability of ∆En (n =1, 2,...., 6), E1rot , E2rot , and E3rot among different side chains, our primary aim is to derive an exact mathematical relationship between ξ and these parameters that would accurately describe the dependance of these parameters on ξ. To derive this relationship, the conformational free energy surface of a side chain is divided into six equal segments of interval

π . 3

The segments are chosen such that each

segment begins (or ends) at an energy minimum and ends (or begins) at the adjacent energy maximum associated with that minimum (see Figure 3). Given the range of φ (−π ≤ φ ≤ π) and the positions of the energy minima (φmin,1 = ±π, φmin,2 =

−π , 3

and φmin,3 = π3 ) on the

energy surface, the range of each segment can be shown to be (n − 4) π3 ≤ φ ≤ (n − 3) π3 (for n =1, 2,..., 6), where n is the index of a given segment. For instance, for the energy surface shown in Figure 3, n =1 corresponds to a segment between φ = −π and φ =

−2π . 3

Each

segment is then fitted using a shifted three-fold symmetric function of the form Fn (φ) = n Enmin + ∆E [1 − cos (3φ − π)], where Enmin and ∆En correspond to the energy of the rotamer 2

state (i.e., the energy minimum) and the associated barrier corresponding to the nth segment of the energy surface (see Figure 3). For instance, the fitting function for the first segment (i.e., n = 1) with −π ≤ φ ≤

−2π 3

would be F1 (φ) = E1min +

∆E1 2

[1 − cos (3φ − π)]. A few

representative free energy profiles fitted with piecewise three-fold symmetric functions are shown in Figure 4 and in Figure S2 (Supporting Information). The piecewise fitted free energy functions so obtained can be used to derive an exact relationship between ξ, ∆En ,

18

ACS Paragon Plus Environment

Page 19 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

E1rot , E2rot , and E3rot . The values of ξ for the individual segments of the free energy profile can be calculated as follows

ξn

Z (n−3) π 3 1 = [Fn,P (φ) − Fn,G (φ)] dφ 2π (n−4) π3    Z (n−3) π  3 ∆En,G ∆En,P min min [1 − cos (3φ − π)] − En,G + [1 − cos (3φ − π)] dφ = En,P + 2 2 (n−4) π3    1 1 min min = En,P − En,G + (∆En,P − ∆En,G ) (18) 6 2

The net ξ of a residue determined by adding ξn of all six segments (i.e., ξ =

P6

n=1 ξn )

can

be written as

ξ =

where h∆EiP =

1 6

P6

n=1

 1 hEPmin i − hEGmin i + (h∆EP i − h∆EG i) 2

∆En,P and h∆EiG =

1 6

P6

n=1

∆En,G are the mean activation barriers

in the protein environment and in the gas phase, respectively. hEPmin i = and hEGmin i =

rot +∆E rot +∆E rot ∆E1,G 2,G 3,G 3

(19)

rot +∆E rot +∆E rot ∆E1,P 2,P 3,P 3

are the mean rotamer energies of a side chain in the protein

environment and in the gas phase, respectively. Eq. 19 is a generalised form of Eq. 17 (see Case-I) and it accurately describes the dependance of relative stabilities of rotamer states and activation barriers on the strength of noncovalent interactions of side chains in proteins. If the rotamer states are isoenergetic both in the gas phase and protein environment with hEPmin i = hEGmin i, we arrive at the result obtained in Case-I. In Eq. 19, hEPmin i − hEGmin i is a measure of the difference between the mean energies of rotamer states in protein environment and in the gas phase. The observation that ξ is directly proportional to hEPmin i − hEGmin i indicates that the relative stabilities of rotamer states can be sensitively tuned by varying the strength of noncovalent interactions of side chains. Eq. 19 also provides a systematic way of determining ξ from experimentally measurable activation barriers of residues in the gas phase and protein environment. To validate the accuracy of Eq. 19, we calculated the mean relative stabilites of conformational 19

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 51

states and activation barriers of the side chains directly from the ABF-derived conformational free energy profiles and determined ξ using Eq. 19. Figure 5 and Figure S3 (Supporting Information) shows the comparison of ξ obtained using the piecewise integration approach (denoted as ξS ) with ξABF obtained directly from free energy surface (using Eq. 11). It is evident from Figure 5 and Figure S3 (Supporting Information) that the values of ξ obtained using Eq. 11 and using Eq. 19 are consistent with each other.

Case-III: Fourier Expansion of Free Energy Profiles The dependance of ξ on activation barriers can also be investigated using the truncated Fourier expansion method. 72,94 In this approach, the conformational energy surface of a side chain is expressed as a combination of one-fold, two-fold, and three-fold symmetric energy functions as follows

F (φ) = ∆E0 +

3 X ∆En n=1

2

[1 − cos n (φ − δn )]

(20)

where δn is the phase angle and ∆E0 and ∆En are the coefficients of the cosine expansion. The absolute values of the coefficients ∆E1 , ∆E2 , and ∆E3 are the measures of the barriers in the one-fold, two-fold, and three-fold energy functions, respectively. The three-fold term is due to steric interactions arising from bond-bond repulsion, the one-fold term is mainly due to dipolar and steric interactions, while the two-fold term is due to charge delocalization between lone pair orbitals and the rotating bonds. 72,94 A few representative side chain free energy profiles fitted using the truncated Fourier expansion are shown in Figure 6. Given the Fourier coefficients and parameters involved in the Fourier expansion of the conformational energy surface of a methyl-containing residue in the protein environment (∆E0,P , ∆En,P ,

20

ACS Paragon Plus Environment

Page 21 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

δn,P ) and in the gas phase (∆E0,G , ∆En,G , δn,G ), ξ can be written as 1 ξ = 2π

Z

1 − 2π

Z

π

∆E0,P + −π π

−π

= (∆E0,P = (∆E0,P f where h∆Eief = P

1 3

P3

n=1

"

"

3 X ∆En,P

2

n=1 3 X

# [1 − cos n (φ − δn,P )] dφ

# ∆En,G ∆E0,G + [1 − cos n (φ − δn,G )] dφ 2 n=1  3  X ∆En,P − ∆En,G − ∆E0,G ) + 2 n=1 h i 3 ef f ef f − ∆E0,G ) + h∆EiP − h∆EiG 2

f ∆En,P and h∆Eief G =

1 3

P3

n=1

(21)

∆En,G are the mean effective bar-

riers obtained by averaging the barriers in the one-fold, two-fold, and three-fold symmetric energy functions in the protein environment and in the gas phase, respectively. It is to be noted that h∆Eief f,P and h∆Eief f,G are different from the mean barriers h∆EiP and h∆EiG , which were obtained by averaging over the six barriers on the side chain conformational energy surface, reported in the piecewise integration approach (Eq. 19). The derived relationship between ξ, h∆Eief f,P and h∆Eief f,G is different from Eq. 19 in that Eq. 21 allows us to probe the relative contributions of specific interactions including dipolar and steric interactions, charge delocalization between lone pair orbitals and the rotating bonds, and steric interactions arising from bond-bond repulsion to the conformational barriers and ξ. A comparison of ξ obtained using the Fourier expansion approach (denoted as ξF ) with ξABF obtained directly from free energy surface (Eq. 11) shown in Figure 7 and Figure S4 (Supporting Information) reveals that the values of ξ obtained using Eq. 11 and using Eq. 21 are consistent with each other.

21

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Residue-specific Variation of Noncovalent Interaction Strengths in Proteins Given FG (φ) and FP (φ), ξ (defined in the Eq. 11) can be calculated for all methyl-containing residues of all proteins. As indicated in the Theory section, ξ provides a site-specific measure of the mean strength of noncovalent energy of the individual residues of proteins. Using the set of ξ so obtained for different residues, the ranges of ξ for different residue types can be examined. The results reveal that ξ varies between -4.3 kcal/mol and 3.73 kcal/mol for all proteins. The calculated minimum (ξmin ), maximum (ξmax ), and mean values (hξi) of ξ for different residue types are provided in Table 1. The results indicate that ξ is positive for all δ-methyls (LEUδ1 , LEUδ2 and ILEδ ) and it is negative for γ-methyls (VALγ1 , VALγ2 , THRγ2 ) and methionine residues. The least value of hξi observed for THR residues (hξi ≈ −3 kcal/mol) indicates that THR residues are likely to be found in most stable environment relative to other methyl-containing residues. Since THR is the only hydrophilic methyl-containing residue, it is mostly located on protein surfaces and exposed to solvent. The excess stability of THR residues in proteins can be attributed to their ability to form hydrogen bonds with solvent and other residues of proteins. hξi of γ-methyls are relatively lower than those of δ-methyls. The order of hξi of γ− methyls is THRγ2 >VALγ1 > ILEγ2 > VALγ2 . A broad range of ξ observed for ILEγ2 and VALγ1 residues suggests that these residues are likely to be located in diverse microenvironments in proteins. Figure 8 shows the magnitudes of ∆ξ (∆ξ = ξ − ξmin ) for various methyl-containing side chains (methyl carbons are color coded based on the ∆ξ values of the corresponding residues) of proteins. It is observed that the methyl-bearing residues of myoglobin with ∆ξ ∼0 are located near internal cavities known as xenon cavities 95,96 with large free volumes indicating that side 2 chains near internal cavities in proteins are likely to be highly flexible with lower Oaxis 2 values (∆ξ versus Oaxis correlation is detailed in the next section). This observation is in

line with an eariler study that demonstrated that internal cavities act as activation centers of fast side chain dynamics in proteins. 97 22

ACS Paragon Plus Environment

Page 22 of 51

Page 23 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Side Chain Flexibility Versus Noncovalent Interaction Strength The correlation between residue-specific noncovalent interactions and side chain dynamics 2 for individual residues of proteins was examined by comparing ∆ξ (∆ξ = ξ − ξmin ) and Oaxis

(Figure 9). Although the range of ξ is different for different residue types (see Table 1), ∆ξ varies from zero to ∼ 4.8 kcal/mol for all residues of all proteins. Figure 9 and Figure S5 2 (Supporting Information) show that Oaxis increases with ∆ξ in a sigmoidal-like fashion and 2 2 reaches a saturation value of Oaxis ∼ 0.75 at ∆ξ ∼ 3 kcal/mol; Oaxis is insensitive to ∆ξ for 2 versus ∆ξ data indicates that the side chains with ∆ξ > 3 kcal/mol. The observed Oaxis

∆ξ greater than 3 kcal/mol are expected to be conformationally restricted within a rotamer 2 2 2 is ∼ 0.75). Although the theoretical upper limit of Oaxis value (Oaxis state with higher Oaxis 2 ∼ 0.75 can be attributed to libration of methyl one, the calculated saturation value of Oaxis 2 versus ∆ξ correlation shown in symmetry axis in the most-probable rotamer state. The Oaxis

Figure 9 indicates that regardless of the topology and molecular weight of a protein, for each residue type, the highly flexible residue of that type posseses the lowest ∆ξ and rigid side 2 chains posses higher ∆ξ. This unique universal correlation between ∆ξ and Oaxis observed

for all proteins reveals that there exists a critical range of ∆ξ for proteins within which the side chain dynamics is highly sensitive to noncovalent interactions and that within this critical range side chain flexibility decreases gradually with a gradual increase in ∆ξ and that the rate of change of flexibility with respect to ∆ξ is same for all residue types of all proteins. 2 The observed ∆ξ versus Oaxis correlation was fitted by using the following sigmoidal

function

2 Oaxis

2 2 Oaxis,max − Oaxis,min 2 = + Oaxis,min 1 + a10 e−a1 ∆ξ

(22)

2 2 and the best-fit parameters are a0 =0.1165, a1 =1.67, Oaxis,max =0.75, Oaxis,min =0.0005. 2 Here, the parameter a1 quantifies the rate of change of Oaxis with ∆ξ and the best-fit value

23

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 51

of a1 suggests that it is approximately equal to the inverse of the thermal energy at a given temperature, T i.e., a1 ≈

1 kB T

2 2 and Oaxis,min where kB is Boltzmann’s constant. Oaxis,max

2 . correspond to the calculated maximum and minimum values of Oaxis

Conformational Entropy Versus Noncovalent Interaction Strength 2 Many experimental and theoretical studies have examined the correlation between Oaxis and

the side chain conformational entropy (Sconf ). 42,60,63,64,72,98 However, the dependence of Sconf on noncovalent interaction strength of residues in proteins has not been studied much. To probe the sensitivity of the conformational entropy associated with sub-nanosecond timescale fast side chain dynamics to noncovalent interactions, the dependence of the calculated Sconf on ∆ξ was examined (Figure 10 and Figure S6 (Supporting Information)). The observed ∆ξ versus Sconf correlation reveals that Sconf gradually decreases with increase in ∆ξ and reaches a saturation at ∆ξ ∼ 3 kcal/mol; For ∆ξ > 3 kcal/mol, Sconf is insensitive to changes in ∆ξ. For all proteins investigated here, Sconf varies between 3.73 cal/mol/K and 9.13 cal/mol/K with MET residues having highest Sconf (range 7.76-9.13 cal/mol/K). The minimum, maximum, mean values and the range of Sconf for different residue types are provided in Table 2. The calculated ranges for different residue types are consistent with results obtained from NMR experiments and earlier MD studies on various proteins. 63,64,68,98 The residues with lower ∆ξ have higher Sconf than those with higher ∆ξ. The comparison of 2 the observed Oaxis versus ∆ξ (Figure 9) and Sconf versus ∆ξ correlations indicate that the

side chains with higher conformational flexibility posses higher conformation entropy and those with restricted dynamics have least conformational entropy. The calculated ∆ξ versus Sconf data can be fitted using the following function

Sconf =

1 + b0 e−b1 min min max (Sconf − Sconf ) + Sconf (b ∆ξ−b ) 2 1 1 + b0 e

(23)

max min and the best-fit parameters are b0 =1.02, b1 =1.26, b2 =1.21, Sconf =6.75, Sconf =3.95. Here,

24

ACS Paragon Plus Environment

Page 25 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

max min the parameter b2 quantifies the rate of change of Sconf with ∆ξ. Sconf and Sconf correspond

to the calculated maximum and minimum values of Sconf

Conclusions It is widely recognized that side chain conformational fluctuations are essential for protein function. 43,49,51,68 Solution NMR spectroscopy has emerged as a powerful technique to probe the diversity, origin, and responses to perturbations (including ligand-binding, pressure, and temperature) of sub-nanosecond fast side chain conformational dynamics of proteins. 62,99,100 By leveraging the advances in isotopic-labelling strategies and the advent of new pulse sequences, site-specific NMR experiments permit more precise measurement of 2 ), which are measures of the degree of side chain fast side chain motional parameters (Oaxis 2 ≤ 1 (highly restricted), of individual residues of flexibility with 0 (highly flexible) ≤ Oaxis 2 proteins. 34–42,49–52,56–59 A broad distribution of NMR-measured Oaxis for many proteins in-

dicated that the heterogenous nature of fast side chain motions is necessary for proteins to maintain a balance between flexibility and stability. Although it is generally believed that 2 ) of a side chain is highly sensitive to the microenvironment surrounding flexibility (i.e., Oaxis

it, NMR experiments have shown a weak or no correlation between side chain flexibility and various microenvironmental factors including the local packing density, solvent exposure, and the degree of burial of side chains. 43,45–53,61,61 In the present study, using strength of noncovalent interaction (ξ) of each residue with its surroundings as a collective variable, the combined effect of various aforementioned mi2 croenvironmental factors on side chain flexibility (i.e., on Oaxis ) and conformational entropy

(i.e., on Sconf ) was investigated. Using the fundamental principles of statistical mechanics, an exact relationship between the conformational free energy surfaces and strengths of noncovalent interaction of protein side chains was derived. Given this mathematical relationship and the conformational energy surfaces of individual residues determined accurately using

25

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the adaptive biasing force (ABF) enhanced sampling method, ξ of individual residues were 2 2 as a , Sconf on ξ was examined. While, using Oaxis calculated and the dependance of Oaxis 2 and Sconf has been proxy for conformational entropy (Sconf ), the correlation between Oaxis

studied extensively for many proteins, the precise correlation between side chain flexibility, conformational entropy and noncovalent interactions of side chains has not been explored much. Our results reveal that the relative stabilities of different side chain conformational states, and the activation barriers separating them are highly sensitive to ξ. An exact mathematical relationship between ξ, activation barriers, and the relative stabilities of different conformational states of side chains derived using three different methods (Fourier expansion method, a piecewise integration approach, and direct integration of conformational free energy surfaces) indicate that the allowed range of ξ for a given residue type is critically determined by the range of the mean activation barriers of residues of that particular type. Given that the mean activation barriers of side chains of proteins are finite, our results unambiguously demonstrate that the allowed range of ξ for all residues of all proteins must also be finite and that the range of ξ is different for different residue types (VALγ1 and ILEγ2 have highest ξ range followed by ILEδ and LEUδ2 whereas MET has least range). The calculated values of ξ for eight proteins investigated here range between -4.3 kcal/mol and 3.73 kcal/mol. 2 The derived correlations between Oaxis , Sconf and ξ using CHARMM and AMBER force

field parameters reveal that side chain flexibility and conformational entropy decrease with increasing ξ in a sigmoidal-like fashion for all residue types. The observed sigmoidal-like 2 Oaxis versus ξ and Sconf versus ξ correlations were found to be common among proteins 2 suggesting a possibility of a universal Oaxis -based (or Sconf -based) noncovalent interaction 2 strength scale (analogous to the Oaxis -based entropy meter 60,63,64 ) to determine ξ directly 2 from NMR-measured Oaxis . These universal parametric relationships also suggest that side 2 chain flexibility of any residue can be reversibly varied from highly flexible (with Oaxis ∼ 0) 2 to highly restricted (with Oaxis ∼ 1) by increasing its noncovalent interaction strength by ∼3

26

ACS Paragon Plus Environment

Page 26 of 51

Page 27 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

kcal/mol. The results also indicate that for each residue type there exists a critical range of ξ determined primarily by the mean side chain conformational barriers, beyond which side chain flexibility and conformational entropy are insensitive to noncovalent interactions. The present results not only provided insights into the influence of noncovalent interactions on side chain dynamics and thermodynamics of proteins but also illustrated the complementary nature of advanced computational techniques in providing microscopic insights into fast side chain motions of proteins. It is anticipated that the established mathematical correlations 2 , Sconf , activation barriers, and ξ and their conbetween experimentally measurable Oaxis

nections with the conformational free energy landscape may be useful for rational design of potential drugs for protein targets. The microscopic connection between side chain conformational dynamics, conformational entropy, and noncovalent interactions of side chains in proteins established here also opens up new avenues to probe perturbation-induced (for example, ligand-binding, temperature, pressure) changes in fast side chain dynamics and thermodynamics by comparing conformational energy surfaces of proteins in the native and perturbed states. The present work is timely given that the noise-filtered electron densities around side chain dihedral angles obtained from X-ray diffraction and DFT-based methods are beginning to provide access to side chain conformational energy surfaces and the nature of noncovalent interactions of proteins. 101–104

Acknowledgement R.R.T. acknowledges Department of Science and Technology (DST) for the INSPIRE fellowship.

Supporting Information Available Additional Figures (Figures S1-S6) are provided.

27

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

References (1) Kannan, N.; Vishveshwara, S. Aromatic Clusters: A Determinant of Thermal Stability of Thermophilic Proteins. Protein Eng. 2000, 13, 753–761. (2) Griffiths-Jones, S. R.; Searle, M. S. Structure, Folding, and Energetics of Cooperative Interactions between the β-Strands of a de Novo Designed Three-Stranded Antiparallel β-Sheet Peptide. J. Am. Chem. Soc. 2000, 122, 8350–8356. (3) Kim, K. S.; Tarakeshwar, P.; Lee, J. Y. Molecular Clusters of π-Systems: Theoretical Studies of Structures, Spectra, and Origin of Interaction Energies. Chem. Rev. 2000, 100, 4145–4186. ˇ (4) Sponer, J.; Leszczynski, J.; Hobza, P. Electronic Properties, Hydrogen Bonding, Stacking, and Cation Binding of DNA and RNA Bases. Biopolymers 2001, 61, 3–31. (5) Lamoureux, J. S.; Maynes, J. T.; Glover, J. M. Recognition of 5’-YpG-3’ Sequences by Coupled Stacking/Hydrogen Bonding Interactions with Amino Acid Residues. J. Mol. Biol. 2004, 335, 399–408. (6) Robertson, A. D.; Murphy, K. P. Protein Structure and the Energetics of Protein Stability. Chem. Rev. 1997, 97, 1251–1268. (7) Gilson, M. K.; Honig, B. Calculation of the Total Electrostatic Energy of a Macromolecular System: Solvation Energies, Binding Energies, and Conformational Analysis. Proteins 1988, 4, 7–18. (8) Warshel, A.; Russell, S. T. Calculations of Electrostatic Interactions in Biological Systems and in Solutions. Q. Rev. Biophys. 1984, 17, 283–422. (9) Frieden, E. Non-covalent Interactions: Key to Biological Flexibility and Specificity. J. Chem. Educ. 1975, 52, 754.

28

ACS Paragon Plus Environment

Page 28 of 51

Page 29 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(10) Angelina, E. L.; Andujar, S. A.; Tosso, R. D.; Enriz, R. D.; Peruchena, N. M. Noncovalent Interactions in Receptor–ligand Complexes. A Study Based on the Electron Charge Density. J. Phys. Org. Chem. 2014, 27, 128–134. (11) Rose, G. D.; Wolfenden, R. Hydrogen Bonding, Hydrophobicity, Packing, and Protein Folding. Annu. Rev. Biophy. Biom. 1993, 22, 381–415. (12) Matthews, C. R. Pathways of Protein Folding. Annu. Rev. Biochem. 1993, 62, 653– 683. (13) Burley, S.; Petsko, G. Weakly Polar Interactions in Proteins. Adv. Protein Chem. 1988, 39, 125–189. (14) Goldenberg, D. P.; Creighton, T. E. Circular and Circularly Permuted Forms of Bovine Pancreatic Trypsin Inhibitor. J. Mol. Biol. 1983, 165, 407–413. (15) Fersht, A. R. Relationships Between Apparent Binding Energies Measured in Sitedirected Mutagenesis Experiments and Energetics of Binding and Catalysis. Biochemistry 1988, 27, 1577–1580. (16) Fersht, A. R.; Shi, J.; Knill, J.; Jones, D.; Lowe, A.; Wilkinson, D.; Blow, P. B.; Carter, P.; Waye, M.; Winter, G. Hydrogen Bonding and Biological Specificity Analysed by Protein Engineering. Nature 1985, 314, 235–238. (17) Varadarajan, R.; Lambright, D. G.; Boxer, S. G. Electrostatic Interactions in Wildtype and Mutant Recombinant Human Myoglobins. Biochemistry 1989, 28, 3771– 3781. (18) Patching, S. G. Surface Plasmon Resonance Spectroscopy for Characterisation of Membrane Protein-Ligand Interactions and its Potential for Drug Discovery. BBABiomembranes 2014, 1838, 43–55.

29

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(19) Loo, J. A. Studying Noncovalent Protein Complexes by Electrospray Ionization Mass Spectrometry. Mass Spectrom. Rev. 1997, 16, 1–23. (20) Daniel, J. M.; Friess, S. D.; Rajagopalan, S.; Wendt, S.; Zenobi, R. Quantitative Determination of Noncovalent Binding Interactions Using Soft Ionization Mass Spectrometry. Int. J. Mass Spectrom. 2002, 216, 1–27. (21) Veenstra, T. D. Electrospray Ionization Mass Spectrometry in the Study of Biomolecular Non-covalent Interactions. Biophys. Chem. 1999, 79, 63–79. (22) Goldenberg, D. P. Dissecting the Roles of Individual Interactions in Protein Stability: Lessons from a Circularized Protein. J. Cell. Biochem. 1985, 29, 321–335. (23) Bloom, J. D.; Labthavikul, S. T.; Otey, C. R.; Arnold, F. H. Protein Stability Promotes Evolvability. P. Natl. Acad. Sci. USA 2006, 103, 5869–5874. (24) Masso, M.; Vaisman, I. I. Accurate Prediction of Stability Changes in Protein Mutants by Combining Machine Learning with Structure Based Computational Mutagenesis. Bioinformatics 2008, 24, 2002–2009. (25) Makhatadze, G. I.; Privalov, P. L. Energetics of Protein Structure. Adv. Protein Chem. 1995, 47, 307–425. (26) Matthews, B. W. Structural and Genetic Analysis of Protein Stability. Annu. Rev. Biochem. 1993, 62, 139–160. (27) Adhikary, R.; Zimmermann, J.; Liu, J.; Dawson, P. E.; Romesberg, F. E. Experimental Characterization of Electrostatic and Conformational Heterogeneity in an SH3 Domain. J. Phys. Chem. B 2013, 117, 13082–13089. (28) Fried, S. D.; Boxer, S. G. Measuring Electric Fields and Noncovalent Interactions Using the Vibrational Stark Effect. Accounts Chem. Res. 2015, 48, 998–1006.

30

ACS Paragon Plus Environment

Page 30 of 51

Page 31 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(29) Park, E. S.; Andrews, S. S.; Hu, R. B.; Boxer, S. G. Vibrational Stark Spectroscopy in Proteins: A Probe and Calibration for Electrostatic Fields. J. Phys. Chem. B 1999, 103, 9813–9817. (30) Chattopadhyay, A.; Boxer, S. G. Vibrational Stark Effect Spectroscopy. J. Am. Chem. Soc. 1995, 117, 1449–1450. (31) Weiss, S. Fluorescence Spectroscopy of Single Biomolecules. Science 1999, 283, 1676– 1683. (32) Ward, L. D. Measurement of Ligand Binding to Proteins by Fluorescence Spectroscopy. Method enzymol. 1985, 117, 400–414. (33) Garcia, R. A.; Pantazatos, D.; Villarreal, F. J. Hydrogen/Deuterium Exchange Mass Spectrometry for Investigating Protein-Ligand Interactions. Assay Drug Dev. Techn. 2004, 2, 81–91. (34) Muhandiram, D.; Yamazaki, T.; Sykes, B. D.; Kay, L. E. Measurement of 2H T1 and T1. rho. Relaxation Times in Uniformly 13C-Labeled and Fractionally 2H-Labeled Proteins in Solution. J. Am. Chem. Soc. 1995, 117, 11536–11544. (35) Chatfield, D. C.; Augsten, A.; D’Cunha, C. Correlation Times and Adiabatic Barriers for Methyl Rotation in SNase. J. Biomol. NMR 2004, 29, 377–385. (36) Chatfield, D. C.; Wong, S. E. Methyl Motional Parameters in Crystalline lAlanine:Molecular Dynamics Simulation and NMR. J. Phys. Chem. B 2000, 104, 11342–11348. (37) Chatfield, D. C.; Szabo, A.; Brooks, B. R. Molecular Dynamics of Staphylococcal Nuclease: Comparison of Simulation with 15 N and 13 C NMR Relaxation Data. J. Am. Chem. Soc. 1998, 120, 5301–5311.

31

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(38) Karplus, M.; Petsko, G. A. Molecular Dynamics Simulations in Biology. Nature 1990, 347, 631–639. (39) Best, R. B.; Clarke, J.; Karplus, M. What Contributions to Protein Side-chain Dynamics are Probed by NMR Experiments? A Molecular Dynamics Simulation Analysis. J. Mol. Biol. 2005, 349, 185–203. (40) Best, R. B.; Clarke, J.; Karplus, M. The Origin of Protein Sidechain Order Parameter Distributions. J. Am. Chem. Soc. 2004, 126, 7734–7735. (41) Case, D. A. Molecular Dynamics and NMR Spin Relaxation in Proteins. Acc. Chem. Res. 2002, 35, 325–331. (42) Lee, A. L.; Wand, A. J. Microscopic Origins of Entropy, Heat Capacity and the Glass Transition in Proteins. Nature 2001, 411, 501–504. (43) Igumenova, T. I.; Frederick, K. K.; Wand, A. J. Characterization of the Fast Dynamics of Protein Amino Acid Side Chains Using NMR Relaxation in Solution. Chem. Rev. 2006, 106, 1672–1699. (44) Vallurupalli, P.; Hansen, D. F.; Stollar, E.; Meirovitch, E.; Kay, L. E. Measurement of Bond Vector Orientations in Invisible Excited States of Proteins. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 18473–18477. (45) Kay, L. E. Protein Dynamics from NMR. Biochem. Cell Biol. 1998, 76, 145–152. (46) Kay, L. E. NMR Studies of Protein Structure and Dynamics. J. Magn. Reson. 2005, 173, 193–207. (47) Wand, A. J. Dynamic Activation of Protein Function: A View Emerging from NMR Spectroscopy. Nat. Struct. Mol. Biol. 2001, 8, 926–931. (48) Br¨ uschweiler, R. New Approaches to the Dynamic Interpretation and Prediction of NMR Relaxation Data from Proteins. Curr. Opin. Struct. Biol. 2003, 13, 175 – 183. 32

ACS Paragon Plus Environment

Page 32 of 51

Page 33 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(49) Boehr, D. D.; Dyson, H. J.; Wright, P. E. An NMR Perspective on Enzyme Dynamics. Chem. Rev. 2006, 106, 3055–3079. (50) Palmer, A. G.; Hochstrasser, R. A.; Millar, D. P.; Rance, M.; Wright, P. E. Characterization of Amino Acid Side Chain Dynamics in a Zinc-Finger Peptide Using Carbon-13 NMR Spectroscopy and Time-Resolved Fluorescence Spectroscopy. J. Am. Chem. Soc. 1993, 115, 6333–6345. (51) Palmer III, A. G. NMR Characterization of the Dynamics of Biomacromolecules. Chem. Rev. 2004, 104, 3623–3640. (52) Sapienza, P. J.; Lee, A. L. Using NMR to Study Fast Dynamics in Proteins: Methods and Applications. Curr. Opin. Pharmacol. 2010, 10, 723–730. (53) Tugarinov, V.; Kay, L. E. Methyl Groups as Probes of Structure and Dynamics in NMR Studies of High-Molecular-Weight Proteins. ChemBioChem 2005, 6, 1567–1577. (54) Lipari, G.; Szabo, A.; Levy, R. M. Protein Dynamics and NMR Relaxation: Comparison of Simulations with Experiment. Nature 1982, 300, 197–198. (55) Ollerenshaw, J. E.; Tugarinov, V.; Skrynnikov, N. R.; Kay, L. E. Comparison of 13 CH3 , 13

CH2 D, and 13 CHD2 Methyl Labeling Strategies in Proteins. J. Biomol. NMR 2005,

33, 25–41. (56) Lee, A. L.; Flynn, P. F.; Wand, A. J. Comparison of 2 H and

13

C NMR Relaxation

Techniques for the Study of Protein Methyl Group Dynamics in Solution. J. Am. Chem. Soc. 1999, 121, 2891–2902. (57) Hansen, D. F.; Neudecker, P.; Vallurupalli, P.; Mulder, F. A.; Kay, L. E. Determination of Leu Side-Chain Conformations in Excited Protein States by NMR Relaxation Dispersion. J. Am. Chem. Soc. 2009, 132, 42–43.

33

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(58) Ishima, R.; Petkova, A. P.; Louis, J. M.; Torchia, D. A. Comparison of Methyl Rotation Axis Order Parameters Derived from Model-Free Analyses of 2 H and 13 C Longitudinal and Transverse Relaxation Rates Measured in the Same Protein Sample. J. Am. Chem. Soc. 2001, 123, 6164–6171. (59) Ruschak, A. M.; Kay, L. E. Methyl Groups as Probes of Supra-Molecular Structure, Dynamics and Function. J. Biomol. NMR 2010, 46, 75–87. (60) Frederick, K. K.; Marlow, M. S.; Valentine, K. G.; Wand, A. J. Conformational Entropy in Molecular Recognition by Proteins. Nature 2007, 448, 325–329. (61) Mittermaier, A.; Kay, L. E.; Forman-Kay, J. D. Analysis of Deuterium RelaxationDerived Methyl Axis Order Parameters and Correlation with Local Structure. J. Biomol. NMR 1999, 13, 181–185. (62) Ishima, R.; Torchia, D. A. Protein Dynamics from NMR. Nat. Struct. Mol. Biol. 2000, 7, 740. (63) Marlow, M. S.; Dogan, J.; Frederick, K. K.; Valentine, K. G.; Wand, A. J. The Role of Conformational Entropy in Molecular Recognition by Calmodulin. Nat. Chem. Biol. 2010, 6, 352–358. (64) Kasinath, V.; Sharp, K. A.; Wand, A. J. Microscopic Insights into the NMR Relaxation-Based Protein Conformational Entropy Meter. J. Am. Chem. Soc. 2013, 135, 15092–15100. (65) Stone, M. J. NMR Relaxation Studies of the Role of Conformational Entropy in Protein Stability and Ligand Binding. Accounts Chem. Res. 2001, 34, 379–388. (66) Prabhu, N. V.; Lee, A. L.; Wand, A. J.; Sharp, K. A. Dynamics and Entropy of a Calmodulin-Peptide Complex Studied by NMR and Molecular Dynamics. Biochemistry 2003, 42, 562–570. 34

ACS Paragon Plus Environment

Page 34 of 51

Page 35 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(67) Li, D.-W.; Showalter, S. A.; Br¨ uschweiler, R. Entropy Localization in Proteins. J. Phys. Chem. B 2010, 114, 16036–16044. (68) DuBay, K. H.; Geissler, P. L. Calculation of Proteins’ Total Side-Chain Torsional Entropy and its Influence on Protein-Ligand Interactions. J. Mol. Biol. 2009, 391, 484–497. (69) Baudry, J. van der Waals Interactions and Decrease of the Rotational Barrier of Methyl-Sized Rotators:A Theoretical Study. J. Am. Chem. Soc. 2006, 128, 11088– 11093. (70) Kitson, D. H.; Hagler, A. T. Theoretical Studies of the Structure and Molecular Dynamics of a Peptide Crystal. Biochemistry 1988, 27, 5246–5257. (71) Krishnan, M.; Smith, J. C. Reconstruction of Protein Side-Chain Conformational Free Energy Surfaces from NMR-Derived Methyl Axis Order Parameters. J. Phys. Chem. B 2012, 116, 4124–4133. (72) Rajitha, R. T.; Smith, J. C.; Krishnan, M. Hidden Regularity and Universal Classification of Fast Side Chain Motions in Proteins. J. Am. Chem. Soc. 2014, 136, 8590–8605. (73) Vijay-Kumar, S.; Bugg, C. E.; Cook, W. J. Structure of Ubiquitin Refined at 1.8 ˚ AResolution. J. Mol. Biol. 1987, 194, 531–544. (74) Perryman, A. L.; Zhang, Q.; Soutter, H. H.; Rosenfeld, R.; McRee, D. E.; Olson, A. J.; Elder, J. E.; David Stout, C. Fragment-Based Screen against HIV Protease. Chem. Biol. Drug Des. 2010, 75, 257–268. (75) Vojtˇechovsk´ y, J.; Chu, K.; Berendzen, J.; Sweet, R. M.; Schlichting, I. Crystal Structures of Myoglobin-Ligand Complexes at Near-Atomic Resolution. Biophys. J. 1999, 77, 2153–2174. 35

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(76) Norledge, B.; Hay, R.; Bateman, O.; Slingsby, C.; Driessen, H. Towards a Molecular Understanding of Phase Separation in the Lens: a Comparison of the X-ray Structures of Two HighTc γ-Crystallins, γE and γF, with Two LowTc γ-Crystallins, γB and γD. Exp. Eye Res. 1997, 65, 609–630. (77) Xu, Z.; Bernlohr, D.; Banaszak, L. The Adipocyte Lipid-binding Protein at 1.6-˚ A Resolution: Crystal Structures of the Apoprotein and with Bound Saturated and Unsaturated Fatty Acids. J. Biol. Chem. 1993, 268, 7874–7884. (78) Noble, M.; Musacchio, A.; Saraste, M.; Courtneidge, S.; Wierenga, R. Crystal Structure of the SH3 Domain in Human Fyn; Comparison of the Three-dimensional Structures of SH3 Domains in Tyrosine Kinases and Spectrin. EMBO J. 1993, 12, 2617– 2624. (79) Hyberts, S.; Goldberg, M.; Havel, T.; Wagner, G. The Solution Structure of Eglin c Based on Measurements of Many NOEs and Coupling Constants and its Comparison with X-ray Structures. Protein Sci. 1992, 1, 736–751. (80) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kal´e, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781–1802. (81) MacKerell, A. D.; Bashford, D.; Bellott,; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S. et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586–3616. (82) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys 1983, 79, 926–935.

36

ACS Paragon Plus Environment

Page 36 of 51

Page 37 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(83) Darve, E.; Pohorille, A. Calculating Free Energies Using Average Force. J. Chem. Phys. 2001, 115, 9169–9183. (84) H´enin, J.; Fiorin, G.; Chipot, C.; Klein, M. L. Exploring Multidimensional Free Energy Landscapes Using Time-Dependent Biases on Collective Variables. J. Chem. Theory Comput. 2009, 6, 35–47. (85) H´enin, J.; Chipot, C. Overcoming Free Energy Barriers Using Unconstrained Molecular Dynamics Simulations. J. Chem. Phys. 2004, 121, 2904–2914. (86) Chipot, C.; Pohorille, A. Free Energy Calculations: Theory and Applications in Chemistry and Biology; Springer GmbH: New York, 2007; Vol. 86; pp 143–149. (87) Li, D.-W.; Br¨ uschweiler, R. NMR-based Protein Potentials. Angew. Chem. 2010, 49, 6778–6780. (88) Lindorff-Larsen, K.; Piana, S.; Palmo, K.; Maragakis, P.; Klepeis, J. L.; Dror, R. O.; Shaw, D. E. Improved Side-Chain Torsion Potentials for the Amber ff99SB Protein Force Field. Proteins 2010, 78, 1950–1958. (89) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Comparison of Multiple Amber Force Fields and Development of Improved Protein Backbone Parameters. Proteins 2006, 65, 712–725. (90) Case, D. A.; Betz, R. M.; Cerutti, D. S.; III, T. E. C.; Darden, T. A.; Duke, R. E.; Giese, T. J.; Gohlke, H.; Goetz, A. W.; Homeyer, N. et al. Ambertools16 ; University of California, San Francisco, 2016. (91) Xue, Y.; Pavlova, M. S.; Ryabov, Y. E.; Reif, B.; Skrynnikov, N. R. Methyl Rotation Barriers in Proteins from 2 H Relaxation Data Implications for Protein Structure. J. Am. Chem. Soc. 2007, 129, 6827–6838.

37

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(92) Roh, J.; Novikov, V.; Gregory, R.; Curtis, J.; Chowdhuri, Z.; Sokolov, A. Onsets of Anharmonicity in Protein Dynamics. Phys. Rev. Lett. 2005, 95, 038101. (93) Krishnan, M.; Smith, J. C. Response of Small-Scale, Methyl Rotors to Protein- Ligand Association: A Simulation Analysis of Calmodulin- Peptide Binding. J. Am. Chem. Soc. 2009, 131, 10083–10091. (94) Radom, L.; Hehre, W. J.; Pople, J. A. Molecular Orbital Theory of the Electronic Structure of Organic Compounds. XIII. Fourier Component Analysis of Internal Rotation Potential Functions in Saturated Molecules. J. Am. Chem. Soc. 1972, 94, 2371– 2381. (95) Tilton Jr, R. F.; Kuntz Jr, I. D.; Petsko, G. A. Cavities in Proteins: Structure of a Metmyoglobin-Xenon Complex Solved to 1.9 ˚ A. Biochemistry 1984, 23, 2849–2857. (96) Brunori, M.; Gibson, Q. H. Cavities and Packing Defects in the Structural Dynamics of Myoglobin. EMBO reports 2001, 2, 674–679. (97) Krishnan, M.; Kurkal-Siebert, V.; Smith, J. C. Methyl Group Dynamics and the Onset of Anharmonicity in Myoglobin. J. Phys. Chem. B 2008, 112, 5522–5533. (98) Glass, D. C.; Krishnan, M.; Smith, J. C.; Baudry, J. Three Entropic Classes of Side Chain in a Globular Protein. J. Phys. Chem. B 2013, 117, 3127–3134. (99) Tzeng, S.-R.; Kalodimos, C. G. Protein Dynamics and Allostery: an NMR view. Curr. Opin. Struct. Biol. 2011, 21, 62–67. (100) LeMaster, D. M. NMR Relaxation Order Parameter Analysis of the Dynamics of Protein Side Chains. J. Am. Chem. Soc. 1999, 121, 1726–1742. (101) Lang, P. T.; Ng, H.-L.; Fraser, J. S.; Corn, J. E.; Echols, N.; Sales, M.; Holton, J. M.; Alber, T. Automated Electron-Density Sampling Reveals Widespread Conformational Polymorphism in Proteins. Protein Sci. 2010, 19, 1420–1431. 38

ACS Paragon Plus Environment

Page 38 of 51

Page 39 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(102) Fraser, J.; van den Bedem, H.; Samelson, A.; Lang, P.; Holton, J.; Echols, N.; Alber, T. Accessing Protein Conformational Ensembles Using Room-Temperature X-ray Crystallography. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 16247–16252. (103) Johnson, E. R.; Keinan, S.; Mori-Sanchez, P.; Contreras-Garcia, J.; Cohen, A. J.; Yang, W. Revealing Noncovalent Interactions. J. Am. Chem. Soc. 2010, 132, 6498– 6506. (104) Contreras-Garc´ıa, J.; Yang, W.; Johnson, E. R. Analysis of Hydrogen-Bond Interaction Potentials from the Electron Density: Integration of Noncovalent Interaction Regions. J. Phys. Chem. A 2011, 115, 12983–12990.

39

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 51

Table 1: Minimum (ξmin ), maximum (ξmax ), mean (hξi) and range of ξ (in kcal/mol units) for different residue types of all systems Residue type ILEδ ILEγ2 LEUδ1 LEUδ2 VALγ1 VALγ2 THRγ2 MET

ξmin ξmax hξi range -0.249 3.737 0.837 3.986 -2.150 2.619 -1.375 4.769 -0.137 1.858 0.644 1.995 -0.33 3.129 0.897 3.459 -2.51 1.943 -1.58 4.453 -2.24 0.76 -1.324 3 -4.30 -0.683 -2.993 3.617 -1.25 -0.531 0.914 0.719

min max Table 2: Minimum (Sconf ), maximum (Sconf ), mean (hSconf i) and range of Sconf (in cal/mol/K units) for different residue types of all systems

Residue type ILEδ ILEγ2 LEUδ1 LEUδ2 VALγ1 VALγ2 THRγ2 MET

min Sconf 4.177 3.797 4.814 4.367 3.728 4.008 3.865 7.762

max Sconf hSconf i 6.939 5.661 6.768 5.828 7.210 6.247 7.224 5.988 6.602 5.671 6.564 5.467 6.601 5.456 9.136 8.644

40

ACS Paragon Plus Environment

range 2.762 2.970 2.395 2.857 2.874 2.553 2.736 1.373

Page 41 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 1: Schematic illustration of the division of the system into sub-system-I and subsystem-II (b) Sub-system-I corresponds to a methyl-containing residue while the remaining residues of the protein, solvents and ions constitute the sub-system-II. The reaction coordinates (θ and φ) that describe methyl dynamics are also shown.

41

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2: Scatter plots of θ/φ obtained from MD trajectories and conformational free energy profiles computed using the ABF method (solid line) shown for representative side chains in the gas phase (red) and in protein (black).

42

ACS Paragon Plus Environment

Page 42 of 51

Page 43 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 3: Schematic illustration of a typical side chain conformational energy surface showing activation barriers and relative rotamer energies.

43

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4: A few representative side chain free energy profiles fitted using piecewise approach are shown.

44

ACS Paragon Plus Environment

Page 44 of 51

Page 45 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 5: Comparison of ∆ξABF obtained directly from ABF-derived free energy profiles with ∆ξS obtained from piecewise integration approach.

45

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 6: A few representative side chain free energy profiles fitted using the truncated Fourier expansion are shown.

46

ACS Paragon Plus Environment

Page 46 of 51

Page 47 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 7: Comparison of ∆ξABF obtained directly from ABF-derived free energy profiles with ∆ξF obtained from truncated Fourier expansion method.

47

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 8: The locations of methyl groups color coded based on their ∆ξ values are shown for all proteins studied.

48

ACS Paragon Plus Environment

Page 48 of 51

Page 49 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2 Figure 9: Dependance of Oaxis on ∆ξ fitted with a sigmoidal function of the form given by Eq. 22.

49

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 10: Dependance of Sconf on ∆ξ fitted with a sigmoidal function of the form given in Eq. 23.

50

ACS Paragon Plus Environment

Page 50 of 51

Page 51 of 51 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 11: TOC Graphic

51

ACS Paragon Plus Environment