The Role of the Hydroxyl Group in Propofol–Protein Target

May 26, 2017 - Propofol (PFL, 1-hydroxyl-2,6-diisopropylbenzene) is currently used widely as one of the most well-known intravenous anesthetics to rel...
4 downloads 11 Views 2MB Size
Article pubs.acs.org/JPCB

The Role of the Hydroxyl Group in Propofol−Protein Target Recognition: Insights from ONIOM Studies Ling Qiu,†,‡ Jianguo Lin,† Qingzhu Liu,† Shanshan Wang,† Gaochao Lv,† Ke Li,† Haiming Shi,† Zhengkun Huang,† and Edward J. Bertaccini*,‡,§ †

Key Laboratory of Nuclear Medicine, Ministry of Health, & Jiangsu Key Laboratory of Molecular Nuclear Medicine, Jiangsu Institute of Nuclear Medicine, Wuxi 214063, P. R. China ‡ Department of Anesthesia, Stanford University School of Medicine, 300 Pasteur Drive, Stanford, California 94305, United States § Palo Alto VA Health Care System, 112A, PAVAHCS, 3801 Miranda Avenue, Palo Alto, California 94304, United States S Supporting Information *

ABSTRACT: Propofol (PFL, 1-hydroxyl-2,6-diisopropylbenzene) is currently used widely as one of the most well-known intravenous anesthetics to relieve surgical suffering, but its mechanism of action is not yet clear. Previous experimental studies have demonstrated that the hydroxyl group of PFL plays a dominant role in the molecular recognition of PFL with receptors that lead to hypnosis. To further explore the mechanism of anesthesia induced by PFL in the present work, the exact binding features and interaction details of PFL with three important proteins, human serum albumin (HSA), the pH-gated ion channel from Gloeobacter violaceus (GLIC), and horse spleen apoferritin (HSAF), were investigated systematically by using a rigorous three-layer ONIOM (M06-2X/6-31+G*:PM6:AMBER) method. Additionally, to further characterize the possible importance of such hydroxyl interactions, a similar set of calculations was carried out on the anesthetically inactive fropofol (FFL, 1-fluoro-2,6-diisopropylbenzene) in which the fluorine was substituted for the hydroxyl. According to the ONIOM calculations, atoms in molecules (AIM) analyses, and electrostatic potential (ESP) analyses, the significance of hydrogen bond, halogen bond, and hydrophobic interactions in promoting proper molecular recognition was revealed. The binding interaction energies of PFL with different proteins were generally larger than FFL and are a significant determinant of their differential anesthetic efficacies. Interestingly, although the hydrogen-bonding effect of the hydroxyl moiety was prominent in propofol, the substitution of the 1-hydroxyl by a fluorine atom did not prevent FFL from binding to the protein via a halogen-bonding interaction. It therefore became clear that multiple specific interactions rather than just hydrogen or halogen bonds must be taken into account to explain the different anesthesia endpoints caused by PFL and FFL. The contributions of key residues in ligand−receptor binding were also quantified, and the calculated results agreed with many available experimental observations. This work will provide complementary insights into the molecular mechanisms of anesthetic action for PFL from a robust theoretical point of view. This will not only assist in interpreting experimental observations but will also help to develop working hypotheses for further experiments and future drug design.



INTRODUCTION

to the Cys-loop superfamily of pentameric ligand-gated ion channels (pLGICs) that also include glycine (GlyR), serotonin (5-HT3), and nicotinic acetylcholine (nAChR) receptors.7−13 However, so far only a few proteins have been proven to have specific interactions with propofol.13−15 On the other hand, quantitative structure activity relationships (QSAR) have also been examined for propofol and its analogues with various substituents on the alkyl phenol backbone. These studies reveal that the formation of an intermolecular hydrogen bond involving the proton of the propofol hydroxyl group is a very

General anesthetics have been employed for nearly two centuries to relieve the suffering due to surgery and other invasive procedures. This has been acknowledged as a cornerstone of modern medicine. In spite of their long and widespread applications in the clinical setting, the mechanism underlying how general anesthetics exert their effects still remains unclear and controversial.1−4 Propofol is known as one of the most commonly used intravenous general anesthetics for induction and maintenance of anesthesia.5 Extensive studies have been performed to explore its mechanism of action.6−11 On one hand, many studies have suggested that membrane proteins might be the most probable target of PFL, such as the GABAA (γ-aminobutyric acid type A) receptors, which belong © 2017 American Chemical Society

Received: March 3, 2017 Revised: May 20, 2017 Published: May 26, 2017 5883

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

binding features and interaction details of propofol and fropofol within a set of well described protein binding sites were studied and compared for the first time using the state-of-the-art ONIOM (Our own N-layered Integrated molecular Orbital and molecular Mechanics) methodology.27 This method has become an important technique for studying very large and complex biological systems by allowing the division of an entire system into several regions, each of which can be studied with a different theoretical level of computational chemistry.27 Hence, it is feasible to implement high-resolution quantum mechanics computations for the most interesting portion while not requiring it for the entire system thereby providing great detail for relevant interaction regions while optimizing computational efficiencies. We have also successfully employed this method to explore the mechanism of action of several inhalational and intravenous anesthetics.26 Based on the ONIOM calculations presented in this manuscript in which we compare the binding features and interaction energies of PFL and FFL, we identify the contributions of key residues for propofol binding within different protein binding sites as well as clarify the role of the 1hydroxyl in propofol−protein target recognition within different proteins. Furthermore, a comparison of the present theoretical modeling results with the previously reported experimental data18 sheds new light on the mechanism of action of propofol. Therefore, the present study provides complementary insights into propofol’s detailed molecular mechanism of action. This will not only assist in the interpretation of experimental observations but also help the development of future experiments and drug design.

important binding interaction that contributes to the specific molecular recognition by target proteins, thereby allowing ligand efficacy.16,17 Recently, to evaluate the role of the 1hydroxyl group in determining the pharmacologically relevant molecular recognition and anesthetic efficacy of propofol, Woll and co-workers designed and synthesized a propofol analogue in which the 1-hydroxyl group is substituted by a fluorine atom to make 1-fluoro-2,6-diisopropylbenzene. This was referred to as fropofol (FFL) to emphasize its similar physicochemical properties to those of propofol (Scheme 1) though it is Scheme 1. Chemical Structures of Propofol (PFL) and Fropofol (FFL)

completely inactive in regard to the production of anesthesia.18 Through a series of in vitro and in vivo experimental assays, they concluded that hydrogen bonding was necessary for the molecular recognition of PFL by the proteins that were necessary for hypnotic activity. However, the interaction details between the ligand and the protein within the binding sites were not studied and compared at an atomic level. Since identifying the binding interactions of general anesthetics with three-dimensional structures of target proteins is of vital importance toward a greater understanding of the mechanisms of action as well as for further drug development, X-ray crystal structure analysis has attracted considerable attention in order to provide atomistic information about different ligand− receptor interactions.12−15 However, despite the increasing number of crystal structures deposited in the Protein Data Bank (PDB), few structures relevant to anesthetic−protein interactions are available because of the difficulties in obtaining proper crystal structures of such complexes. Therefore, to make use of what is available and expand such analyses, it is desirable to employ more efficient computational chemistry tools to explore and identify interactions between the anesthetics and the proteins. With the rapid development of computer technology and theoretical methods, theoretical modeling has become an effective way to explore the mechanism of action of general anesthetics as a complement to experimental analyses.12,19−26 For example, the binding sites and interactions of various anesthetics within proteins have been determined using a variety of theoretical approximations incorporated into in silico modeling tools.12,19−22 The binding sites and common binding features within some well-characterized anesthetic−protein complexes have also been studied by our group using theoretical modeling methods.23−26 The results from such works have shown that important anesthetic−protein interactions often involved several common hydrophobic and hydrophilic types of interactions within an amphiphilic binding cavity, including hydrogen bonds, halogen bonds, van der Waals forces, and other weak polar interactions that many previous studies have implicated as necessary for reversible anesthetic binding. For the present work, in order to further distinguish the reason for propofol potency over that of fropofol, the exact



