Modeling Protein SâAromatic Motifs Reveals Their Structural and

Subscriber access provided by - Access paid by the | UCSB Libraries

B: Biophysical Chemistry and Biomolecules

Modeling Protein S-Aromatic Motifs Reveals Their Structural and Redox Flexibility Esam A. Orabi, and Ann M. English J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b00089 • Publication Date (Web): 13 Mar 2018 Downloaded from http://pubs.acs.org on March 14, 2018

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 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 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.

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 13 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

Modeling Protein S-Aromatic Motifs Reveals Their Structural and Redox Flexibility Esam A. Orabi[a,†] and Ann M. English*[a,b] a.

Department of Chemistry and Biochemistry, Concordia University, 7141 Sherbrooke Street West, Montréal, Québec H4B 1R6, Canada b. Center for Research in Molecular Modeling (CERMM) and Quebec Network for Research on Protein Function, Engineering, and Applications (PROTEO). ABSTRACT: S-aromatic motifs are important noncovalent forces for protein stability and function but remain poorly understood. Hence, we performed quantum calculations at the MP2(full)/6-311++G(d,p) level on complexes between Cys (H2S, MeSH) and Met (Me2S) models with models of Phe (benzene, toluene), Trp (indole, 3-methylindole), Tyr (phenol, 4-methylphenol) and His (imidazole, 4-methylimidazole). The most stable gas-phase conformers exhibit binding energies of −2 to −6 kcal/mol and the S atom lies perpendicular to the ring plane. This reveals preferential interaction with the ring π-system except in the imidazoles where S binds edge-on to an N atom. Complexation tunes the gas-phase vertical ionization potentials of the ligands over as much as 1 eV and strong σ- or π-type H-bonding supports charge transfer to the H-bond donor, rendering it more oxidizable. When the S atom acts as an H-bond acceptor (N/O‒Har···S), calibration of the CHARMM36 force field (by optimizing pair-specific Lennard-Jones parameters) is required. Implementing the optimized parameters in molecular dynamics simulations in bulk water, we find stable Saromatic complexes with binding free energies of −0.6 to −1.1 kcal/mol at ligand separations up to 8 Å. The aqueous S-aromatics exhibit flexible binding conformations but edge-on conformers are less stable in water. Reflecting this, only 0.3 to 10% of the Sindole, S-phenol and S-imidazole structures are stabilized by N/O‒Har···S or S‒H···Oar/Nar σ-type H-bonding. The wide range of energies and geometries found for S-aromatic interactions and their tunable redox properties expose the versatility and variability of the S-aromatic motif in proteins and allow us to predict a number of their reported properties.

INTRODUCTION S-aromatic interactions refer to noncovalent complexation between S-containing and aromatic groups. These interactions are common among small chemical structures1–3 and in proteins.4–14 Protein-based S-aromatic interactions are known as S-aromatic motifs, and they form on association of the Cys and Met side chains with those of Phe, Tyr, Trp or His. Close to 40 years ago it was noted in 21 proteins that contacts between S and aromatic C atoms are more common than expected.4,5 Since then a number of bioinformatic analyses (Table 1) have exposed the high frequency of S-aromatic motifs,6– 14 including those involving cystine residues. A stabilizing role for these motifs has been suggested7–9,13,14 and confirmed in a number of small peptides.15–17 Thus, while studies to date demonstrate that S-aromatic motifs are commonplace and contribute to stabilizing polypeptide secondary structure, the data in Table 1 reveal that they are poorly described. Critical unknowns include factors that influence their geometry (Figure 1) and the stabilities of S-His and S-Tyr interactions in water or proteins have not been assessed although we have modeled S-phenolate and S-imidazolium.18 In this work we focus on modeling S-interactions with the neutral aromatics to delimit their properties and to resolve reported inconsistencies (Table 1). Specifically, the affinity of Cys and Met (modeled with H2S, MeSH and Me2S) for Phe, Trp, Tyr and His (modeled with benzene/toluene, indole/3-

methylindole, phenol/4-methylphenol and imidazole/4methylimidazole) are investigated in the gas phase and in water to cover the range of microenvironments possible in proteins. The S-aromatics were investigated in water by molecular dynamics simulations following calibration of the CHARMM36 all-atom additive force field (FF) where necessary.19

Figure 1. Geometrical parameters for the S-aromatic complexes of (a) benzene, toluene, phenol and 4-methylphenol, (b) indole and 3-methylindole and (c) imidazole and 4-methylimidazole (A = H, CH3; B = H, OH). The S-atom-to-ring-centroid (X) distance is (C6 ring for the indoles). The angle between S, X and any point Y on the vector normal to the ring plane that passes through X is . En-face, edge-on and intermediate refer to binding of the S atom to the aromatic face ( < 30°), edge ( > 60°) or in between (30° ≤ ≤ 60°), respectively. In the H2O-aromatic complexes O replaces S to give and .

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

09 program.35 Over 1400 structures, generated for the H2Oand S-aromatic complexes by placing the interacting molecules at various relative orientations without constrains, converged to 347 stable and structurally different conformers (no imaginary frequencies). Binding energies are corrected for basis set superposition error (BSSE) using the counterpoise (CP) procedure of Boys and Bernardi,36 and corrected (ECP) and uncorrected (E) values are reported. Electrostatic potentials Vs(r) on the surface of H2S and MeSR and their O counterparts (H2O, MeOR) were computed at the same level of theory and the resultant maps shed light on the differences in their aromatic complexes. How complexation alters the redox properties of the interacting ligands in the gas phase also was investigated. Vertical ionization potentials (IPV) were calculated for the optimized geometries of the free ligands and their complexes. The IPV represents the energy difference between a neutral molecule and its radical cation without considering geometry relaxation following electron removal. To assess the site of oxidation in the complexes, the sum of the Mulliken atomic spin densities were calculated for the non-aromatic ligand in each radical cation complex, and to assess the extent of charge transfer in the neutral complexes, we performed natural bond orbital (NBO) analysis on all optimized structures and report the net charge on the non-aromatic ligands. Ab initio rigid potential energy curves (PECs) were generated at the MP2(full)/6-311++G(d,p) level of theory for two conformers of the MeSR complexes with benzene, indole, phenol and imidazole. PECs were calculated by moving the Sligand from 3 to 10 Å in 0.1-Å increments from the point of projection of its S atom on the ring plane or from the ring centroid. The interacting molecules are maintained at their ab initio optimized geometry during the scan and the PECs corrected for BSSE were used to calibrate the FF for the Saromatic interactions. Molecular mechanics (MM) calculations. These were performed with CHARMM,37 using the CHARMM36 all-atom additive protein FF for MeSR and the aromatic compounds,19 the additive model of Kamath et al. for H2S38 and the TIP3P model for H2O.39 MM binding energies using both the default (EMM) and optimized (EMM,opt) FF are calculated as the difference in energy between the complex and the isolated ligands. The geometry of each complex and the relative orientation of the two ligands are constrained in their ab initio gasphase MP2(full)/6-311++G(d,p) optimized structures. FF parameterization. FF calibration is required for the reliable description of certain S-aromatic interactions.10,18,40,41 We reported that the CHARMM36 FF yields binding energies for MeSR-imidazolium and MeSRphenolate in poor agreement with the ab initio values,18 and find here that the FF also underestimates the binding energies of the MeSR complexes of the neutral aromatic with N/O‒ Har···S σ-type H-bonding. Hence, we calibrated the FF for these complexes using our reported procedures18 as outlined in the Supporting Information (SI). Reliable determination of the stability of the S-aromatics in water requires FFs that also correctly describe S-H2O and H2O-aromatic interactions. The default FF gives hydration free

Table 1. Bioinformatics analyses of protein S-aromatic motifs No of proteinsa

Sulfur residue residue

Aromatic Aromatic matic residue residue

CutCut-offb (Ǻ) Ǻ)

rSX or or rSC (Å) (Å)

Geometryc

Ref

36 36

C,M

F,Y,W

11

5.3

edge-on

6

49 49

C

F,H,Y,W

3.8

3.6c

en-face

7

60 60

C-Cd

W

‒

~4

en-face

8

393 393

M

F

3.8

3.6c

en-face

9

4.5‒6.5 5.5 4.5 3.5‒4.0 5.0‒5.5

edge-on/ intermed

10

edge-on

11

en-face edge-on

12

5.0

NRe

13

609 609

C

F

9

604 604

M C-Cd

F,Y

10

753 753

C

F,Y,W

12

~80,000 ~80,000

M

F,Y,W

20

a

Page 2 of 13

b

Number of protein structures analysed. The cut-off is the maximum rSX (see Figure 1) or rSC value considered. cThe rSC is the average separation between S and the aromatic C atoms in the protein S-aromatic motifs examined. dEn-face geometry between the S‒S midpoint and the indole benzene ring. eNot reported.