THEORETICAL METHODS Setup of Models. The X-ray structures of propofol bound to human serum albumin (HSA, PDB ID 1E7A),14 the pHgated ion channel GLIC (pH-gated bacterial homologue from Gloeobacter violaceus, PDB ID 3P50),13 and horse spleen apoferritin (HSAF, PDB ID 3F33)15 were obtained from the RCSB Protein Data Bank and employed as the starting points for the present study. The corresponding starting structures for fropofol bound to GLIC, HSA, and HSAF were obtained by substituting the hydroxyl group in propofol within 3P50, 1E7A, and 3F33, which were referred to as 1E7A′, 3P50′, and 3F33′ to emphasize their similar structural features. All these structure files underwent a similar preparatory sequence prior to theoretical modeling: side chain corrections, adding hydrogen atoms, setting molecular mechanics (MM) parameters, optimizing hydrogen geometries, and defining subsets for subsequent ONIOM calculations. Hydrogen atoms were added and optimized using the Discovery Studio 4.1 software suite (Accelrys Inc., San Diego, CA), and conventional protonation states for all amino acids at pH 7.4 were adopted. Since the proteins HSA, GLIC, and HSAF are organized as a dimer, pentamer, and 24-mer, respectively, there are two, five, and 24 binding sites in which the ligand can be located in each. To improve the efficiency of theoretical calculations, therefore, the active sites located in one subunit of HSA (1E7A, 1E7A′) and GLIC (3P50, 3P50′) or at the interface between subunits 1 and 22 of HASF (3F33, 3F33′) were considered as the starting models. Moreover, these models were further truncated to minimize any unnecessary computational burdens according to the partition scheme noted below. To maintain the main scaffold of surrounding regions of the active site as they were observed in the X-ray structures, all residues with at least one atom interacting with any atom of the ligand within the 5884

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B Table 1. Summary of the ONIOM Models Used in the Present Studya PDB ID

protein

ligand

ligand ID

ONIOM-1

ONIOM-2

RMSD-1/Å

RMSD-2/Å

1E7A 1E7A′ 3P50 3P50′ 3F33 3F33′

HSA HSA GLIC GLIC HSAF HSAF

propofol fropofol propofol fropofol propofol fropofol

PFL4001 FFL4001 PFL319 FFL319 PFL2001 FFL2001

31:339:1197 30:338:1196 31:359:1247 30:358:1246 31:285:1441 30:284:1440

132:339:1197 131:338:1196 90:359:1247 89:358:1246 98:285:1441 97:284:1440

0.7926 0.7596 0.7559 0.8845 1.6678 1.2613

0.7877 0.7510 0.7718 0.8142 2.3792 2.3811

a

ONIOM-1 and ONIOM-2 models with the number of atoms studied by the high level of theory, the medium level of theory, and the low level of theory, respectively. RMSD-1 and RMSD-2 denote the root-mean-square deviation of optimized structures within ONIOM-1 and ONIOM-2 models relative to the original structures, respectively.

Figure 1. Structures of propofol bound with HSA (left), subunit A of HSA with propofol bound to an inner 10 Å protein cavity shown in sticks (middle), and three-layer ONIOM-1 model for propofol bound to an inner 10 Å protein cavity with the high-level layer illustrated in ball-and-stick, the medium-level layer in yellow line, and the low-level layer in gray line (right).

system is a part cut from the real system, and the model system is a part cut from the middle system. The hybrid meta DFT method M06-2X/6-31+G* was employed to describe the model system (the default convergence criteria of Gaussian0929 were used throughout), since its capability of providing accurate descriptions for both hydrogen bond and dispersive interactions has been verified.31 A diffuse function was included in order to treat the lone electron pairs properly and has been demonstrated as a requirement for adequately describing hydrogen-bonded systems. The PM632 semiempirical QM method and AMBER33 molecular mechanics were used to treat the middle system and the real system, respectively. Accordingly, the three-layer ONIOM method notation (M062X/6-31+G*:PM6:AMBER) describes the system in the present work. Nonstandard residues in the whole system were manually assigned GAFF34 and AMBER (Parm99)33 parameters. Atomic charges were also assigned manually to each atom in the nonstandard residues and determined by a restrained fit to the electrostatic potential (ESP) via the Merz− Kollman procedure,35 which was obtained from the HF/631G* calculations as implemented in the Gaussian09 program suite.29 The electronic embedding (EE) scheme28 was used in the ONIOM calculations to include the polarization effect of the protein environment on the QM region. The interactions of propofol and fropofol with different proteins were studied and compared by using two ONIOM schemes, which were referred to as ONIOM-1 and ONIOM-2. In the ONIOM-1 model, the high-level layer only included the atoms of the active ligand, propofol or fropofol, and the medium-level layer included all the atoms from the residues within a 4 Å sphere of the ligand. The latter region sufficiently considers the possible hydrogen- or halogen-bonding inter-

interatomic distance threshold of 10 Å were taken into account. The final models span a total of 1196−1441 atoms (Table 1), including any amino acids and waters that may interact directly with the ligand and are required to account for the effects of different protein environments. ONIOM Calculations. The ONIOM methodology28 implemented in the Gaussian09 program package (Gaussian Inc., Wallingford, CT, USA)29 was chosen to study the interactions of propofol and fropofol bound to different proteins. The three-layer ONIOM method combined three different methods: a high-level quantum mechanics (HQ) region, a low-level quantum mechanics (LQ) region, and the largest and most distant region from the ligand described by molecular mechanics (MM) methods. This subset scheme has been proven to efficiently describe large biological systems in previous works.30 According to the ONIOM methodology, the energy calculations for the system were divided into three overlapping layers, designated as the high-level, medium-level, and low-level layers. The energy for the real system was computed using the follow equation:30 high low ONIOM medium low medium Ereal = Emodel + Emiddle + Ereal − Emodel − Emiddle

(1)

EONIOM real

Ehigh model,

Emedium middle ,

Elow real,

Emedium model ,

Elow middle

where and denote the energy of the real system predicted by the ONIOM method, the energy of the model system predicted at the high level of theory, the energy of the middle system predicted at the medium level of theory, the energy of the real system predicted at the low level of theory, the energy of the model system predicted at the medium level of theory, and the energy of the middle system predicted at the low level of theory, respectively. The entire system is referred to as the real system, the middle 5885

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 2. Structures of propofol bound with GLIC (left), subunits A and B of GLIC with propofol bound to an inner 10 Å protein cavity shown in sticks (middle), and three-layer ONIOM-1 model for propofol bound to an inner 10 Å protein cavity with the high-level layer illustrated in ball-andstick, the medium-level layer in yellow line, and the low-level layer in gray line (right).

Figure 3. Structures of propofol bound with HSAF (left), subunits 1 and 22 of HSAF with propofol bound to the interfacial cavity shown in sticks (middle), and three-layer ONIOM-1 model for propofol bound to an inner 10 Å protein cavity with the high-level layer illustrated in ball-and-stick, the medium-level layer in yellow line, and the low-level layer in gray line (right).

determined via the standard counterpoise (CP) method of Boys and Bernardi37 as implemented in the Gaussian09 program.29 In the ONIOM-2 scheme, the model system was improved by incorporating key residues into the high layer. These key residues were those which displayed strong interactions with the ligand in the ONIOM-1 calculations. These residues were further optimized within the ONIOM-2 partition scheme in order to determine accurate binding interactions by taking into account the relaxation of important environmental residues. The detailed information about the layer divisions and the PDB residues within the ONIOM-1 and ONIOM-2 models was listed in Tables S1 and S2.

actions between the ligand and the surrounding protein. The remainder constitutes the low-level layer (Figures 1−3). During ONIOM-1 optimization, the geometries of propofol and fropofol were fully optimized while the atomic positions of the protein were frozen with hydrogen atoms serving as link atoms to saturate the dangling bonds spanning between two layers. The minimum nature of the stationary point was checked through harmonic vibrational analysis (no imaginary frequency). The direct binding interaction energy (BE) between the ligand and the receptor was calculated according to the following formula:36 BE = Ecomplex − (E ligand + Eprotein)

(2)



where Ecomplex, Eligand, and Eprotein represented the energy of the ligand−protein complex, the ligand, and the protein after ONIOM-1 optimization, respectively. Moreover, to evaluate individual contributions of amino acid residues to the ligand binding to the protein, the pairwise interaction energies (ΔEAB) between the ligand and the individual residues within the binding site of 4 Å were also calculated using the following approach:36 ΔEAB = EAB − EA − E B + BSSE

RESULTS AND DISCUSSION As is well-known, molecular recognition in biological systems is driven by specific intermolecular interactions between a small molecule and its target receptor. The strength of the recognition is usually characterized by the binding affinity, which in turn determines the biological activity of a compound.38 Similar molecular recognition processes are essential for determining many similar biological functions. To investigate and compare the molecular recognition processes of PFL and FFL, available high-resolution crystal structures of propofol bound to the membrane protein GLIC (PDB 3P50)13 and the water-soluble proteins (HSA, PDB

(3)

where EAB was the energy of the ligand−residue pair, EA and EB were the energies of the isolated individual fragments, and BSSE was the correction for the basis set superposition error 5886

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 4. Overlay of experimental (blue) and theoretical (ONIOM-1, red) structures, showing the change in the ligand structure (top, PFL; bottom, FFL) upon relaxation within the binding site.

1E7A;14 and HSAF, PDB 3F3315) were used as the models in the present work. GLIC, HSA, and HSAF are three important proteins, not only because they are amenable to X-ray structure determination but, in particular, GLIC belongs to the family of pLGICs that are considered as the principal targets of general anesthetics;1−4 HSA is one of the most abundant proteins in the circulatory system and plays a key role in determining the distribution and pharmacokinetics of general anesthetics,14 and HSAF has the highest affinity for anesthetics among any protein studied to date.15 Additionally, the affinity of HSAF for a broad range of general anesthetics correlates well with their anesthetic potency.39 ONIOM-1 Model. The binding features of PFL and FFL within three different proteins were first studied and compared by using the ONIOM-1 scheme. The molecular overlay of the optimized structures and the experimental complexes is presented in Figure 4, and the root-mean-square deviation (RMSD) between the optimized structure and the crystal structure was applied as a criterion to evaluate any distortion of the ligand binding modes (Table 1). Nonbonded interactions between PFL or FFL and the proteins HSA, GLIC, and HSAF are shown in Figure 5 and summarized in Table S3. From Table 1, one can see that the RMSD values for the complexes 1E7A, 1E7A′, 3P50, and 3P50′ were all less than 1 Å, indicating that relaxation of the ligand within the protein HSA and GLIC produced little alteration in the structures relative to the crystal structures. However, the RMSD values for both 3F33 and 3F33′ were larger than 1 Å, indicating that there might exist more strong interactions between the ligand and the protein HSAF than expected from the initial coordinates. This will be further analyzed and discussed below. Comparison between PFL and FFL bound to the same protein showed that the difference in RMSDs was small when PFL and FFL interacted similarly with the protein HSA (−0.033 Å), but the differences in RMSDs of PFL and FFL were large when they interacted with different amino acids of the protein GLIC (0.1286 Å) and HSAF (−0.4065 Å). This can be observed from Figure 5 and Table S3. This also indicated that when the OH group was

changed to F, nonbonded interactions between the ligand and the proteins would obviously change to allow a better fit of the ligand to the binding cavity. Many crystal structures deposited in the Protein Data Bank have revealed that ligands are bound to different receptors by means of different kinds of noncovalent interactions, such as hydrogen bonds (HB), halogen bonds (XB), electrostatic interactions, van der Waals forces, hydrophobic interactions, and so on.40 Previous studies have also reported that the pharmacological effect of an individual anesthetic might be dependent on distinctive chemical features giving rise to a relatively higher affinity for a particular protein target.13−18 Hence, understanding the molecular interactions between anesthetics and their protein binding sites is essential for achieving a greater mechanistic understanding of anesthetic action and for the future design of novel agents with desirable effects at specific targets. From Figure 5, it is observed that there are hydrogen bonds or halogen bonds between the ligand, PFL or FFL, and the proteins in addition to several hydrophobic interactions. For instance, in the complex of PFL bound to HSA (1E7A), besides the hydrophobic interactions between PFL and the amino acids LEU387, ILE388, CYS392, PHE403, LEU407, VAL433, CYS438, ALA449, LEU453, and ARG485, there existed one conventional hydrogen bond between the PFL and the amino acid LEU430, which agreed well with experimental observation.14 Here, the “conventional hydrogen bond” refers to that between the classical hydrogen bond donor atoms N, O, P, and S or the hydrogen atoms connected to such atoms and the hydrogen bond acceptor atoms N, O, P, and S. Additionally, atoms of element type F, Cl, Br, and I can also serve as the hydrogen bond acceptors. Similarly, one halogen bond formed between FFL and the amino acid LEU430 besides similar hydrophobic interactions between FFL and the amino acids LEU387, ILE388, CYS392, PHE403, LEU407, VAL433, CYS438, ALA449, LEU453, and ARG48. This clearly showed that the loss of the hydrophilic feature (OH group) did not 5887

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 5. Nonbonded interactions between PFL (right) or FFL (left) and the protein HSA, GLIC, and HSAF within the ONIOM-1 calculation scheme.

interactions between PFL and FFL, the calculated hydrogen bond or halogen bond details between the ligand and the proteins have been summarized in Table 2. The hydrogen bond or halogen bond strengths have also been evaluated by the atoms in molecules (AIM) theory45 using the AIM2000 2.0 program,46 which is a very useful tool to characterize atomic and molecular interactions. The topological and energetic properties for the H-bond critical point are often used as a measure of the H-bonding strength. Here, only the dissociation energies of H-bond or X-bond (EHB, in kJ/mol) are listed in Table 2.

prevent FFL binding to the protein via the formation of a halogen-bonding interaction. It is known that the hydrogen bond is one of the most important classes of intermolecular interactions in biology,38,40,41 which usually provides a directionality and specificity of interaction that is a fundamental aspect of molecular recognition for a ligand within its protein binding site. Recently, the role of halogen bonds in the fields of medicinal chemistry and molecular design have been increasingly recognized as significant.38,40,42−44 In order to analytically compare the hydrogen bond and halogen bond 5888

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B Table 2. Hydrogen/Halogen Bond Details Calculated with the ONIOM-1 and ONIOM-2 Schemesa

a

PDB

HB donor

HB acceptor

dD−H/Å

1E7A 1E7A′ 3F33′

PFL4001:H1 FFL4001:F1 22-ARG59:HD2

LEU430:O LEU430:O 1-FFL2001:F1

0.973 1.358 1.090

1E7A 1E7A′ 3P50 3P50′ 3F33 3F33′

PFL4001:H1 FFL4001:F1 TYR254:HH TYR254:HH 22-ARG59:HD2 22-ARG59:HD2