Reversible and irreversible oxidation of Cys, Met, Tyr, Trp and His is vitally import in protein catalysis, cellular oxidative stress and acceleration of electron transfer through proteins.20–27 Complexation will modulate the redox activity of these residues but this has been largely neglected and seemingly contradictory observations remain unexplained. For example, S···π intramolecular thioether-arene interaction renders the S moiety more prone to oxidation in small molecules28–30 but a proteomics analysis suggests that Met residues in Saromatic motifs are less susceptible to oxidation.31 To address such contradictions, we calculated the ionization potentials of the free ligands and their complexes. Our results provide a quantitative and predictive measure of how complexation modifies the redox reactivity of the ligands and, by extrapolation, of the residues in S-aromatic motifs. In sum, we delimit the binding geometries, stabilities and redox tuning of Saromatic motifs.

COMPUTATIONAL METHODS Gas-phase geometries and binding energies are reported for all stable conformers of the 16 MeSR-aromatics (R = H, CH3) and 16 H2S-/H2O-aromatics examined. High-level ab initio quantum mechanical (QM) calculations accurately predict the strength and directionality of S-aromatic interactions in the gas phase but are computationally prohibitive for studies in bulk water. The performance of the somewhat computationally less-demanding AM132 and PM333 semi-empirical methods were tested on gas-phase H2S-benzene and Me2S-benzene but both methods significantly underestimate the ab initio binding energies of the en-face conformers.34 Adding an empirical dispersion correction term yields the correct H2Sbenzene binding energy but doubles that of Me2S-benzene,34 suggesting that semi-empirical methods inadequately describe S-aromatic interactions. FFs provide the best alternative to the computational cost shortcoming of ab initio QM methods and are used here following calibration to obtain S-aromatic binding free energies in water. Ab initio calculations. These were performed at the MP2(full)/6-311++G(d,p) level of theory using the Gaussian

2


Page 3 of 13 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


energies for the free aromatics in good agreement with experiment42 and only slightly overestimates the binding energy of the global minimum conformers of the H2O-aromatics (see below). Since both the default and our previously optimized FF for MeSR-H2O interactions18 do not reproduce the experimental hydration free energy of MeSR, the FF is readjusted in this work (see SI). Molecular dynamics (MD). Simulations were performed in a box of 1000 water molecules plus one MeSR and one aromatic with cubic periodic boundary conditions in an isothermal-isobaric ensemble (NpT) at 298.15 K and 1.0 atm as described elsewhere.18 Potential of mean force (PMF). Umbrella sampling is used to derive PMFs between MeSR and the aromatics in bulk water as described previously.18 The reaction coordinate (Figure 1) is sampled over 2.5 to 12.0 Å in 0.5-Å increments and the system is simulated for 2.5 ns at each value. The binding affinities for the MeSR-aromatics are calculated using both the default and our optimized FF.

Ab initio electrostatic potential surfaces for numerous aromatics were calculated by Mecozzi et al.45 The benzene and phenol rings exhibit similar negative electrostatic centres but the indole C6 ring possesses a larger and more intense negative electrostatic potential than benzene. In imidazole, the electronegative potential is centered over the pyridine N and is relatively weak over the rest of the ring.45

Investigation of S-aromatic properties in bulk water. The complexation geometry and the bonding interactions that stabilize MeSR-aromatics in bulk water are found from the MD simulations. We investigated N/O‒Har···S and S‒H···Nar/Oar σ-type H-bonding in the complexes using a cutoff distance of 2.8 Å between H and its acceptor atom. This value is slightly larger than the gas-phase H···O/S distances observed for the H-bonded S‒H···O and O‒H···S interactions in the H2S-H2O pair (2.2‒2.6 Å).43 The structures of each complex during a 140-ns MD trajectory (280,000 frames) were analyzed, the values were binned (0.25-Å bin width), and the percentage of H-bonded structures in each bin is plotted vs . The directionalities of the S-aromatic interactions are found by additionally binning the angles (Figure 1) observed in the structures. The number of structures in each bin =( , ) were counted and a plot of the 2D PMF function (>?@ lnA=( , )/2π sin C) vs sin and cos reveals the preference of the S atom for the face or the edge of the aromatic moiety.18

Figure 2. Electrostatic potential maps of (A) H2S, MeSH and Me2S; and (B) H2O, MeOH and Me2O. The geometry of each molecule is presented next to its surface. Two orientations are shown with the S or O atom pointing toward (left side of each panel) and away from the viewer (middle of each panel). Note the much larger potential range (kcal/mol) in the O vs. S ligands. Table 2. Vs,min and Vs,max (kcal/mol) in the electrostatic potential maps of H2S, H2O, MeSR and MeORa Molecule

Vs,min b

Location c

Vs,max b

Location d

H2 S

−15.7 (2)

HSH

18.2 (2)

H

H2 O

−33.8 (1)

HOH

36.3 (2)

H

12.2 (3) 15.3 (1)

methyl H sulfhydryl H Trans methyl H OC bond f hydroxyl H In-plane methyl H g CSC plane h In-plane methyl H i Out-of-plane methyl H

MeSH

−18.4 (2)

CSH

MeOH

−33.6 (1)

COH

Me2S

−20.1 (2)

CSC

7.8 (1) 10.1 (1) 34.3 (1) 10.1 (2) 11.2 (2)

Me2O

−31.7 (1)

COC

9.0 (2) 9.2 (4)

RESULTS Ab initio calculations Electrostatic potential maps. These are useful in analyzing and predicting molecular reactivity.44 The distinctive properties of the O- and S-containing ligands discussed in the sections below arise from considerably larger negative (Vs,min) and positive (Vs,max) electrostatic potentials Vs(r) on the OH vs. SH atoms (Table 2). These values clearly reveal why O is a much better H-bond acceptor and donor than S. Nonetheless, the methyl H atoms of Me2S have the capacity to form stronger H-bonds and electrostatic interactions since their Vs,max values are 1‒2 kcal/mol more positive than those of Me2O. Notably, Vs,min is relatively insensitive to H2O methylation but becomes more negative on H2S methylation (Table 2). Also, it is observed that the S lone pairs are spatially isolated compared to the overlapping electron density of the O lone pairs (Figure 2).

e

g

a

Potentials computed at the MP2(full)/6-311++G(d,p) level. bThe number of Vs,min and Vs,max in each map is in parentheses. cThe two Vs,min in H2S and MeSR are above and below the RSR plane; the single Vs,min of H2O and MeOR is in the ROR plane. dWhere Vs,max is centered on H, this atom is bolded and underlined. eMethyl H in the COH plane. Note no Vs,max is resolved for the other two methyl Hs due to overlap with Vs,max of the hydroxyl H. f Vs,max is along the OC bond, between the three methyl Hs. gThe two methyl Hs in the CSC or COC plane. hThe two Vs,max are above and below the CSC plane. iThe four methyl Hs out of the COC plane.

Overview of geometries and binding energies for the gas-phase S- and H2O-aromatic complex-

3



Page 4 of 13

es. Table 3 summarizes the properties of the 347 stable structures found Table 3. Number (n) of energy-minima en-face, edge-on and intermediate conformers plus binding-energy range (−ECP, kcal/mol) for complexes of H2O, H2S, MeSH and Me2S with the aromaticsa,b,c Aromatic Benzene Toluene Indole 3-Me-indole Phenol d 4-Me-phenol Imidazole 4-Me-imidazole TOTALS

H2S En-face n −ECP

Edge-on n −ECP

En-face n −ECP

Edge-on n −ECP

Intermediate n −ECP

En-face n −ECP

Edge-on n −ECP

Intermediate n −ECP

En-face n −ECP

Edge-on n −ECP

1 5 4 3 3 4 3 2 25

2 8 13 13 10 10 7 6 69

2 3 4 2 5 4 2 3 25

0 3 8 3 1 3 5 2 25

1 4 7 7 9 8 2 5 43

0 0 2 1 1 0 0 1 5

1 5 14 12 9 7 9 8 65

2 7 6 7 5 12 3 3 45

1 1 1 1 1 3 0 0 8

2 3 8 8 5 6 2 2 36

2.2 2.4–2.6 3.1–3.9 3.3–4.2 1.6–2.4 2.7–2.8 2.5–2.8 2.9–3.0

MeSH

0.3–0.6 0.1–0.6 0.5–2.5 0.2–2.5 0.3–2.9 0.2–2.8 0.7–3.8 0.6–4.0

Me2S

2.2–2.6 2.8–3.1 3.1–4.2 4.3 2.1–3.6 3.0–3.6 2.9–3.4 2.9–3.6

1.0–1.1 0.4–1.6 1.0–1.2 1.1 0.4–1.2 1.5–3.8 3.8

2.4 2.4–3.0 2.6–4.7 3.7–5.0 2.0–5.1 2.4–5.1 3.1–3.5 2.4–3.7

H2O

3.6 3.9 2.4

3.3

1.4 1.2–1.6 1.2–4.5 0.9–4.4 0.9–2.5 0.9–2.6 0.5–4.7 0.6–4.7