LEU430:O LEU430:O PFL319:O1 FFL319:F1 PFL2001:O1 FFL2001:F1

0.968 1.351 0.976 0.970 1.093 1.091

dH···A/Å ONIOM-1 2.789 3.602 2.476 ONIOM-2 2.609 2.985 1.899 1.925 3.097 2.878

dD···A/Å

∠D−H···A/deg

∠H···A−Y/deg

EHB/kJ·mol−1

3.537 4.848 3.515

134.3 152.7 159.0

137.6 127.9 157.6

−2.63 −0.92 −6.55

3.025 4.315 2.805 2.828 3.526 3.332

106.1 168.0 153.2 153.8 104.1 105.1

164.4 163.9 100.9 131.7 114.8 120.3

−3.74 −5.81 −33.83 −30.59 −2.13 −11.26

In the complex 1E7A′, dD−H, dH···A, dD···A, ∠D−H···A, and ∠H···A−Y correspond to dC−F, dF···O, dC···O, ∠C−F···O, and ∠F···O−C.

Table 3. Binding Interaction Energies (BE) for PFL/FFL with Different Proteins Computed with the ONIOM-1 and ONIOM-2 Schemesa model

complex

Ecomplex/au

Eligand/au

Eprotein/au

BE/kJ·mol−1

ONIOM-1

1E7A−PFL 1E7A′-FFL 3P50−PFL 3P50′−FFL 3F33−PFL 3F33′−FFL 1E7A−PFL 1E7A′−FFL 3P50−PFL 3P50′−FFL 3F33−PFL 3F33′−FFL

−545.7168928 −569.7252702 −546.1137946 −570.1299264 −358.8727209 −382.8711498 −3683.689255 −3707.697448 −1833.025645 −1857.037348 −1788.617435 −1812.614111

−543.1097202 −567.1275431 −543.1081408 −567.1278908 −543.1101182 −567.1275363 −543.1061233 −567.1249752 −543.1093826 −567.1272508 −543.1045303 −567.1235436

−2.579904465 −2.579904465 −2.981238393 −2.981238393 184.2777936 184.2777936 −3140.542721 −3140.542780 −1289.868682 −1289.869511 −1245.350169 −1245.436386

−71.59 −46.79 −64.10 −54.60 −106.06 −56.20 −106.10 −77.96 −124.92 −106.56 −427.26 −142.25

ONIOM-2

a

Ecomplex, Eligand, and Eprotein denoted the energy of the complex, the ligand, and the protein after ONIOM optimization, respectively.

similar ligands, propofol and fropofol, to the same protein cavity by comparing the calculated BE values because the difference in entropy of propofol and fropofol is little and the temperature is identical. Hence, according to the formula ΔG = ΔH − TΔS, the binding Gibbs free energy difference between two different ligands, Δ(ΔG), is approximately equal to the difference in enthalpy change Δ(ΔH) between the complex and the sum of its constituent parts. In other words, the binding affinity can be approximately reflected by the energy difference between the complex and the sum of its constituent parts [Ecomplex − (Eligand + Eprotein)] as calculated. Accordingly, one can infer that PFL might bind with higher affinity for HSA than FFL since it contains both hydrophilic and hydrophobic features. This was consistent with the previous experimental study in which PFL showed higher affinity than FFL for a protein model containing hydrogen-bonding interactions.18 The ONIOM-1 calculations also showed that there were no hydrogen-bonding interactions in the complex of PFL bound to GLIC (3P50), but did demonstrate hydrophobic interactions between PFL and the amino acids PRO120, ILE202, MET205, LEU206, ILE258, and ILE259. Similarly, there was no halogenbonding interaction between FFL and GLIC; yet there were several hydrophobic interactions between FFL and the amino acids PRO120, PHE121, ILE202, MET205, LEU206, ILE258, and ILE259. The binding interaction energy between PFL and GLIC was calculated to be −64.10 kJ·mol−1, while that between FFL and GLIC was −54.60 kJ·mol−1, suggesting that PFL might bind to GLIC with higher affinity than FFL. However, this differed from the previous study in which Woll et al.

According to the Jeffrey classification of hydrogen bonds (HB),41 the intermolecular distances between the proton donor and the proton acceptor are as follows: dH···A within 1.2−1.5 Å or dD···A within 2.2−2.5 Å characterizes a strong HB; dH···A within 1.5−2.2 Å or dD···A within 2.5−3.2 Å characterizes a moderate HB, and dH···A longer than 2.2 Å or dD···A longer than 3.2 Å is classified as a weak HB. Therefore, it can be inferred that the hydrogen-bonding interaction between PFL and HSA is weak. This can also be deduced from the calculated bond dissociation energy (−2.63 kJ·mol−1). Comparison of the bond dissociation energy of HB in PFL and XB in FFL showed that the hydrogen-bonding interaction between PFL and HSA was stronger than the halogen-bonding interaction between FFL and HSA. Since the hydrophobic interactions between PFL and FFL with HSA are approximately identical, it is therefore inferred that the binding affinity of PFL for HSA might be larger than that of FFL, but as is discussed below the interactions are even more complicated than that. Moreover, the binding interaction energies (BE) between the ligands and the proteins calculated according to the previous method36 also showed the same trend (Table 3). The BE for PFL with HSA (−71.59 kJ·mol−1) is obviously larger than that for FFL with HSA (−46.79 kJ·mol−1). Although the calculated BE values cannot directly correspond to the binding constant Kd that can be derived from the “binding Gibbs free energy (ΔG)” (ΔG = −RT ln Kd), they still can reflect the binding affinity of the ligand for the receptors to a certain extent especially when comparing two ligands in the same binding site. It also seems feasible to compare the binding affinity of such 5889

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 6. Individual contributions of each amino acid residue within the 4 Å protein cavity for the binding of PFL or FFL to the proteins HSA, GLIC, and HSAF, respectively. The pairwise interaction energies were calculated using the ONIOM-1 and ONIOM-2 schemes with M06-2X/631+G* and M06-2X/6-311++G** methods, respectively. For clarity, the basis sets 6-31+G* and 6-311++G** were denoted by BS1 and BS2, respectively.

reported that FFL had higher affinity than PFL for a protein site in the absence of hydrogen bond interactions.18 Such a discrepancy between the computational and experimental results might stem in large part from the underestimation of the ligand hydration effects in protein−ligand binding, which would lead to overestimation of the biochemical potency of a ligand with many hydrophilic moieties.47 In the complex of PFL bound to HSAF (3F33), only hydrophobic interactions existed between PFL and the amino acids LEU24, TYR28, LEU31, and ARG59, consistent with the experimental observation that the propofol found a highly hydrophobic pocket in the HSAF binding pocket.15 The closest distance between PFL and the amino acid ARG59 of subunit 22 was 3.353 Å (PFL2001:O1···22-ARG59:HA), which met the criterion for a hydrogen bond and agreed well with the experimental observation that hydrogen bonds were not apparent in the crystal structure of PFL with HSAF.15 However, one carbon hydrogen bond formed between FFL and the amino acid ARG59 of subunit 22 in the complex 3F33′ (FFL2001:F1···22-ARG59:HD2:2.476 Å) in addition to hydrophobic interactions between FFL and the amino acids LEU24,

TYR28, LEU31, ALA55, and ARG59. The carbon−hydrogen bond interactions are weaker hydrogen bonds where the C−H group acts as the hydrogen bond donor.48 There will be weak polarization of the C−H bond in the presence of a suitable acceptor atom, and the carbon atom can serve as a donor if it is either in an acetylene group or adjacent to an oxygen or nitrogen atom. This further demonstrated that the substitution of the OH group by F did not prevent the ligand FFL from binding with the protein HSAF via the formation of a hydrogen bond. However, detailed inspection of protein structures suggested that although this type of interaction was of some significance, a C−H bond will have only minor contributions to the overall stability of a protein structure.48 According to the calculated binding interaction energies, one could infer that PFL can bind with almost 2-fold higher affinity than FFL (−106.06 kJ·mol−1 vs −56.20 kJ·mol−1) to the protein HSAF. This might be the main reason for explaining the different anesthetic efficacy induced by PFL and FFL to some extent. In short, ONIOM-1 calculations suggested that the substitution of 1-hydroxyl by the fluorine atom would lead to weaker binding 5890

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 7. Overlay of experimental (blue) and theoretical (ONIOM-2, red) structures, showing the structural changes of the ligand (top, PFL; bottom, FFL) and key residues involved in strong nonbonded interactions upon relaxation within the binding site.

kJ/mol. In essence, the agreement between M06-2X/6-31+G* and M06-2X/6-311++G** calculated results was prominent and the contribution of each amino acid residue to ligand− protein binding was independent of the basis set. Therefore, it can be concluded that the method M06-2X/6-31+G* is suitable for the description of PFL− or FFL−protein complexes. The pairwise interaction energies were further compared to determine the relative contributions of various amino acid residues in the global PFL/FFL−protein binding process. According to the calculated results in Figure 6, the amino acids ASN391, CYS392, ARG410, CYS438, and ARG485 in HSA are important contributors to the total binding interactions of PFL and FFL to HSA (1E7A and 1E7A′), with respective interaction energies of −13.21 (−13.19), −21.29 (−19.94), −21.06 (−15.62), −20.18 (−15.81), and −12.48 (−6.34) kJ/mol. As for the binding of PFL or FFL to the protein GLIC (3P50 and 3P50′), the amino acids ILE202, TYR254, and ILE258 are the most important contributors to the total interactions, with respective interaction energies of −10.82 (−14.50), −12.77 (−4.60), and −3.98 (−14.97) kJ/mol. With respect to the binding of PFL to the protein HSAF (3F33), ARG59 in the subunit 22 was the most important contributor to the total binding interaction although it was not hydrogen-bonded to the ligand PFL after ONIOM-1 optimization. The pairwise interaction energy calculated at the M06-2X/6-311+ +G**level was −30.2 kJ/mol, which was 4-fold larger than the total interaction energy (−7.3 kJ/mol). However, the amino acids LEU24 and ARG59 in subunit 1 also contributed more to the total binding interactions of FFL with HSAF in addition to ARG59 in subunit 22. The pairwise interaction energies were −13.16, −14.57, and −21.90 kJ/mol, respectively, which accounted for 31%, 34%, and 51% of the total interaction energies, respectively. ONIOM-2 Model. To improve both the simulated geometries and interaction energy calculations while attempting to minimize large increases in computational time, key residues involved in strong interactions with the ligands derived from the ONIOM-1 calculations were also incorporated into the high layer of the ONIOM-2 models. In addition, since the crystal structures and ONIOM-1 calculations both indicated that the

interactions to the specific propofol binding sites, despite the presence of hydrogen bond or halogen bond interactions. A purely geometric analysis of such weak interactions was not sufficient to quantify or discriminate the contribution of each single residue to the ligand−receptor molecular recognition or binding process. Therefore, further calculations were necessary for revealing the molecular recognition of the ligands PFL and FFL in different proteins. The pairwise interaction energies involved in the PFL− or FFL−protein interactions were calculated at the M06-2X/6-31+G* level to show the individual contributions of each amino acid to the binding of PFL or FFL within each receptor. For clarity and comparison, the calculations only took into account the contributions of the amino acid residues within a 4 Å sphere of the ligand binding cavity. To check the reliability of the basis set (6-31+G*) used, the pairwise interaction energies were also calculated with the 6-311++G** basis set. The calculated results have been summarized in Figure 6, which clearly shows the trends in pairwise interaction energies as a function of method and basis set, and displays which residues play more important roles in ligand−receptor binding. In general, the negative interaction energies indicate attractive interactions for the binding of the ligand to the protein, while the positive values indicate repulsive interactions. As shown in Figure 6, M06-2X/6-31+G* predicted negative interaction energies for most of the residues, indicating that they were attractive for the binding of the ligand to the protein. As a whole, the differences among the pairwise interaction energies caused by two different basis sets were small (less than 2 kJ/mol) although the basis set 6-311++G** produced a slightly more negative interaction energy than 6-31+G*. The discrepancy was slightly larger for the ligand−amino acid pair involved in an important hydrogen bond or π-stacking interaction. For instance, the interaction energy between PFL and the amino acid ARG59 in the complex 3F33 computed at the M06-2X/6-31+G* level was −28.6 kJ/mol, while that at the M06-2X/6-311++G** level was ca. −30.2 kJ/mol. Similarly, the interaction energy between FFL and ARG59 in the complex 3F33′ computed at the M06-2X/6-31+G* level was −20.9 kJ/ mol, while that at the M06-2X/6-311++G** level was −21.9 5891

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 8. Overlay of PFL (top) and FFL (bottom) in different receptors HSA, GLIC, and HSAF: experimental (blue), ONIOM-1 (red), and ONIOM-2 (gray) geometries.

the complex 3P50 in which the hydroxyl group of propofol could form a hydrogen bond with the amino acid residue TYR254.13 Compared with the ONIOM-1 calculations, ONIOM-2 calculations produced some new nonbonded interactions. This was mainly due to the fact that only the flexibility of the ligand PFL or FFL was considered in the ONIOM-1 calculations while the flexibility of key residues involved in strong interactions with the ligand were also considered in the ONIOM-2 calculations. Therefore, nonbonded interactions between the ligand and the proteins could change significantly to allow a better fit of the ligand to the binding cavity. The hydrogen bond and halogen bond details are also listed in Table 2 for comparison with ONIOM-1 calculations. In general, it is thought that the main determinant of the strength of a hydrogen bond or halogen bond is its length. The shorter the hydrogen or halogen bond distance, the stronger the hydrogen bond or halogen bond strength. From Table 2, one can find that the hydrogen bonds and halogen bonds within the ONIOM-2 models are stronger than those within the ONIOM1 models. The strength of hydrogen-bonding or halogenbonding interactions (EHB) really depended on the H/X···A or D···A distance. This also further revealed that the ligand and key residues in the ONIOM-2 calculations would reorient dramatically to optimize various favorable interactions between the ligand and nearby functional groups. It is acknowledged that the binding interaction energy usually depends on the size of the binding pocket considered in the calculations. By increasing the radius of the binding pocket, more residues will be included which may contribute more to the binding process. As anticipated, the ONIOM-2 calculations produced larger binding interaction energies (BE) between the ligand and the receptor than the ONIOM-1 calculations (Table 3). Moreover, one can observe from Table 3 that the calculated binding interaction energies between the ligand PFL and the proteins in ONIOM-2 models were all larger than those between the ligand FFL and the proteins, further implying that PFL might bind with higher affinity to the proteins HSA, GLIC, and HSAF than FFL. This further confirmed the experimental study in which PFL showed higher affinity than FFL for a protein model containing hydrogen bond interactions.18 Accordingly, it appears that the binding interaction energies could explain the differences between PFL anesthetic activity and the lack thereof in FFL.