2.4 2.5–3.3 2.6–5.2 3.9–5.4 2.2–5.9 2.4–5.9 2.7–3.8 2.5–3.8

2.4 3.0 4.0 4.1 2.7 2.7–3.8

1.1–1.3 1.2–1.4 0.3–5.4 0.4–5.2 1.3–6.0 1.2–5.9 5.9–6.7 5.7–6.9

a Properties of the global and local minimum conformers were calculated in the gas phase at the MP2(full)/6-311++G(d,p) level and are reported in Tables 4 and 5 plus S1‒S32. bEn-face ( < 30°), edge-on ( > 60°) or intermediate geometry (30° ≤ ≤ 60°) (Figure 1). cThe –ECP value (corrected for basis set superposition error as outlined in the Computational Methods) for the global minimum conformer of each complex (Tables 4 and 5) is bolded and underlined. dThe H2S-phenol global minimum conformer (omitted from the table due to space constrains) has intermediate geometry with ECP = −3.0 kcal/mol.

(local plus global minimum confers) and their coordinates are listed in the SI. The en-face complexes are mainly stabilized by π-type H-bonding with the S-ligand or H2O acting as an Hbond donor to the ring’s π system (S/O/C‒H···πar) and by S···πar interactions. Edge-on complexes are stabilized by σtype H-bonding with the aromatic accepting (O/N/C‒ Har···S/O) or donating an H-bond (S/O‒H···Nar/Oar). Any Saromatic conformer with O/N‒Har···S or S‒H···Nar/Oar bonding has a complexation energy (-ECP) > 2.0 kcal/mol.

Notably, the lone-pair directionality seen in H2O and H2S (Figure 2) is retained in their complexes. For example, when accepting an H-bond from indole, H2O (4c) and H2S (S3k‒m) form a Z‒O···Har angle of 160° vs a Z‒S···Har angle of 110‒ 120° (Z is a point in the bisector of the HOH or HSH angle). Our ECP values for the H2S-aromatics agree well with the literature values (kcal/mol) for 3a (−2.34, −2.64 and −2.74; CCSD(T)/aug-cc-pVXZ levels, X = D, T, Q),46 3c (−4.49; CCSD(T) level at the complete basis set limit),47 and 3f (−2.78 and −3.68; aug-cc-pVDZ basis set at the B3LYP and MP2 levels, respectively).48

Interactions in the H2S- and H2O-aromatics. Structural data and binding energies for the global minimum conformers of the H2S- and H2O-aromatics are given in Figure 3 and Table 4. Because of the larger atomic radius of S compared to O, rSX is 0.2‒0.6 Å longer than rOX while the larger electrostatic potentials on H2O vs H2S (Table 2) result in stronger H-bonds in the H2O- vs H2S-aromatics. Hence, for structurally similar complexes, the H2O-aromtics are 0.3‒3.0 kcal/mol more stable than the H2S-aromatics. H2S and H2O bind en face to benzene (3a, 4a) and toluene (3b, 4b) via π-type (S/O‒H···πar) H-bonding. Additional N‒Har···S σ-type H-bonding contributes to the greater stability of the distorted en-face H2S-indoles (3c, 3d) whereas N‒ Har···O bonding alone stabilizes the edge-on H2O-indoles (4c, 4d). Both H2S and H2O accept an H-bond from the phenolic O atom (O–Har···S/O) (3e, 3f, 4e, 4f). However, H2S-phenol (3e) exhibits intermediate geometry with both σ- and π-type H-bonding. In contrast, the ligands preferentially donate an Hbond to the imidazoles, giving edge-on conformers with S/O‒ H···Nar binding (3g, 3h, 4g, 4h). Donation of an H-bond to phenolic O (S5j, S6i, S13e, S14h) or acceptance of a H-bond from imidazole N (S7f, S7h, S8c, S8d, S15a, S16a) yields complexes that are 1‒2 kcal/mol less stable than the global minimum complexes. The preference of the aromatic to donate (phenols) or accept (imidazoles) an H-bond reflects its relative basicity, and which prevails in a protein will depend on the local environment.

Interactions in the MeSR- and MeORaromatics. The S‒H···πar bonded en-face MeSH-benzene (5a) and MeSH-toluene (5b) conformers (Figure 4) are ~0.5 kcal/mol more stable (Table 5) than their H2S congeners (3a, 3b). This arises because the more polarizable MeSH ligand exerts larger dispersion forces despite its less positive sulfhydryl H compared to H2S (Vs,max; Table 2). To accommodate N/O‒Har···S σ-type bonding, 5c‒5f exhibit intermediate geometry and again dispersion forces likely contribute to their ~1‒2 kcal/mol higher stability over 3c‒3e (Table 3). The edgeon MeSH-imidazoles with N‒Har···S (5g) and S‒H···Nar (5h) bonding are essentially isoenergetic as are the MeSHimidazole conformers, S23g/S23h, which underscores the flexibility in H-bonding between the S‒H and imidazole moieties. The global minimum Me2S-aromatics (Figure 4) possess intermediate (6a‒6f) or distorted edge-on geometry (6g, 6h). With one exception (6a), these conformers are 0.2‒1.0 kcal/mol more stable than their MeSH analogues (Table 5). This arises from stronger H-bond donation (O/N/C‒Har···S) to Me2S, which has the most negative S atom (Vs,min, Table 2) as well as larger dispersion forces in the Me2S vs MeSH complexes. The MeOR-aromatics were not modeled here. However, the reported ECP values (kcal/mol; MP2/aug-cc-pVDZ level with B3LYP/augg-cc-pVDZ geometries)49,50 for edge-on conformers of the MeOH (−7.25), MeSH (−5.16), Et2O (–9.25)

4


Page 5 of 13 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


gen (blue), oxygen (red) and sulfur (yellow). Black and green dotted lines designate intermolecular σ- and π-type H-bonds, respectively. The atomic coordinates are provided in the SI.

and Et2S (–7.30) complexes of 4-methylphenol, which are stabilized by O‒Har···O/S H-bonding, reflect the greater negative potential on the O vs. S atom (Vs,min; Table 2). Since the methyl Hs are more positive in Me2S vs. Me2O (Vs,max; Table 2), en-face Me2S-benzene (–4.13) stabilized by C‒H···πar Hbonding is more stable than Me2O-benzene (–3.19) (E, kcal/mol at the MP2/6-311+G(d,p) level).13 Vertical ionization potentials and charge transfer. The IPV measures the energy required to remove an electron from an atom or molecule and thus its ease of oxidation. The IPVs calculated for the ligands are close to the experimental values51–60 (Table 6), and are on average within 3.2% of the calculated61 and experimental62 IPVs of the amino acids,

revealing both the reliability of the model chemistry and the suitability of the ligands as models. Methylation lowers the ligand IPV as it increases molecular polarizability and facilitates electron removal. Electron-donation by the methyl also helps stabilize the radical cation formed on oxidation. The data in Tables 6, S1‒S32 indicate that complexation modulates the IPV of the ligands by up to 1eV. For example the IPV of the isolated 3-methylindole is 7.88 eV compared to 8.82 eV in structure 3d of the H2S-3-methylindole complex (Table 6). The high IPV for H2O ensures that the aromatic is the center ionized in all the H2O-aromatics (Table 6, S9‒S16). Similarly the indole ring is the site of oxidation in the H2Sindoles (Table 6, S3, S4). H2S has the highest IPV of the Sligands and is less frequently ionized than MeSH, which in turn is less frequently ionized than Me2S with the lowest IPV (Table 6). However, in the majority of S-aromatics the center ionized depends not only on the ∆IPV of the interacting ligands but also on their complexation mode. For example, when Me2S donates an H-bond in the conformers of Me2S-3methylindole, which model the common Met-Trp motif in proteins, it also gives up the electron (S28k, l, n, o, p, r, s; Figure S28) but in the remaining conformers (including the global minimum conformer 6c; Figure 4), the indole donates

Table 4. Ab initio equilibrium distance ( / , Å), angle ( / degrees) and binding energy (ECP, kcal/mol) of the global minimum conformers of the H2S- and H2O-aromaticsa Conformer/ Conformer/aromatic romatic 3a,4 3a,4a ,4a b benzene enzene 3b,4 3b,4b ,4b t toluene oluene 3c,4 3c,4c ,4c i indole ndole 3d,4 3d,4d ,4d 3 3-methylindole ethylindole 3e,4 3e,4e ,4e p phenol henol 3f,4f 3f,4f 4 f 4-methylphenol ethylphenol 3g,4 3g,4g ,4g i imidazole midazole 3h,4 3h,4h 4-methylimidazole ethylimidazole ,4h 4

H2S

rSXa

θSXY

ECP

3.70 3.68 3.68 3.68 4.25 5.25 4.51 4.54

3 12 22 22 49 90 85 83

–2.17 –2.61 –3.92 –4.18 –3.01 –2.82 –3.84 –3.99

rOX 3.35 3.33 5.19 5.17 4.66 4.66 4.05 4.03