amino acid LEU430 was hydrogen-bonded to the ligand PFL in 1E7A and halogen-bonded to the ligand FFL in 1E7A′, respectively, it was also included in the high layer of ONIOM-2 models. Detailed information about the layer division and PDB residues within the ONIOM-2 models is listed in Table S2. Overlay of experimental and ONIOM-2 optimized structures is shown in Figure 7. These calculations produced the structural changes of the ligand PFL or FFL and key residues involved in strong nonbonded interactions upon relaxation within different protein binding sites. Compared with the experimental structures, the most significant changes in the ONIOM-2 optimizations occurred for the structural optimizations in PFL and FFL bound to HSAF (3F33 and 3F33′), respectively. This could be observed clearly from Figure 8, which shows the overlay of experimental, ONIOM-1, and ONIOM-2 optimized structures of PFL and FFL in the different receptors HSA, GLIC, and HSAF. The RMSD values for the complexes 1E7A, 1E7A′, 3P50, and 3P50′ were all less than 1 Å (Table 1), indicating that the optimized structures deviated only slightly relative to the experimental structures, while those for 3F33 and 3F33′ were larger than 2 Å, indicating that the optimized structures changed significantly in comparison with the experimental structures. This was mainly due to the fact that the flexibilities of key residues involved in strong interactions with the ligand were considered in the ONIOM-2 calculations, while only the flexibility of the ligand PFL or FFL was considered in the ONIOM-1 calculations. Hence, the ligand and key residues would adjust their positions significantly to allow for better molecular recognition and binding within the protein environment. Taken together, the protein HSAF has a larger influence on the recognition of PFL and FFL than the proteins HSA and GLIC in spite of ONIOM-1 or ONIOM-2 calculations. On the other hand, one can find that the discrepancies in RMSDs of PFL and FFL within the same protein were very small, which might be attributed to the fact that ONIOM-2 optimizations produced similar nonbonded interactions between PFL and FFL within the same protein (Figure 9, Table S4). Besides forming similar hydrophobic interactions with the same residues in each protein, PFL and FFL were all hydrogen-bonded to the amino acid TYR254 in GLIC and ARG59 in HSAF, respectively, except that PFL was hydrogenbonded to the amino acid LEU430 in HSA while FFL was halogen-bonded to LEU430 in HSA, respectively. This supported the previous experimental hypothesis concerning 5892

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 9. Nonbonded interactions between PFL (right) or FFL (left) and the proteins HSA, GLIC, and HSAF within the ONIOM-2 calculation scheme.

Based on the ONIOM-2 optimized structures, the pairwise interaction energies between the ligand PFL or FFL and the residues within a 4 Å sphere of the protein cavity were also calculated at the level of M06-2X/6-31+G* and M06-2X/6311++G**, respectively. These calculations revealed the individual contributions of the amino acids to the ligand− receptor binding processes. The calculated results have also been summarized in Figure 6 for better comparison with those obtained from the ONIOM-1 calculations. From Figure 6, it

was clearly observed that the basis set had little effect on the pairwise interaction energies, but the ONIOM model had relatively large effects on the pairwise interaction energies. On the whole, the relaxation of key residues produced more negative pairwise interaction energies between the ligand and the residues in the ONIOM-2 models than those in the ONIOM-1 models. This was again mainly due to the fact that only the flexibility of the ligand PFL or FFL was considered in the ONIOM-1 calculations while the flexibility of key residues 5893

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B

Figure 10. ESP of PFL (top) and FFL (bottom) within different proteins (HSA, HSAF, and GLIC) calculated at the M06-2X/6-31+G* level and mapped onto the isosurface of the electron density (0.001 electron per Å3). Blue and red regions represent positive and negative potential areas, respectively.

negativity. A halogen atom serving as the hydrogen bond acceptor is a typical example of such interactions. As shown in Figure 9 and Table 2, ONIOM-2 optimizations showed the formation of carbon−hydrogen bonding between the fluorine atom of FFL and the hydrogen atom of the amino acid TYR254 in the protein GLIC, as well as between the fluorine atom of FFL and the hydrogen atom of the amino acid ARG59 in the protein HSAF. A halogen atom can also function as a halogen bond donor according to the unique IUPAC definition of the halogen bond (R−X···Y),51 where X is any halogen atom with an electronpoor region, R is a group covalently bound to X, and Y is typically a molecular entity possessing at least one electron-rich region. This was demonstrated here by the existence of the halogen-bonding interaction between FFL and the amino acid LEU430 of the protein HSA. While seemingly counterintuitive, this formation was characterized by the fact that although the electrostatic potential at the fluorine atom was more negative than that at the hydrogen atom, it was still less negative than that at the oxygen atom. Furthermore, comparison of the ESPs of propofol and fropofol also showed that although the two had similar molecular volumes (192 Å3 vs 189 Å3),18 the F substitution altered the electron distribution in the entire ligand differently than the OH. Not only was the ESP at the OH group in propofol different from that at the fluorine atom in fropofol, but the ESP at the benzene ring in propofol within the three different proteins was more negative than that in fropofol. This suggested that the benzene ring in propofol could form more attractive interactions with electron-poor partners than fropofol. This differential electron distribution further explains why the calculated binding interaction energies of propofol with three different proteins were all larger than those of fropofol. This also verified the hypothesis that was recently proposed by Woll and co-workers that the substitution of the 1hydroxyl by the fluorine atom could lead to a change in the electron distribution, which would affect the binding modes of the ligands to the proteins.18 Accordingly, the ESP analysis further demonstrated that propofol could bind with higher affinity than fropofol to the same receptor. The molecular recognition of target receptors did not only depend on the hydrogen bond interaction since the simultaneous cooperativity of noncovalent interactions could be responsible for the

involved in strong interactions with the ligand were also allowed to optimize in the ONIOM-2 calculations. Although the interaction energies were different between the ONIOM-1 and ONIOM-2 models, the trends were consistent on the whole, with most key residues still playing important roles in ligand−protein binding via large pairwise interaction energies. This could prove to be quite instructive for the future design of novel anesthetics as a means to improve ligand−protein binding affinity and efficacy. In particular, the interaction energies between PFL and key residues were stronger than those between FFL and corresponding residues within the proteins HSA and GLIC, while those between PFL and FFL with HSAF deviated relatively remarkably. This would illustrate that despite the dominant effect of hydrogen bonds involving the hydroxyl group for the specific recognition of target proteins, the total binding processes are due to more interactions than the hydrogen-bonding interaction alone. ESP Analysis. The molecular electrostatic potential (ESP)49 usually acts as a useful parameter to understand sites for electrophilic or nucleophilic attack reactions as well as nonbonded interactions. The blue regions of the surface illustration represent the positive electrostatic potential relating to electrophilic reactivity, while the red regions represent the negative electrostatic potential correlating with nucleophilic reactivity. In order to investigate whether the substitution of the 1-hydroxyl by the fluorine atom would result in electronic changes that would modulate ligand binding, ESP surfaces of propofol and fropofol within different proteins were analyzed and shown in Figure 10. As expected there is a positive electrostatic potential around the proton of the 1-hydroxyl group of propofol, while there is a negative electrostatic potential located at the oxygen atom of the 1-hydroxyl group and the benzene ring of propofol. Therefore, the hydrogen and the oxygen atom of the 1-hydroxyl group in PFL could act as the hydrogen bond donor and hydrogen bond acceptor, respectively. This can be observed from Table 2, which lists the hydrogen bond details. On the contrary, the electrostatic potential at the fluorine atom of FFL is negative, indicating that it could serve as a hydrogen bond acceptor. According to the definition of the hydrogen bond,50 halogen atoms can function as electron-rich (nucleophilic) sites and form net attractive interactions with electron-poor (electrophilic) partners due to their relatively high electro5894

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

The Journal of Physical Chemistry B



additional affinity and specificity of ligand binding. Together, it is the combination of subtle interactions which ensures that correct molecular recognition could occur. Based on this work, it can be more definitively stated that the different anesthetic endpoints caused by propofol and fropofol may be attributed to multiple specific interactions rather than merely the intermolecular hydrogen-bonding interactions alone.

AUTHOR INFORMATION

Corresponding Author

* Phone: +16504935000. Fax: +16508523423. Email:[email protected]. ORCID

Edward J. Bertaccini: 0000-0002-9062-0566



Notes

CONCLUSIONS To explore the significance of the hydroxyl group in determining the anesthetic efficacy of propofol in the present study, binding features of propofol and its fluorine substituted analogue, fropofol, within three important proteins (GLIC, HSA, and HSAF) have been investigated and compared in detail using a rigorous three-layer ONIOM (M06-2X/631+G*:PM6:AMBER) method. First, only the flexibility of the ligand was considered in the high level quantum mechanics layer within the ONIOM-1 scheme to take into account the effect of different protein environments on the binding features. Then, the flexibility of key residues involved in strong interactions with the ligand derived from ONIOM-1 calculations was incorporated into the high layer of the ONIOM-2 scheme. Both ONIOM-1 and ONIOM-2 calculations showed that substitution of the 1-hydroxyl in propofol by the fluorine atom in fropofol did not prevent binding to the protein. However, in general propofol had higher binding affinities than fropofol to the same protein cavities according to the calculated binding interaction energies regardless of hydrogen or halogen bond interactions. As a whole, the more robust ONIOM-2 calculations seem to better explain why propofol can lead to the desired anesthesia endpoint while fropofol cannot. The dominant role of key residues (ASN391, CYS392, ARG410, CYS438, and ARG485 in HSA; ILE202, TYR254, and ILE258 in GLIC; LEU24, ARG59, and GLU63 in HSAF) in contributing to the binding of propofol or fropofol within the three different proteins was highlighted by the calculations of pairwise interaction energies. ESP analysis also provided some complementary insights that the substitution of the 1-hydroxyl by a fluorine atom changed the electron distribution within the entire ligand allowing, for instance, the benzene ring in propofol to form more attractive interactions with electronpoor partners than fropofol. This clearly would result in different binding affinities for the same protein. Therefore, the molecular binding and subsequent anesthesia induced by propofol cannot be attributed only to the binding contribution by the hydroxyl group, but several other nonbonded interactions should also be taken into account collectively. Such vital information will provide new insights into the mechanism of action of propofol at a truly molecular level and offer new perspectives for the development of novel anesthetics targeting specific binding sites.



Article

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the National Natural Science Foundation of China (21371082, 21501074), Natural Science Foundation of Jiangsu Province (BK20141102, BK20151118), and 333 Project of Jiangsu Province (BRA2016518). The authors would also like to thankfully acknowledge the support of the Stanford University Department of Anesthesia and the United States Department of Veterans Affairs.



REFERENCES

(1) Franks, N. P.; Lieb, W. R. Molecular and cellular mechanisms of general anaesthesia. Nature 1994, 367, 607−614. (2) Hemmings, H. C.; Akabas, M. H.; Goldstein, P. A.; Trudell, J. R.; Orser, B. A.; Harrison, N. L. Emerging molecular mechanisms of general anesthetic action. Trends Pharmacol. Sci. 2005, 26, 503−510. (3) Chau, P. L. New insights into the molecular mechanisms of general anaesthetics. Br. J. Pharmacol. 2010, 161, 288−307. (4) Bertaccini, E. J. The Molecular Mechanisms of Anesthetic Action: Updates and Cutting Edge Developments from the Field of Molecular Modeling. Pharmaceuticals 2010, 3, 2178−2196. (5) Kotani, Y.; Shimazawa, M.; Yoshimura, S.; Iwama, T.; Hara, H. The experimental and clinical pharmacology of propofol, an anesthetic agent with neuroprotective properties. CNS Neurosci. Ther. 2008, 14, 95−106. (6) Trapani, G.; Altomare, C.; Sanna, E.; Biggio, G.; Liso, G. Propofol in anesthesia. Mechanism of action, structure-activity relationships, and drug delivery. Curr. Med. Chem. 2000, 7, 249−271. (7) Lingamaneni, R.; Krasowski, M. D.; Jenkins, A.; Truong, T.; Giunta, A. L.; Blackbeer, J.; MacIver, M. B.; Harrison, N. L.; Hemmings, H. C. Anesthetic properties of 4-iodopropofol: Implications for mechanisms of anesthesia. Anesthesiology 2001, 94, 1050− 1057. (8) Yip, G. M. S.; Chen, Z. W.; Edge, C. J.; Smith, E. H.; Dickinson, R.; Hohenester, E.; Townsend, R. R.; Fuchs, K.; Sieghart, W.; Evers, A. S.; et al. A propofol binding site on mammalian GABAA receptors identified by photolabeling. Nat. Chem. Biol. 2013, 9, 715−720. (9) Jayakar, S. S.; Dailey, W. P.; Eckenhoff, R. G.; Cohen, J. B. Identification of Propofol Binding Sites in a Nicotinic Acetylcholine Receptor with a Photoreactive Propofol Analog. J. Biol. Chem. 2013, 288, 6178−6189. (10) Jayakar, S. S.; Zhou, X. J.; Chiara, D. C.; Dostalova, Z.; Savechenkov, P. Y.; Bruzik, K. S.; Dailey, W. P.; Miller, K. W.; Eckenhoff, R. G.; Cohen, J. B. Multiple Propofol-binding Sites in a γAminobutyric Acid Type A Receptor (GABAAR) Identified Using a Photoreactive Propofol Analog. J. Biol. Chem. 2014, 289, 27456− 27468. (11) Chiara, D. C.; Gill, J. F.; Chen, Q.; Tillman, T.; Dailey, W. P.; Eckenhoff, R. G.; Xu, Y.; Tang, P.; Cohen, J. B. Photoaffinity Labeling the Propofol Binding Site in GLIC. Biochemistry 2014, 53, 135−142. (12) Laurent, B.; Murail, S.; Shahsavar, A.; Sauguet, L.; Delarue, M.; Baaden, M. Sites of Anesthetic Inhibitory Action on a Cationic LigandGated Ion Channel. Structure 2016, 24, 595−605. (13) Nury, H.; Van Renterghem, C.; Weng, Y.; Tran, A.; Baaden, M.; Dufresne, V.; Changeux, J. P.; Sonner, J. M.; Delarue, M.; Corringer, P. J. X-ray structures of general anaesthetics bound to a pentameric ligand-gated ion channel. Nature 2011, 469, 428−431.

ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.7b02079. Tables of layer divisions and PDB residues within the ONIOM-1 and -2 models and nonbonded interactions (PDF) Coordinate and output files for the binding sites (PDB1, PDB2, PDB3) 5895

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896

Article