H2O

θOXY

ECP

3.7 10 90 89 90 90 80 77

>2.42 >2.96 >5.37 >5.20 >6.02 >5.86 >6.71 >6.88

a

See structures in Figure 3. ECP values are corrected for basis set superposition error (BSSE) as outlined in the Computational Methods.

Figure 3. Optimized gas-phase geometry of the global minimum conformer of the H2S and the H2O complex with benzene (3a ,4a), toluene (3b, 4b), indole (3c, 4c), 3-methylindole (3d, 4d), phenol (3e, 4e), 4-methylphenol (3f, 4f), imidazole (3g, 4f) and 4methylimidazole (3h, 4h) at the MP2(full)/6-311++G(d,p) level of theory. Atom color code: hydrogen (white), carbon (gray), nitro-

Figure 4 Optimized gas-phase geometry of the global minimum conformer of the MeSH and MeSR complexes of benzene (5a, 6a) toluene (5b, 6b), indole (5c, 6c), 3-methylindole (5d, 6d) phenol

5



a

See structures in Figure 4. ECP values corrected for basis set superposition error (BSSE) as outlined in the Computational Methods.

(5e, 6e), 4-methylphenol (5f, 6f), imidazole (5g, 6g) and 4methylimidazole (5h, 6h). See caption of Figure 3 for details. Table 5. Ab initio equilibrium distance (KLM , Å), angle (NLMO , degrees) and binding energy (ECP, kcal/mol) of the global minimum conformers of the MeSH- and Me2S-aromatic complexesa Conformer/Aromatic Conformer/Aromatic /Aromatic 5a, 6a , 6a Benzene Benzene 5b, 6b , 6b Toluene Toluene 5c, 6c , 6c Indole Indole 5d, 6d , 6d 3-Methylindole Methylindole 5e, 6e , 6e Phenol Phenol 5f, 6f , 6f 4-Methylphenol Methylphenol 5g, 6g , 6g Imidazole Imidazole 5h, 6h , 6h 4-Methylimidazole Methylimidazole

MeSH MeSH rSX 3.71 3.66 4.02 4.01 4.15 4.11 4.51 4.42

θSXY 1.5 13 32 32 47 46 76 88

ECP >2.59 >3.14 >4.74 >5.04 >5.05 >5.12 >3.77 >3.80

Me2S rSX 4.16 4.10 4.13 4.13 4.52 4.49 4.43 4.41

θSXY 33 31 41 40 56 54 73 71

Page 6 of 13

ECP >2.41 >3.33 >5.22 >5.38 >5.89 >5.92 >4.74 >4.68

Table 6. Calculated IPVs (eV) for the free ligands and the H2O- and S-aromatics a Ligand Ligand

Free ligand ligand

Expt b

Ref Ref

H2O

12.60

12.65

51

H2S

9.99

10.46

52

H2O-aromatic c

H2S-aromatic c

MeSHMeSH-aromatic c

Me2S-aromatic c

MeSH MeSH

9.17

9.45

53

Me2S

8.53

8.70

54

Benzene Benzene

9.51

9.24

55

9.97 (9.18 ‒ 9.97) 9.97

9.38 (9.35 ‒ 9.67)

8.68 (8.68 ‒ 8.83)

8.19 (8.19 ‒ 8.49)

Toluene Toluene

9.45

8.83

56

9.78 (8.97 ‒ 9.78) 9.78

9.37 (9.15 ‒ 9.59)

8.69 (8.65 ‒ 9.15)

8.18 (8.14 ‒ 8.51)

Indole Indole

8.19

7.76

57

7.85 (7.85 7.85 ‒ 8.52)

8.40 (8.04 8.40 ‒ 8.86)

8.36 (8.20 8.36 ‒ 8.94)

8.03 (7.89 8.03 ‒ 8.76)

3-Methylindole Methylindole

7.88

7.51

57

7.53 (7.53 ‒ 8.24) 7.53

8.82 (7.74 ‒ 8.82) 8.82

9.02 (8.13 ‒ 9.43) 9.02

8.62 (7.78 ‒ 8.62) 8.62

Phenol Phenol

9.20

8.51

58

8.73 (8.73 ‒ 9.53) 8.73

9.22 (8.91 ‒ 10.04) 9.22

9.21 (8.65 ‒ 9.21) 9.21

9.11 (8.10 ‒ 9.11) 9.11

4-Methylphenol Methylphenol

8.77

8.34

59

8.36 (8.36 ‒ 9.09) 8.36

8.53 (8.53 ‒ 9.60) 8.53

8.81 (8.58 ‒ 8.99) 8.81

8.67 (8.07 ‒ 8.81) 8.67

Imidazole Imidazole

9.18

8.81

60

9.52 (8.69 ‒ 9.52) .52

9.10 (8.87 ‒ 9.59)

8.91 (8.40 ‒ 9.22) 8.91

8.91 (7.97 ‒ 8.91) 8.91

4-Methylimidazole Methylimidazole

8.78

NRb

9.10 (8.32 ‒ 9.10) 9.10

9.01 (8.50 ‒ 9.46)

8.46 (8.26 ‒ 9.14)

8.55 (7.91 ‒ 8.75) 8.55

b

a

Vertical ionization potentials (IPVs) calculated at the MP2(full)/6-311++G(d,p) level. Experimental IPVs reported for the free ligands (NR, not reported). c IPV for the global minimum conformer (Figure 3 and 4) and IPV range for the energy-minima conformers (SI Figures S1-S32 and Tables S1-S32) of each complex is in parentheses. The bolded red and green font indicates removal of 1.0 electron and 0.3-0.4 electron from the aromatic ligand (see text).

Table 7. The rSX (Å) values and PMF minima (kcal/mol) for the MeSR-aromatics in bulk water Aromatic Benzene Toluene Indole 3-Methylindole Phenol 4-Methylphenol Imidazole 4-Methylimidazole

MeSH Optimized FF b KUV 4.8–5.3 4.8–5.2 4.9–5.3 4.9–5.3 4.7–5.1 4.7–5.2 4.3–4.6 4.3–4.4

PMF −0.7 −0.7 −0.6 −0.6 −0.8 −1.0 −0.8 −0.9

Default FF b KUV 4.7–5.3 4.9–5.2 4.8–5.2 4.9–5.2 4.9–5.1 4.6–5.3 4.3–4.7 4.3–4.7

PMF −0.6 −0.7 −0.5 −0.7 −0.7 −0.6 −0.7 −0.5

Me2S Optimized FF b PMF KUV 4.3–5.2 −0.6 4.6–5.1 −0.6 3.8–5.1 −0.9 4.6–5.3 −0.8 3.6–5.3 −0.9 3.4–5.1 −1.1 4.3–4.6 −1.0 4.3–4.7 −1.0

Default FF b KUV 4.8–5.4 5.1–5.3 4.9–5.3 4.9–5.4 4.8–5.4 4.8–5.4 4.4–4.8 4.4

PMF −0.8 −0.6 −0.8 −0.7 −0.6 −0.7 −0.6 −0.6

a

Data from Figure S35. The PMF values give the binding free energies (kcal/mol). bValues calculated using the optimized and default FF parameters, respectively (see text).

an H-bond to the S atom and gives up the electron (Table S28). Thus, if ∆IPV is not too large, the H-bond-donor generally becomes the electron donor. In certain MeSR-aromatics (e.g., S18a, S21a, S26a, S29a) ionization of the aromatic would strengthen C‒Har···S bonding but would repel the elec-

tropositive MeSR methyl Hs so the MeSR fragment, which has the lower IPV, is oxidized. The extent of charge transfer in the complexes is small as revealed by the low net charge on the non-aromatic ligand calculated from NBO analyses (Tables S1‒S32). Complexes that exhibit the most charge transfer have strong N/O‒

6


Page 7 of 13 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


6).18 Although the indoles become unstable at > 70° (Figure 6B,F), the remaining MeSR-aromatics exhibit full conformational flexibility, ranging from en-face to edge-on geometries. Nonetheless, en-face geometry is preferred with a few exceptions: the distorted en-face/intermediate geometry of Me2S-phenol (Figure 6G) and the intermediate geometry of the MeSR-imidazole (Figure 6D,H). The broad minimum in the 1D PMF of the indoles and phenols (Figure S35) splits into two minima at rSX < 4 Å and ~ 5 Å in the 2D PMF map (Figure 6). The presence of a minimum at rSX < 4 Å suggests that S···πar interaction contributes to stabilizing these aqueous complexes. Interestingly, en-face binding is not favored at rSX ~6 Å (Figure 6) since the gap between the ligands is too small to accommodate a water molecule. Nonetheless, the aqueous MeSR-aromatics persist up to rSX ~ 8 Å and a weak shoulder is observed in some PMFs at ~ 8.5 Å (Figure S35), corresponding to a solvent-separated pair. Hence, S-aromatic interactions are long-range in water despite being weakened by competing interactions between the ligands and solvent.