The Journal of Physical Chemistry B (14) Bhattacharya, A. A.; Curry, S.; Franks, N. P. Binding of the General Anesthetics Propofol and Halothane to Human Serum Albumin. J. Biol. Chem. 2000, 275, 38731−38738. (15) Vedula, L. S.; Brannigan, G.; Economou, N. J.; Xi, J.; Hall, M. A.; Liu, R.; Rossi, M. J.; Dailey, W. P.; Grasty, K. C.; Klein, M. L.; et al. A unitary anesthetic binding site at high resolution. J. Biol. Chem. 2009, 284, 24176−24184. (16) Krasowski, M. D.; Jenkins, A.; Flood, P.; Kung, A. Y.; Hopfinger, A. J.; Harrison, N. L. General anesthetic potencies of a series of propofol analogs correlate with potency for potentiation of gammaaminobutyric acid (GABA) current at the GABA(A) receptor but not with lipid solubility. J. Pharmacol. Exp. Ther. 2001, 297, 338−351. (17) Krasowski, M. D.; Hong, X.; Hopfinger, A. J.; Harrison, N. L. 4D-QSAR analysis of a set of propofol analogues: Mapping binding sites for an anesthetic phenol on the GABA(A) receptor. J. Med. Chem. 2002, 45, 3210−3221. (18) Woll, K. A.; Weiser, B. P.; Liang, Q. S.; Meng, T.; McKinstryWu, A.; Pinch, B.; Dailey, W. P.; Gao, W. D.; Covarrubias, M.; Eckenhoff, R. G.Role for the Propofol Hydroxyl in Anesthetic Protein Target Molecular Recognition. ACS Chem. Neurosci. 2015, 6, 927− 935. (19) Vemparala, S.; Domene, C.; Klein, M. L. Computational Studies on the Interactions of Inhalational Anesthetics with Proteins. Acc. Chem. Res. 2010, 43, 103−110. (20) Murail, S.; Wallner, B.; Trudell, J. R.; Bertaccini, E.; Lindahl, E. Microsecond Simulations Indicate that Ethanol Binds between Subunits and Could Stabilize an Open-State Model of a Glycine Receptor. Biophys. J. 2011, 100, 1642−1650. (21) Liu, R.; Perez-Aguilar, J. M.; Liang, D.; Saven, J. G. Binding Site and Affinity Prediction of General Anesthetics to Protein Targets Using Docking. Anesth. Analg. 2012, 114, 947−955. (22) Leon, I.; Millan, J.; Cocinero, E. J.; Lesarri, A.; Castano, F.; Fernandez, J. A. Mimicking anaesthetic-receptor interaction: a combined spectroscopic and computational study of propofol···phenol. Phys. Chem. Chem. Phys. 2012, 14, 8956−8963. (23) Bertaccini, E. J.; Trudell, J. R.; Franks, N. P. The common chemical motifs within anesthetic binding sites. Anesth. Analg. 2007, 104, 318−324. (24) Bertaccini, E. J. The Molecular Mechanisms of Anesthetic Action: Updates and Cutting Edge Developments from the Field of Molecular Modeling. Pharmaceuticals 2010, 3, 2178−2196. (25) Bertaccini, E. J.; Yoluk, O.; Lindahl, E. R.; Trudell, J. R. Assessment of Homology Templates and an Anesthetic Binding Site within the γ-Aminobutyric Acid Receptor. Anesthesiology 2013, 119, 1087−1095. (26) Qiu, L.; Lin, J.; Bertaccini, E. J. Insights into the Nature of Anesthetic-Protein Interactions: An ONIOM Study. J. Phys. Chem. B 2015, 119, 12771−12782. (27) Chung, L. W.; Sameera, W. M. C.; Ramozzi, R.; Page, A. J.; Hatanaka, M.; Petrova, G. P.; Harris, T. V.; Li, X.; Ke, Z.; Liu, F.; et al. The ONIOM Method and Its Applications. Chem. Rev. 2015, 115, 5678−5796. (28) Vreven, T.; Byun, K. S.; Komáromi, I.; Dapprich, S.; Montgomery, J. A.; Morokuma, K.; Frisch, M. J. Combining Quantum Mechanics Methods with Molecular Mechanics Methods in ONIOM. J. Chem. Theory Comput. 2006, 2, 815−826. (29) 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 09 (Revision D.01); Gaussian, Inc.: Wallingford, CT, 2013. (30) Morokuma, K.; Wang, Q.; Vreven, T. Performance Evaluation of the Three-Layer ONIOM Method: Case Study for a Zwitterionic Peptide. J. Chem. Theory Comput. 2006, 2, 1317−1324. (31) Hohenstein, E. G.; Chill, S. T.; Sherrill, C. D. Assessment of the Performance of the M05-2X and M06-2X Exchange-Correlation Functionals for Noncovalent Interactions in Biomolecules. J. Chem. Theory Comput. 2008, 4, 1996−2000. (32) Stewart, J. J. P. Application of the PM6Method to Modeling Proteins. J. Mol. Model. 2009, 15, 765−805.

(33) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. A Second Generation Force Field for the Simulation of Proteins, Nucleic Acids, and Organic Molecules. J. Am. Chem. Soc. 1995, 117, 5179−5197. (34) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A.Development and Testing of a General Amber Force Field. J. Comput. Chem. 2004, 25, 1157−1174. (35) Singh, U. C.; Kollman, P. A. An approach to computing electrostatic charges for molecules. J. Comput. Chem. 1984, 5, 129− 145. (36) Kitisripanya, N.; Saparpakorn, P.; Wolschann, P.; Hannongbua, S. Binding of huperzine A and galanthamine to acetylcholinesterase, based on ONIOM method. Nanomedicine 2011, 7, 60−68. (37) 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. (38) Bissantz, C.; Kuhn, B.; Stahl, M. A medicinal chemist’s guide to molecular interactions. J. Med. Chem. 2010, 53, 5061−5084. (39) Butts, C. A.; Xi, J.; Brannigan, G.; Saad, A. A.; Venkatachalan, S. P.; Pearce, R. A.; Klein, M. L.; Eckenhoff, R. G.; Dmochowski, I. J. Identification of a fluorescent general anesthetic, 1-aminoanthracene. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 6501−6506. (40) Mahadevi, A. S.; Sastry, G. N. Cooperativity in Noncovalent Interactions. Chem. Rev. 2016, 116, 2775−2825. (41) Jeffrey, G. A. An Introduction to Hydrogen Bonding; Oxford University Press: Oxford, 1997. (42) Cavallo, G.; Metrangolo, P.; Milani, R.; Pilati, T.; Priimagi, A.; Resnati, G.; Terraneo, G. The Halogen Bond. Chem. Rev. 2016, 116, 2478−2601. (43) Wilcken, R.; Zimmermann, M. O.; Lange, A.; Joerger, A. C.; Boeckler, F. M. Principles and applications of halogen bonding in medicinal chemistry and chemical biology. J. Med. Chem. 2013, 56, 1363−1388. (44) Xu, Z.; Yang, Z.; Liu, Y.; Lu, Y.; Chen, K.; Zhu, W. Halogen bond: its role beyond drug-target binding affinity for drug discovery and development. J. Chem. Inf. Model. 2014, 54, 69−78. (45) Bader, R. F. W. Atoms in Molecules: A Quantum Theory; Oxford University Press: Oxford, 1990. (46) Biegler-König, F.; Schönbohm, J. Update of the AIM2000 -Program for Atoms in Molecules. J. Comput. Chem. 2002, 23, 1489− 1494. (47) Shoichet, B. K.; Leach, A. R.; Kuntz, I. D. Ligand solvation in molecular docking. Proteins: Struct. Proteins: Struct., Funct., Genet. 1999, 34, 4−16. (48) Derewenda, Z. S.; Lee, L.; Derewenda, U. The occurrence of C−H···O hydrogen bonds in proteins. J. Mol. Biol. 1995, 252, 248− 262. (49) Johnson, B. G.; Gill, P. M. W.; Pople, J. A. Computing Molecular Electrostatic Potentials with the PRISM Algorithm. Chem. Phys. Lett. 1993, 206, 239−246. (50) Arunan, E.; Desiraju, G. R.; Klein, R. A.; Sadlej, J.; Scheiner, S.; Alkorta, I.; Clary, D. C.; Crabtree, R. H.; Dannenberg, J. J.; Hobza, P.; et al. Defining the hydrogen bond: An account. Pure Appl. Chem. 2011, 83, 1619−1636. (51) Desiraju, G. R.; Ho, P. S.; Kloo, L.; Legon, A. C.; Marquardt, R.; Metrangolo, P.; Politzer, P.; Resnati, G.; Rissanen, K. Definition of the halogen bond. Pure Appl. Chem. 2013, 85, 1711−1713.

5896

DOI: 10.1021/acs.jpcb.7b02079 J. Phys. Chem. B 2017, 121, 5883−5896