Har···S/O or S/O‒H···Nar/Oar σ-type H-bonding (e.g., 3e-h, 5eh, 6e-h). Besides strengthening complexation, significant charge transfer to the H-bond donor render it the site oxidized. Stability of the S-aromatics in water. The preference of an aromatic for MeSR vs. H2O can be predicted from the gas-phase binding energies of the global minimum conformers of the MeSR-aromatics and those of the competing MeSR-H2O and H2O-aromatic complexes (Table S33). However, such predictions neglect ligand hydration and do not consider the formation of local-minima conformers. Thus, the stability of the MeSR-aromatics was rigorously examined in bulk water using MD simulations. We calibrated the FF for the MeSR-aromatics using the optimized LJ parameters in Table S34, which yields binding energies (EMM,opt) closer to the ab initio values (ECP) than the default FF (EMM) (Tables S17‒S32; Figure S33 and S34). The improvement in calculating the interaction energy using the optimized FF is seen by comparing the mean unsigned and signed errors (MUE and MSE) calculated with both the default and optimized FFs for each S- and H2O-aromatic complex (Table S35). Since our previously optimized LJ parameters for MeSRH2O18 overestimate the MeSR hydration free energy, they were readjusted to produce values that agree better with experiment (Table S36). The calculated PMFs for the complexes in bulk water with the optimized FF (Figure S35) give binding free energies of −0.6 to −1.1 kcal/mol for the MeSR-aromatics (Table 7). For most of the complexes, the optimized FF predicts larger binding free energies compared to the default FF (Table 7). The Me2S complexes of the polar aromatics are 0.1‒0.3 kcal/mol more stable than their MeSH congeners and ring methylation alters complex stability by ≤ 0.1 kcal/mol except it strengthens MeSR-phenol by 0.2 kcal/mol (Table 7).

DISCUSSION Properties of the gas-phase S-aromatics. The energetic and structural properties of S-aromatic complexes in the gas phase have been the subject of numerous computational investigations.10–13,18,34,40,41,46–49,63–71 However, these studies have been limited in scope, reporting on global-minimum structures only and many focus on H2S-benzene.12,41,46–48,65–68 To understand the chemistry of protein S-aromatic motifs, we extend the study to MeSR interactions with benzene, indole, phenol, imidazole and their methylated forms, the latter being better models for residue side chains. Our extensive search in the configurational space of the H2S and MeSR complexes with the eight aromatics uncovers 303 stable conformers plus 44 conformers for the H2O-aromatics. The more stable conformers are characterized by en-face or intermediate geometry except when edge-on binding is stabilized by H-bonding to an aromatic heteroatom. Our calculated IPVs of the free ligands and their complexes predict how S-aromatic interactions will influence their redox properties. Few experimental studies address this issue but spectroscopic and electrochemical data reveal that constraining a thioether near an aromatic ring depresses its oxidation potential.28–30 In agreement, Me2S is the center ionized in Me2S-benzene and Me2S-toluene and with S close to the aromatic ring (rSX ~ 4 Å), Me2S exhibits IPVs 0.29‒0.39 eV lower than free ligand (Table 6, S25, S26). Most notably, our analysis of Mulliken atomic spin densities in the S-aromatic radical cations reveals that the H-bond donor is typically the electron donor when the ligands possess similar IPVs (e.g., 3f, S1c, S2a, S5m, S7f, S8c). Properties of the aqueous MeSR-aromatics. Affinities in the range of −0.6 to −1.1 kcal/mol are found for the aqueous MeRS-aromatics from the MD simulations (Table 7). The optimized parameters result in up to 0.4 kcal/mol larger binding affinities for the complexes in water but the difference will be much greater for the complexes in hydrophobic environments, where H-bonding interactions are much more sta-

Geometry and interactions in the MeSRaromatics in water. The aqueous complexes exhibit rSX values centered at ~4.5 or ~5 Å (Table 7), which indicate direct complexation between the ligands. The broad minimum displayed by the Me2S-indoles and Me2S-phenols is striking but it is not replicated by the default parameters (Figure S35C,D). Since N/O‒Har···S and S‒H···Nar/Oar σ-type H-bonding stabilizes the gaseous complexes of the polar aromatics (Figure 4), their persistence in water is of interest. Calculating the fraction of H-bonded MeSR-aromatics reveals that ~3% of MeSR-indole (rSX ~ 5.2 Å), ~10% of MeSR-phenol (rSX ~ 4.9 Å) and ~10% of MeSR-imidazole (rSX ~ 4.4 Å) structures exhibit N/O–Har···S H-bonding (Figure 5A,C). The default parameters underestimate the stability of this type of H-bonding (Tables S21‒S24 and S29‒S32) and hence the percentage of structures stabilized by it in water (Figure 5B,D). The presence of S‒H···Nar/Oar H-bonding only in 4% of MeSH-imidazole and 0.3% of MeSH-phenol structures (Figure 5A) reveals the relatively weak H-bond donating capability of MeSH. Figure 5 also reveals that H-bonding in MeSRimidazole is restricted to a narrow range of rSX (4‒5 Å). Bonding directionality is visualized in 2D plots that separately show the PMF dependence on and (Figure

7



ble. The complexes exhibit conformational flexibility (Figure 6)

Page 8 of 13

Comparing our results to the available experimental values, we find that the calculated binding free energies for MeSR-toluene in water (−0.7 and −0.6 kcal/mol; Table 7) are in good agreement with those reported for the Cys-Phe (–0.6 to –2.0 kcal/mol)15 and Met-Phe (–0.5 to –0.8 kcal/mol)16 interaction at positions i, i + 4 of α-helices. In contrast, the MetPhe or Met-Trp interaction that stabilizes a β-hairpin peptide is reported to be only –0.3 kcal/mol,17 which is a half to a third of the values found here (Table 7). Our MD simulations indicate that the S atom in MeSH preferentially acts an H-bond acceptor in water (Figure 5). In agreement, solvent-exposed Cys is reported to be selective toward H-bond donors in over 9,000 nonredundant protein structures.72 In contrast, from the analysis of 500 highresolution protein structures, only 151 Cys are stated to be Hbond acceptors vs. 465 H-bond donors. However, H-bonding of Cys to aromatic N (His and Trp) and hydroxyl O (Tyr, Ser, Thr) was rarely found and solvent exposure was not considered in this study.73 A survey of human proteins to explore the structural determinants of Met oxidation revealed that a buried or solventexposed Met close to Trp or Tyr is less prone to oxidation.31 This led the authors to propose that the S-aromatic motif is a key determinant of Met redox status.31 If, as seen here for Me2S-3-methylindole and Me2S-4-methylphenol (Figure 4), N/O–Har···S H-bonding stabilizes protein-based Met-Trp and Met-Tyr interactions, then the aromatic should be the center oxidized (Table 6) since it is the H-bond donor. Although they were not surveyed,31 the structures of the S-aromatic motifs are of interest to establish the predictive power of our model chemistry. Also, since N/O–Har···S H-bonding is stronger in the gas phase than in water (Table 5 vs. 7), the protein microenvironment should control how well Met is protected in its Saromatic motifs. Recently, it has been proposed that the Trp/Tyr chains abundant in many oxidoreductases protect their active sites by transporting holes to the protein surface.20–27,74 Hence, the ability of S-aromatic motifs to direct hole-hopping to specific Trp and Tyr residues in proteins is of interest. For example, it is well established in cytochrome c peroxidase that interactions with M230 and/or M231 stabilize a cation radical on W191, which is essential for hole-hopping to ferrocytochrome c and hence for peroxidase activity.75 S-aromatic motifs additionally promote protein-protein or protein-ligand interactions. The former is exemplified by the Met-Trp motif that stabilizes calmodulin binding to its target proteins76 and the latter by > 25,000 PDB structures of protein-ligand complexes where Cys and Met interact with ligand aromatic rings.77 Indeed, numerous studies13,14,78–84 report on protein-ligand complexes stabilized by S-aromatic interactions where the aromatic is protein or ligand based. Met oxidation to the sulfoxide (MetO) increases its affinity for aromatic residues both in the gas phase and in a model peptide scaffold.85 The structural and energetic properties of MetO-aromatics have been examined for benzene, phenol and indole rings but the geometry and redox properties of these complexes have not been characterized.

Figure 5. Percentages of N/O–Har···S (solid lines) and S‒ H···Nar/Oar (dotted lines) σ-type H-bonded structures in (A,B) MeSH and (C,D) Me2S complexes of indole, phenol or imidazole in bulk water at 298.15 K vs. rSX. Plots are constructed from 140ns MD simulations of the MeSR-aromatics using the (A,C) optimized and (B,D) default FF parameters; rSX is defined in Figure 1.

with rSX values centered around 4.5 to 5 Å but they persist up to rSX ~ 8 Å in water (Figure S35). The edge-on conformers are relatively weak (Figure 6) since H2O strongly competes with MeSR for H-bonding to the aromatic heteroatoms (Table S33). Hence, only ~10% of MeSR-imidazole and MeSRphenol structures and as low as 3% of MeSR-indole structures are stabilized by H-bond acceptance by MeSR (N/O–Har···S). MeSH donates an H-bond (S‒H···Nar/Oar) even less frequently, stabilizing only 4% of MeSH-imidazole and 0.3% of MeSHphenol structures. The minima at rSX < 4 Å in the 2D PMF maps (Figure 6) indicate that S···πar interactions additionally contribute to stabilizing the aqueous complexes. We reported previously that MeSR-phenolate and MeSR-imidazolium preferentially adopt en-face and intermediate geometry, respectively.18 Again, the geometry of the complexes of these charged aromatics is dictated by the higher affinity of the heteroatoms of the charged aromatics for H2O vs. MeSR.

Insights into S-aromatic motifs in proteins. Our evaluation of S-aromatic binding in the gas phase and bulk water delimits the properties of S-aromatic motifs in hydrophobic and polar protein microenvironments. The binding affinities found for the MeSR-aromatics range from 0.6−6 kcal/mol between water and the gas-phase (Table 5 and 7) so a 10-fold variation in stability is possible for S-aromatic motifs depending on the protein microenvironment alone. H-bonding to neighboring residues and/or ligand ionization18 could further enhance the stability of these motifs. Also, since Saromatic interactions persist over ~4−8 Å (Figure S35), Saromatic motifs likely contribute to protein stability over a wide range of distances as well as bonding configurations.

8


Page 9 of 13 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


Our comprehensive calculations underscore the distinctive properties of the S- vs. O-containing ligands, which account for the special behavior of S-aromatic motifs. We delimit the physical properties of S-aromatic interactions in the gas phase and in water as a basis for further examination of the role of S-aromatic motifs in proteins. Finally, it should be stressed that the default CHARMM36 FF underestimates by 25‒55% the gas-phase binding energies of complexes with N/O‒Har···S bonding whereas our optimized FF offers a reliable tool for exploring S-aromatic motifs in proteins, especially those in nonpolar environments where H-bonding becomes more dominant.

ASSOCIATED CONTENT Electronic Supplementary Information (ESI) available: 32 figures showing the optimized geometry of the local minima conformers of all complexes; 32 tables listing the structural and energetic parameters as well as the IPV, and the net charge and sum of Mulliken atomic spin densities on the non-aromatic ligand of the global plus local minima conformers of all complexes; tables of data for predicting the stability of the S-aromatics in water, with the optimized force-field parameters, with mean unsigned and signed errors for FF-calculated binding energies, and with the components of the hydration free energy; figures showing the potential energy curves for the gaseous MeSR-aromatics and the PMFs for the aqueous complexes; and the atomic coordinates of all (local plus global) optimized geometries. This material is available free of charge via the Internet at http://pubs.acs.org.

AUTHOR INFORMATION Corresponding Author * E-mail: [email protected]. Tel.: +1-514-848-2424, extension 3338; Fax: +1-514-848-2868

Notes † On leave from the Department of Chemistry, Faculty of Science, Assiut University, Assiut 71516, Egypt.

Figure 6. 2D maps of the potential of mean force (PMF, kcal/mol) for the MeSH (left panels) and Me2S complexes (right panels) of (A,E) benzene, (B,F) indole, (C,G) phenol and (D,H) imidazole. Maps are constructed from 140-ns MD simulations of the MeSR-aromatics in bulk water at 298.15 K. The PMF is the function: -?@ lnA=( , )/2π sin C, where and are defined in Figure 1. Note that the S atom is located en-face to the aromatic ring when = 0° (x-axis) and edge-on to the ring plane when = 90° (y-axis).

ACKNOWLEDGMENTS This work was supported by a research grant from the Natural Sciences and Engineering Research Council (NSERC) of Canada awarded to AME. Computations were enabled by the support provided by the Centre for Research in Molecular Modeling (CERMM) at Concordia and by Compute Canada (www.computecanada.ca).

REFERENCES (1) Zauhar, R. J.; Colbert, C. L.; Morgan, R. S.; Welsh, W. J. Evidence for a Strong Sulfur-Aromatic Interaction Derived from Crystallographic Data. Biopolymers 2000, 53, 233–248. (2) Wan, C.-Q.; Han, J.; Mak, T. C. W. Intermolecular S···π interactions in Crystalline Sulfanyl-Triazine Derivatives. New J. Chem. 2009, 33, 707–712. (3) Chylewska, A.; Sikorski, A.; Ogryzek, M.; Makowski, M. Attractive S···π and π-π Interactions in the Pyrazine-2Thiocarboxamide Structure: Experimental and Computational Studies in the Context of Crystal Engineering and Microbiological Properties. J. Mol. Struct. 2016, 1105, 96–104. (4) Warme, P. K.; Morgan, R. S. A Survey of Atomic Interactions in 21 Proteins. J. Mol. Biol. 1978, 118, 273–287.

9



(5) Warme, P. K.; Morgan, R. S. A Survey of Amino Acid Side-Chain Interactions in 21 Proteins. J. Mol. Biol. 1978, 118, 289– 304. (6) Reid, K. S. C.; Lindley, P. F.; Thornton, J. M. SulphurAromatic Interactions in Proteins. FEBS Lett. 1985, 190, 209–213. (7) Pal, D.; Chakrabarti, P. Different Types of Interactions Involving Cysteine Sulfhydryl Group in Proteins. J. Biomol. Struct. Dyn. 1998, 15, 1059–1072. (8) Ioerger, T. R.; Du, C.; Linthicum, D. S. Conservation of Cys-Cys Trp Structural Triads and Their Geometry in the Protein Domains of Immunoglobulin Superfamily Members. Mol. Immunol. 1999, 36, 373–386. (9) Pal, D.; Chakrabarti, P. Non-Hydrogen Bond Interactions Involving the Methionine Sulfur Atom. J. Biomol. Struct. Dyn. 2001, 19, 115–128. (10) Duan, G.; Smith Jr, V. H.; Weaver, D. F. Characterization of Aromatic-Thiol π-Type Hydrogen Bonding and PhenylalanineCysteine Side Chain Interactions Through Ab Initio Calculations and Protein Database Analyses. Mol. Phys. 2001, 99, 1689–1699. (11) Iwaoka, M.; Takemoto, S.; Okada, M.; Tomoda, S. Weak Nonbonded S···X (X = O, N, and S) Interactions in Proteins. Statistical and Theoretical Studies. Bull. Chem. Soc. Jpn. 2002, 75, 1611–1625. (12) Ringer, A. L.; Senenko, A.; Sherrill, C. D. Models of S/π Interactions in Protein Structures: Comparison of the H2S‒Benzene Complex with PDB Data. Protein Sci. 2007, 16, 2216–2223. (13) Valley, C. C.; Cembran, A.; Perlmutter, J. D.; Lewis, A. K.; Labello, N. P.; Gao, J.; Sachs, J. N. The Methionine-Aromatic Motif Plays a Unique Role in Stabilizing Protein Structure. J. Biol. Chem. 2012, 287, 34979–34991. (14) Cordomí, A.; Gómez-Tamayo, J. C.; Gigoux, V.; Fourmy, D. Sulfur-Containing Amino Acids in 7TMRs: Molecular Gears for Pharmacology and Function. Trends Pharmacol. Sci. 2013, 34, 320– 331. (15) Viguera, A. R.; Serrano, L. Side-Chain Interactions between Sulfur-Containing Amino Acids and Phenyalanine in αHelices. Biochemistry 1995, 34, 8771–8779. (16) Stapley, B. J.; Rohl, C. A.; Doig, A. J. Addition of Side Chain Interactions to Modified Lifson-Roig Helix-Coil Theory: Application to Energetics of Phenylalanine-Methionine Interactions. Protein Sci. 1995, 4, 2383–2391. (17) Tatko, C. D.; Waters, M. L. Investigation of the Nature of the Methionine-π Interaction in β-Hairpin Peptide Model Systems. Protein Sci. 2004, 13, 2515–2522. (18) Orabi, E. A.; English, A. M. Sulfur-Aromatic Interactions: Modeling Cysteine and Methionine Binding to Tyrosinate and Histidinium Ions to Assess Their Influence on Protein Electron Transfer. Isr. J. Chem. 2016, 56, 872–885. (19) Best, R. B.; Zhu, X.; Shim, J.; Lopes, P. E. M.; Mittal, J.; Feig, M.; Mackerell, A. D.; Jr. Optimization of the Additive CHARMM All-Atom Protein Force Field Targeting Improved Sampling of the Backbone Φ, ψ and Side-Chain χ(1) and χ(2) Dihedral Angles. J. Chem. Theory Comput. 2012, 8, 3257–3273. (20) Morgan, R. S.; Tatsch, C. E.; Gushard, R. H.; Mcadon, J. M.; Warme, P. K. Chains of Alternating Sulfur and π-Bonded Atoms in Eight Small Proteins. Int. J. Pept. Protein Res. 1978, 11, 209–217. (21) Lawrence, J. D.; Swenson, R. P. Role of Methionine 56 in the Control of the Oxidation−Reduction Potentials of the Clostridium Beijerinckii Flavodoxin: Effects of Substitutions by Aliphatic Amino Acids and Evidence for a Role of Sulfur−Flavin Interactions. Biochemistry 1998, 37, 9668–9678. (22) Pfister, T. D.; Gengenbach, A. J.; Syn, S.; Lu, Y. The Role of Redox-Active Amino Acids on Compound I Stability, Substrate Oxidation, and Protein Cross-Linking in Yeast Cytochrome c Peroxidase. Biochemistry 2001, 40, 14942–14951.

Page 10 of 13

(23) Shih, C.; Museth, A. K.; Abrahamsson, M.; BlancoRodriguez, A. M.; Di Bilio, A. J.; Sudhamsu, J.; Crane, B. R.; Ronayne, K. L.; Towrie, M.; Vlček, A. et al. Tryptophan-Accelerated Electron Flow Through Proteins. Science 2008, 320, 1760–1762. (24) Wang, M.; Gao, J.; Müller, P.; Giese, B. Electron Transfer in Peptides with Cysteine and Methionine as Relay Amino Acids. Angew. Chemie Int. Ed. 2009, 48, 4232–4234. (25) Dempsey, J. L.; Winkler, J. R.; Gray, H. B. Proton-Coupled Electron Flow in Protein Redox Machines. Chem. Rev. 2010, 110, 7024–7039. (26) Takematsu, K.; Williamson, H.; Blanco-Rodríguez, A. M.; Sokolová, L.; Nikolovski, P.; Kaiser, J. T.; Towrie, M.; Clark, I. P.; Vlček Jr, A.; Winkler, J. R. et al. Tryptophan-Accelerated Electron Flow Across a Protein‒Protein Interface. J. Am. Chem. Soc. 2013, 135, 15515–15525. (27) Gray, H. B.; Winkler, J. R. Hole Hopping Through Tyrosine/Tryptophan Chains Protects Proteins from Oxidative Damage. Proc. Natl. Acad. Sci. U. S. A. 2015, 112, 10920–10925. (28) Chung, W. J.; Ammam, M.; Gruhn, N. E.; Nichol, G. S.; Singh, W. P.; Wilson, G. S.; Glass, R. S. Interactions of Arenes and Thioethers Resulting in Facilitated Oxidation. Org. Lett. 2009, 11, 397–400. (29) Ammam, M.; Zakai, U. I.; Wilson, G. S.; Glass, R. S. Anodic oxidation of m-Terphenyl Thio-, Seleno- and Telluroethers: Lowered Oxidation Potentials Due to Chalcogen⋯π Interaction. Pure Appl. Chem. 2010, 82, 555–563. (30) Monney, N. P.-A.; Bally, T.; Bhagavathy, G. S.; Glass, R. S. Spectroscopic Evidence for a New Type of Bonding between a Thioether Radical Cation and a Phenyl Group. Org. Lett. 2013, 15, 4932–4935. (31) Aledo, J. C.; Cantón, F. R.; Veredas, F. J. Sulphur Atoms from Methionines Interacting with Aromatic Residues Are Less Prone to Oxidation. Sci. Rep. 2015, 5, 16955. (32) Dewar, M. J. S.; Zoebisch, E. G.; Healy, E. F.; Stewart, J. J. P. Development and Use of Quantum Mechanical Molecular Models. 76. AM1: A New General Purpose Quantum Mechanical Molecular Model. J. Am. Chem. Soc. 1985, 107, 3902–3909. (33) Stewart, J. J. P. Optimization of Parameters for Semiempirical Methods I. Method. J. Comput. Chem. 1989, 10, 209– 220. (34) Morgado, C. A.; ; McNamara, J. P.; Hillier, I. H.; Burton, N. A.; Vincent, M. A. Density Functional and Semiempirical Molecular Orbital Methods Including Dispersion Corrections for the Accurate Description of Noncovalent Interactions Involving SulfurContaining Molecules. J. Chem. Theory Comput. 2007, 3, 1656–1664. (35) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A. et al. Gaussian, Inc., Wallingford CT, Gaussian 09, Revision E.01, 2013. (36) Boys, S. F.; Bernardi, F. The Calculation of Small Molecular Interactions by the Differences of Separate Total Energies. Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553– 566. (37) Brooks, B. R.; Brooks, C. L.; Mackerell, A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S. et al. CHARMM: The Biomolecular Simulation Program. J. Comput. Chem. 2009, 30, 1545–1614. (38) Kamath, G.; Lubna, N.; Potoff, J. J. Effect of Partial Charge Parametrization on the Fluid Phase Behavior of Hydrogen Sulfide. J. Chem. Phys. 2005, 123, 124505. (39) 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.

10


Page 11 of 13 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


(40) Pranata, J. Sulfur–Aromatic Interactions: A Computational Study of the Dimethyl Sulfide–Benzene Complex. Bioorg. Chem. 1997, 25, 213–219. (41) Sherrill, C. D.; Sumpter, B. G.; Sinnokrot, M. O.; Marshall, M. S.; Hohenstein, E. G.; Walker, R. C.; Gould, I. R. Assessment of Standard Force Field Models against High-Quality Ab Initio Potential Curves for Prototypes of π–π, CH/π, and SH/π Interactions. J. Comput. Chem. 2009, 30, 2187–2193. (42) Macias, A. T.; MacKerell, A. D. CH/π interactions involving aromatic amino acids: Refinement of the CHARMM Tryptophan Force Field. J. Comput. Chem. 2005, 26, 1452–1463. (43) Orabi, E. A.; Lamoureux, G. Simulation of Liquid and Supercritical Hydrogen Sulfide and of Alkali Ions in the Pure and Aqueous Liquid. J. Chem. Theory Comput. 2014, 10, 3221–3235. (44) Politzer, P.; Laurence, P. R.; Jayasuriya, K. Molecular Electrostatic Potentials: An Effective Tool for the Elucidation of Biochemical Phenomena. Environ. Health Perspect. 1985, 61, 191– 202. (45) Mecozzi, S.; West, A. P.; Dougherty, D. A. Cation-π Interactions in Aromatics of Biological and Medicinal Interest: Electrostatic Potential Surfaces as a Useful Qualitative Guide. Proc. Natl. Acad. Sci. U. S. A. 1996, 93, 10566–10571. (46) Tauer, T. P.; Derrick, M. E.; Sherrill, C. D. Estimates of the ab Initio Limit for Sulfur−π Interactions: The H2S−Benzene Dimer. J. Phys. Chem. A 2005, 109, 191–196. (47) Cabaleiro-Lago, E. M.; Rodríguez-Otero, J.; Peña-Gallego, Á. Computational Study of the Interaction of Indole-Like Molecules with Water and Hydrogen Sulfide. J. Chem. Phys. 2011, 135, 134310. (48) Biswal, H. S.; Shirhatti, P. R.; Wategaonkar, S. O−H···O versus O−H···S Hydrogen Bonding I: Experimental and Computational Studies on the p-Cresol·H2O and p-Cresol·H2S Complexes. J. Phys. Chem. A 2009, 113, 5633–5643. (49) Biswal, H. S.; Shirhatti, P. R.; Wategaonkar, S. O−H···O versus O−H···S Hydrogen Bonding. 2. Alcohols and Thiols as Hydrogen Bond Acceptors. J. Phys. Chem. A 2010, 114, 6944–6955. (50) Biswal, H. S.; Wategaonkar, S. O−H···O versus O−H···S Hydrogen Bonding. 3. IR−UV Double Resonance Study of Hydrogen Bonded Complexes of p-Cresol with Diethyl Ether and Its Sulfur Analog. J. Phys. Chem. A 2010, 114, 5947–5957. (51) Snow, K. B.; Thomas, T. F. Mass Spectrum, Ionization Potential, and Appearance Potentials for Fragment Ions of Sulfuric Acid Vapor. Int. J. Mass Spectrom. Ion Process. 1990, 96, 49–68. (52) Prest, H. F.; Tzeng, W.-B.; Brom, J. M.; Ng, C. Y. Molecular Beam Photoionization Study of H2S. Int. J. Mass Spectrom. Ion Phys. 1983, 50, 315–329. (53) Nourbakhsh, S.; Norwood, K.; Yin, H.-M.; Liao, C.-L.; Ng, C. Y. Vacuum Ultraviolet Photodissociation and Photoionization Studies of CH3SH and SH. J. Chem. Phys. 1991, 95, 946–954. (54) Carnovale, F.; Livett, M. K.; Peel, J. B. Indentification of the Gas-Phase Trimer (CH3)2S.(HF)2 by Photoelectron Spectroscopy. J. Am. Chem. Soc. 1983, 105, 6788–6790. (55) Chewter, L. A.; Sander, M.; MüllerDethlefs, K.; Schlag, E. W. High Resolution Zero Kinetic Energy Photoelectron Spectroscopy of Benzene and Determination of the Ionization Potential. J. Chem. Phys. 1987, 86, 4737–4744. (56) Lu, K. T.; Eiden, G. C.; Weisshaar, J. C. Toluene Cation: Nearly Free Rotation of the Methyl Group. J. Phys. Chem. 1992, 96, 9742–9748. (57) Hager, J. W.; Wallace, S. C. Two-Laser Photoionization Supersonic Jet Mass Spectrometry of Aromatic Molecules. Anal. Chem. 1988, 60, 5–10. (58) Lipert, R. J.; Colson, S. D. Accurate Ionization Potentials of Phenol and P-(H2O) from the Electric Field Dependence of the Pump–Probe Photoionization Threshold. J. Chem. Phys. 1990, 92, 3240–3241.

(59) Johnstone, R. A. W.; Mellon, F. A. Effects of Induction and Resonance in the Calculation of Ionization Potentials of Substituted Benzenes by Perturbation Molecular Orbital Theory. J. Chem. Soc. Faraday Trans. 2 1973, 69, 36–42. (60) Main-Bobo, J.; Olesik, S.; Gase, W.; Baer, T.; Mommers, A. A.; Holmes, J. L. The Thermochemistry and Dissociation Dynamics of Internal-Energy-Selected Pyrazole and Imidazole Ions. J. Am. Chem. Soc. 1986, 108, 677–683. (61) Close, D. M. Calculated Vertical Ionization Energies of the Common α-Amino Acids in the Gas Phase and in Solution. J. Phys. Chem. A 2011, 115, 2900–2912. (62) Cannington, P. H.; Ham, N. S. He(I) and He(II) Photoelectron Spectra of Glycine and Related Molecules. J. Electron Spectros. Relat. Phenomena 1983, 32, 139–151. (63) Cheney, B. V.; Schulz, M. W.; Cheney, J. Complexes of Benzene with Formamide and Methanethiol as Models for Interactions of Protein Substructures. Biochim. Biophys. Acta Protein Struct. Mol. Enzymol. 1989, 996, 116–124. (64) Yan, S.; Lee, S.-J.; Kang, S.-W.; Choi, K.-H.; Rhee, S.-K.; Lee, J.-Y. Attractive Sulfur···π Interaction between Fluorinated Dimethyl Sulfur (FDMS) and Benzene. Bull. Korean Chem. Soc. 2007, 28, 959–964. (65) Wang, Y.; Paulus, B. A Comparative Electron Correlation Treatment in H2S-Benzene Dimer with DFT and Wavefunction-Based ab Initio Methods. Chem. Phys. Lett. 2007, 441, 187–193. (66) Cabaleiro-Lago, E. M.; Carrazana-García, J. A.; Rodríguez-Otero, J. Study of the Interaction between Water and Hydrogen Sulfide with Polycyclic Aromatic Hydrocarbons. J. Chem. Phys. 2009, 130, 234307. (67) Sherrill, C. D.; Takatani, T.; Hohenstein, E. G. An Assessment of Theoretical Methods for Nonbonded Interactions: Comparison to Complete Basis Set Limit Coupled-Cluster Potential Energy Curves for the Benzene Dimer, the Methane Dimer, Benzene−Methane, and Benzene−H2S. J. Phys. Chem. A 2009, 113, 10146–10159. (68) Biswal, H. S.; Wategaonkar, S. Sulfur, Not Too Far Behind O, N, and C: SH···π Hydrogen Bond. J. Phys. Chem. A 2009, 113, 12774–12782. (69) Mintz, B. J.; Parks, J. M. Benchmark Interaction Energies for Biologically Relevant Noncovalent Complexes Containing Divalent Sulfur. J. Phys. Chem. A 2012, 116, 1086–1092. (70) Zhou, F.; Liu, R.; Li, P.; Zhang, H. On the Properties of S⋯O and S⋯π Noncovalent Interactions: The Analysis of Geometry, Interaction Energy and Electron Density. New J. Chem. 2015, 39, 1611–1618. (71) Senćanski, M.; Došen-Mićović, L.; Šukalović, A. V.; Kostić-Rajačić, S. Theoretical Insight into Sulfur–Aromatic Interactions with Extension to D2 Receptor Activation Mechanism. Struct. Chem. 2015, 26, 1139–1149. (72) Mazmanian, K.; Sargsyan, K.; Grauffel, C.; Dudev, T.; Lim, C. Preferred Hydrogen-Bonding Partners of Cysteine: Implications for Regulating Cys Functions. J. Phys. Chem. B 2016, 120, 10288–10296. (73) Zhou, P.; Tian, F.; Lv, F.; Shang, Z. Geometric Characteristics of Hydrogen Bonds Involving Sulfur Atoms in Proteins. Proteins Struct. Funct. Bioinforma. 2009, 76, 151–163. (74) Chen, X.; Tao, Y.; Li, J.; Dai, H.; Sun, W.; Huang, X.; Wei, Z. Aromatic Residues Regulating Electron Relay Ability of S-Containing Amino Acids by Formations of S∴π Multicenter Three-Electron Bonds in Proteins. J. Phys. Chem. C 2012, 116 (37), 19682–19688.

(75) Miller, M. A.; Vitello, L.; Erman, J. E. Regulation of Interprotein Electron Transfer By Trp 191 of Cytochrome c Peroxidase. Biochemistry 1995, 34, 12048–12058. (76) Yuan, T.; Weljie, A. M.; Vogel, H. J. Tryptophan Fluorescence Quenching by Methionine and Selenomethionine

11



Residues of Calmodulin: Orientation of Peptide and Protein Binding. Biochemistry 1998, 37, 3187–3195. (77) Imai, Y. N.; Inoue, Y.; Yamamoto, Y. Propensities of Polar and Aromatic Amino Acids in Noncanonical Interactions: Nonbonded Contacts Analysis of Protein-Ligand Complexes in Crystal Structures. J. Med. Chem. 2007, 50, 1189–1196. (78) Gaßel, M.; Breitenlechner, C. B.; Rüger, P.; Jucknischke, U.; Schneider, T.; Huber, R.; Bossemeyer, D.; Engh, R. A. Mutants of Protein Kinase A that Mimic the ATP-binding Site of Protein Kinase B (AKT). J. Mol. Biol. 2003, 329, 1021–1034. (79) Meyer, E. A.; Furler, M.; Diederich, F.; Brenk, R.; Klebe, G. Synthesis and In Vitro Evaluation of 2-Aminoquinazolin-4(3H)One-Based Inhibitors for tRNA-Guanine Transglycosylase (TGT). Helv. Chim. Acta 2004, 87, 1333–1356. (80) Kawatkar, S. P.; Kuntz, D. A.; Woods, R. J.; Rose, D. R.; Boons, G.-J. Structural Basis of the Inhibition of Golgi αMannosidase II by Mannostatin A and the Role of the Thiomethyl Moiety in Ligand-Protein Interactions. J. Am. Chem. Soc. 2006, 128, 8310–8319. (81) Wogulis, M.; Wheelock, C. E.; Kamita, S. G.; Hinton, A. C.; Whetstone, P. A.; Hammock, B. D.; Wilson, D. K. Structural Studies of a Potent Insect Maturation Inhibitor Bound to the Juvenile Hormone Esterase of Manduca Sexta. Biochemistry 2006, 45, 4045– 4057. (82) Mileni, M.; Garfunkle, J.; DeMartino, J. K.; Cravatt, B. F.; Boger, D. L.; Stevens, R. C. Binding and Inactivation Mechanism of a Humanized Fatty Acid Amide Hydrolase by α-Ketoheterocycle Inhibitors Revealed from Cocrystal Structures. J. Am. Chem. Soc. 2009, 131, 10497–10506. (83) Ritschel, T.; Kohler, P.; Neudert, G.; Heine, A.; Diederich, F.; Klebe, G. How to Replace the Residual Solvation Shell of Polar Active Site Residues to Achieve Nanomolar Inhibition of tRNAGuanine Transglycosylase. ChemMedChem 2009, 4, 2012–2023. (84) Salonen, L. M.; Ellermann, M.; Diederich, F. Aromatic Rings in Chemical and Biological Recognition: Energetics and Structures. Angew. Chemie Int. Ed. 2011, 50, 4808–4842. (85) Lewis, A. K.; Dunleavy, K. M.; Senkow, T. L.; Her, C.; Horn, B. T.; Jersett, M. A.; Mahling, R.; McCarthy, M. R.; Perell, G. T.; Valley, C. C. et al. Oxidation Increases the Strength of the Methionine-Aromatic Interaction. Nat. Chem. Biol. 2016, 12, 860–866.

12


Page 12 of 13

Page 13 of 13 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


TOC Graphic

13


Modeling Protein SâAromatic Motifs Reveals Their Structural and

Recommend Documents

Modeling Protein SâAromatic Motifs Reveals Their Structural and