Exploring Conformational Preferences of Leu-enkephalin Using the Conformational Search and Double-Hybrid DFT Energy Calculations.

The conformational preferences of Leu-enkephalin (Leu-Enk) were explored by the conformational search and density functional theory (DFT) calculations. By a combination of low-energy conformers of each residue, the initial structures of the neutral Leu-Enk were generated and optimized using the ECEPP3 force field in the gas phase. These structures were reoptimized at the HF/3-21G(d) and M06-2X levels of theory with 6-31G(d) and 6-31+G(d) basis functions. We finally located the 139 structures with the relative energy <10 kcal mol-1 in the gas phase, from which the structures of the corresponding zwitterionic Leu-Enk were generated and reoptimized at the M06-2X/6-31+G(d) level of theory using the implicit solvation model based on density (SMD) in water. The conformational preferences of Leu-Enk were analyzed using Gibbs free energies corrected by single-point energies calculated at the double-hybrid DSD-PBEP86-D3BJ/def2-TZVP level of theory in the gas phase and in water. The neutral Leu-Enk dominantly adopted a folded structure in the gas phase stabilized by three H-bonds with a βII'-bend-like motif at the Gly3-Phe4 sequence and a close contact between the side chains of Phe4 and Leu5. The zwitterionic Leu-Enk exhibited a folded structure in water stabilized by three H-bonds with double β-bends such as a βII' bend at the Gly2-Gly3 sequence and a βI bend at the Gly3-Phe4 sequence. The calculated ensemble-averaged distance between CGly2 α and CLeu5 α of the zwitterionic Leu-Enk in water is consistent with the value estimated from the simulated annealing using the distance constraints derived from nuclear Overhauser effect spectroscopy (NOESY) spectra in water. Interestingly, the preferred conformations of the neutral and zwitterionic Leu-Enk are new folded structures not predicted by earlier computational studies. According to the refined model of the zwitterionic Leu-Enk bound to δ-opioid receptor (δOR), there were favorable interactions of the terminal charged groups of Leu-Enk with the side chains of charged residues of δOR as well as a favorable CAryl···H interaction of the Phe4 residue of Leu-Enk with Trp284 of δOR. Hence, these favorable interactions would induce the folded structure of the zwitterionic Leu-Enk with double β-bends isolated in water into the "bioactive conformation" like an extended structure when binding to δOR.

Introduction

Enkephalins are endogenous opioid pentapeptides that are generated from the precursor proenkephalin (proEnk) in central and peripheral nervous systems and multiple organ systems.[1−5] The proEnk is proteolytically cleavaged to produce Met-enkephalin (Met-Enk, H-Tyr-Gly-Gly-Phe-Met-OH) and Leu-enkephalin (Leu-Enk, H-Tyr-Gly-Gly-Phe-Leu-OH). Both enkephalins have a higher binding affinity for the δ-opioid receptor (δOR), followed by the μ-opioid receptor and the ability to elicit signaling through G proteins and β-arrestin.[2,3]

Leu-Enk adopted three different conformations, namely, single β-bend,[6,7] double β-bend,[8,9] and extended[10−12] structures in the crystalline state. On the other hand, the peptide chain of Met-Enk was fully extended in the crystal.[12] For four decades, considerable NMR[13−26] and laser Raman/IR[15,27−31] spectroscopic measurements have been performed to investigate the preferred conformations of Leu-Enk in solution and in the solid state. It was known that Leu-Enk exists as a mixture of uncharged and zwitterionic species in dimethyl sulfoxide (DMSO), whereas it is a zwitterionic form in water.[15,21,30,31] According to the results from laser Raman/IR spectroscopic measurements,[27−29,31] Leu-Enk exists in different conformations in different solvents; i.e., it adopts folded β-bend structures in DMSO, whereas it exists as an ensemble of different conformations in water. In particular, it was suggested that the dominant conformation of Leu-Enk is a folded structure with a type I or II′ β-bend at the Gly3–Phe4 sequence with a C10 H-bond between the C=O of the Gly2 residue and the H–N of the Leu5 residue in DMSO by NMR experiments.[13,14,16,18] Furthermore, the side chains of Phe4 and Leu5 residues were predominantly populated in the gauche– (g–) configuration in DMSO.[14,16] It was suggested that this preferred conformation of the Leu5 residue brings the hydrophobic leucyl side chain toward the backbone far from the Tyr1 residue and hence plays a crucial role in the selective interaction of Leu-Enk with δOR.[16] In particular, the MD simulations using distance restraints obtained from the combination of nuclear Overhauser effect spectroscopy (NOESY) and rotating-frame nuclear Overhauser effect spectroscopy (ROESY) experiments supported that Leu-Enk exists as an ensemble of single-bent, double-bent, and extended structures in water, where the bent structures are β-bend-like structures without any H-bonds.[24] In addition, there have been several experimental works focused on the conformational preference of Leu-Enk in micelles and bicelles as membrane-mimicking environment in the absence and/or presence of ganglioside GM1.[23,24,32−38] The structures with single β-bend or double β-bend were found to be predominantly populated in phospholipid micelles and bicelles.

There have been many computational studies performed to identify the most preferred (or lowest-energy) conformation of Leu-Enk based on the conformational search and energy minimization[39−51] as well as the molecular dynamics (MD) simulations[23,36,51−57] using molecular force fields. Based on the conformational search and energy minimization with the ECEPP3[58] force field, three hairpin structures with a βII′ bend centered at the Gly3–Phe4 or Gly2–Gly3 sequences were suggested as the most preferred conformation of Leu-Enk in the absence of solvent.[40,43,48,50] In DMSO, MD simulations indicated that the zwitterionic Leu-Enk forms a clear close contact between N- and C-terminal ends, but there is no preferred conformation for the neutral form.[52] However, the simulations of the zwitterionic Leu-Enk in water suggested that the peptide backbone exists as a mixture of folded and unfolded forms[52−54] and the folded form possesses a N–H(Gly2) → O=C(Leu5) H-bond and a close head-to-tail separation.[53,54] In addition, the water molecules were suggested to play as bridges in the formation of β-turns at the Gly2–Gly3 and Gly3–Phe4 sequences of the folded structure in water.[55]

Limited studies have been focused on the structures, vibrational frequencies, vibrational absorption (VA), vibrational circular dichroism (VCD), and Raman intensities of the neutral Leu-Enk using density functional theory (DFT) methods.[59−62] From starting structures generated from three conformational motifs such as the extended, single β-bend, and double β-bend structures taken from X-ray data,[6,9,11,12] three probable structures with a single β-bend were located at the B3LYP/6-31G(d) level of theory, of which VA and VCD spectra are consistent with those of experiments in DMSO.[60,61] It was suggested that one (i.e., the conformer s4 in refs (59−61)) with a βI′ turn at the Gly2–Gly3 sequence of the three single β-bend structures is the most probable in DMSO.[61] Watson and Hirst identified 23 local minima of Leu-Enk at the EDF1/6-31+G(d) level of theory, of which initial structures were generated from the conformational search and energy minimizations using the CHARMM force field and from three X-ray structures.[62] They found that the lowest-energy structure is a kinked-extended structure with a γ-turn at the Gly2 residue, for which the amide I vibrational frequencies were calculated at the same level of theory. Recently, Frau et al. optimized the extended structure of the neutral Leu-Enk at eight DFT levels of theory with the def2-SVP basis set and the implicit solvation model in water and calculated the chemical descriptors.[63]

There are two crystal structures of δOR complexed with enkephalin-like ligands. The first one is the crystal structure of the human δ-opioid 7TM receptor complexed with the subtype-selective ligand naltrindole at 1.8 Å resolution (PDB ID 4N6H).[64] The second one is the active δOR crystal structure with the potent opioid agonist peptide KGCHM07 at 2.8 Å resolution (PDB ID 6PT2).[65] Sanfelice and Temussi[66] performed the molecular docking of δOR complexed with Leu-Enk using the crystal structure of PDB 4N6H, in which two conformations of Leu-Enk derived from NOE data in a DMSO/H2O mixture (ref (33)) were used starting structures for docking. Recently, Sharma et al. reported the molecular docking of δOR complexed with the Leu-Enk derivative bearing the Tyr1-ψ[(Z)CF=CH]-Gly2 substitution using the crystal structure of PDB 6PT2.[67] However, the “bioactive conformations” of Leu-Enk bound to δOR were different from each other and will be discussed in the Results and Discussion section.

As described above, considerable computational studies have been performed to investigate the conformational preferences of Leu-Enk in the absence and presence of solvent but no consensus conformation has been obtained, though all suggested conformations have similar structural features to those deduced from experiments. In particular, there was not any progress in high-resolution X-ray diffraction studies and quantum-mechanical computational studies at the high level of theory for Leu-Enk even in the 2000s to date. In the present work, we explored the conformational preferences of the neutral form of Leu-Enk in the absence of solvent (i.e., in the gas phase) and the zwitterionic form of Leu-Enk in water using DFT methods. In particular, the energies of preferred conformers of Leu-Enk for both forms were calculated at the double-hybrid DSD-PBEP86 level of theory with dispersion corrections,[68] which is one of the fifth-rung functionals. The conformational preferences were analyzed by conformational Gibbs free energies in the gas phase and in water.

Computational Methods

Chemical structures and torsion angles for the neutral and zwitterionic forms of Leu-Enk are defined in Figure . At the DFT level of theory, all geometry optimizations were performed using the M06-2X[69] functional and the implicit solvation model based on density (SMD)[70] implemented in the Gaussian 09 program.[71] Although the M06-2X functional is a hybrid-meta-generalized gradient approximation (GGA) functional with the improved medium-range correlation energy, it performs well in predicting noncovalent interactions of small molecules and relative stabilities of biological compounds.[72]

Figure 1

Chemical structure and torsion angles of the neutral Leu-Enk. The same definition of torsion angles was applied to the zwitterionic Leu-Enk.

As pointed out by Scheraga,[73,74] a large number of local minima are expected to exist on the multidimensional conformational energy surface of the flexible oligopeptide like Leu-Enk. However, the energy minimization procedure usually leads only to the local minimum closest to the starting structure, rather than to the global (or lowest-energy) minimum.[73,74] It is, therefore, necessary to have efficient methods for searching conformational space to locate the lowest-energy minimum.[73,74] As one of the methods for conformational search, the Scheraga group developed the buildup procedure and applied it to the conformational search of Met-Enk.[75] The buildup procedure starts with low-energy conformations of single residues and uses the combinations of these to build up feasible structures of dipeptides, tripeptides, and so on, carrying out energy minimization at each stage. However, some strategies were applied to eliminate the high-energy conformations due to the huge storage required for many low-energy conformations as longer peptides are built up.

In the present work, first, we performed the conformational search of the neutral form of Leu-Enk in the gas phase, in which initial structures were generated by the combinations of low-energy conformers of each terminally blocked residue without any eliminations. For Ac-X-NHMe, there were 31, 7, 18, and 27 local minima for X = Tyr, Gly, Phe, and Leu, respectively, with the relative energy (ΔEECEPP3) optimized by the ECEPP3 force field less than 5 kcal mol–1.[76] Hence, we generated 738,234 (31 × 7 × 7 × 18 × 27) initial structures for the neutral Leu-Enk by all combinations of these single-residue minima, which were followed by optimizations using the ECEPP3 force field and we obtained 2,633 local minima with ΔEECEPP3 < 5 kcal mol–1. Then, these 2,633 structures were reoptimized at the HF/3-21G(d) level of theory, which resulted in 984 local minima with the relative energy (ΔE) < 20 kcal mol–1. At the M06-2X/6-31G(d) level of theory, we reoptimized these 984 structures and identified 883 local minima with ΔE < 10 kcal mol–1. Then, it was followed by optimizations at the M06-2X/6-31+G(d) level of theory, which finally resulted in 139 structures with ΔE < 10 kcal mol–1.

From the 883 local minima of the neutral Leu-Enk with ΔE < 10 kcal mol–1 at the M06-2X/6-31G(d) level of theory, the corresponding initial structures of the zwitterionic Leu-Enk were built by the transfer of a proton from the C-terminal −COOH to the N-terminal −NH2 using GaussView.[77] For these initial structures of the zwitterionic Leu-Enk, we calculated the single-point energies at the SMD M06-2X/6-31G(d) level of theory in water and we obtained the 114 structures with ΔE < 10 kcal mol–1 in water. Then, these 114 structures were reoptimized at the SMD M06-2X/6-31+G(d) level of theory in water, which finally resulted in the 98 local minima with ΔE < 10 kcal mol–1 in water.

For all local minima of Leu-Enk with ΔE < 10 kcal mol–1 in the gas phase and in water, vibrational frequencies were calculated at the M06-2X/6-31+G(d) level of theory in the gas phase and the SMD M06-2X/6-31+G(d) level of theory in water at 25 °C and 1 atm. The scale factor used is 0.9440 that was chosen to reproduce the experimental frequency of 1707 cm–1 for the amide I band of N-methylacetamide in Ar and N2 matrixes.[78] The zero-point energy correction and the thermal energy corrections were employed in calculating the Gibbs free energy of each conformation. Here, the ideal gas, rigid rotor, and harmonic oscillator approximations were used for the translational, rotational, and vibrational contributions to the Gibbs free energy, respectively.[79]

In addition, single-point energies were calculated at the M06-2X/def2-TZVP level of theory for all local minima with ΔE < 10 kcal mol–1 in the gas phase and in water. We identified the 37 structures of the neutral Leu-Enk with the relative Gibbs free energy (ΔG) < 5 kcal mol–1 at the M06-2X/def2-TZVP//M06-2X/6-31+G(d) level of theory in the gas phase and the 46 structures of the zwitterionic Leu-Enk with ΔG < 5 kcal mol–1 at the M06-2X/def2-TZVP//SMD M06-2X/6-31+G(d) level of theory in water. Then, it was followed by the additional single-point energy calculations at the double-hybrid DSD-PBEP86 functional[68] with D3BJ dispersion corrections[80] and the def2-TZVP basis set for these low-free-energy structures in the gas phase and in water. The DSD-PBEP86-D3BJ level of theory is one of the highest fifth-rung functionals on the Jacob’s Ladder.[81] The DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) and M06-2X/def2-TZVP//M06-2X/6-31+G(d) levels of theory are abbreviated as DSD-dTZ and M062X-dTZ hereafter in the present work.

Recently, the DSD-dTZ level of theory exhibited good performance for relative energies of the 24 representative conformations of the alanine tetrapeptide in the gas phase[82] and the 45 local minima for the cationic nylon-3 dipeptide in water[83] with RMSD = 0.31 and 0.19 kcal mol–1, respectively, against the corresponding CCSD(T)/CBS-limit energies. In particular, the DSD-dTZ level of theory reproduced well relative conformational energies of six conformers of Ac-Ala-NHMe and seven conformers of Ac-Pro-NHMe with RMSD = 0.05 kcal mol–1 against the corresponding CCSD(T)/CBS-limit energies.[84] All three-dimensional (3D) graphics of optimized structures of Leu-Enk were prepared using PyMOL.[85]

Results and Discussion

Neutral Leu-Enk in the Gas Phase

There were the 29 local minima for the neutral Leu-Enk with the relative Gibbs free energy (ΔG) < 5 kcal mol–1 at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) (abbreviated as DSD-dTZ hereafter) level of theory in the gas phase. The torsion angles of these 29 local minima optimized at the M06-2X/6-31+G(d) level of theory are listed in Table S1 in the Supporting Information. The types of H-bond, β-bend, thermodynamic properties, and population of the corresponding structures calculated at the DSD-dTZ level of theory are shown in Table S2 in the Supporting Information, in which the thermodynamic properties and population at the M06-2X/def2-TZVP//M06-2X/6-31+G(d) (abbreviated as M062X-dTZ hereafter) level of theory are also included. The corresponding absolute electronic energies, enthalpies, Gibbs free energies, and single-point energies are listed in Table S3 in the Supporting Information.

Conformational Preferences of Neutral Leu-Enk

We identified the 11 structures for the neutral Leu-Enk with ΔG < 3 kcal mol–1 at the DSD-dTZ level of theory in the gas phase, whose torsion angles optimized at the M06-2X/6-31+G(d) level of theory are listed in Table . The types of H-bond and thermodynamic properties of the corresponding structures are shown in Table . The 3D structures of the first six preferred conformers with ΔG < 2 kcal mol–1 are depicted in Figure , whose Cartesian coordinates are also listed in the Supporting Information.

Table 1

Torsion Angles of 11 Local Minima for the Neutral Leu-EnK Optimized at the M06-2X/6-31+G(d) Level of Theory in the Gas Phasea

	Tyr1		Gly2		Gly3		Phe4			Leu5
conf.	ψ1	χ11	ϕ2	ψ2	ϕ3	ψ3	ϕ4	ψ4	χ41	ϕ5	ψ5	χ51
n01	–3	64	115	–18	84	–73	–74	–15	66	–161	–42	177
n02	–2	66	110	–4	77	–90	–84	–13	61	–154	65	–175
n03	–1	66	113	–15	81	–76	–66	–33	–59	–143	–36	–64
n04	–3	65	115	–16	82	–74	–68	–23	–56	–158	–40	177
n05	–1	64	113	–15	82	–71	–97	–2	–58	–172	8	43
n06	–2	64	112	–13	82	–76	–106	10	61	–151	–20	63
n07	–2	64	112	–13	84	–71	–107	5	70	–166	6	46
n08	–2	64	114	–14	83	–70	–109	14	–57	–166	–20	61
n09	–2	65	114	–15	81	–76	–59	–37	176	–151	–31	–179
n10	–2	64	115	–17	84	–73	–76	–13	67	–160	–45	178
n11	2	67	114	–24	81	–75	–82	101	178	56	34	–49

Torsion angles (°) are defined in Figure . Only 11 local minima with the relative Gibbs free energy (ΔG) < 3 kcal mol–1 are listed.

Table 2

Type of H-bond, β-Bend, Thermodynamic Properties, and Population of 11 Local Minima for the Neutral Leu-EnK Calculated at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) Level of Theory in the Gas Phasea

conf.	H-bond typeb	β-bendc	ΔEd	ΔHe	ΔGf	wg
n01	OH1 → 3, 4 → 2, OH5 → 1		0.00	0.00	0.00	63.6
n02	OH1 → 3, 4 → 2, 5 → 2, OH5 → 1	βII′34	2.48	2.55	1.18	8.6
n03	OH1 → 3, 4 → 2, OH5 → 1		1.95	2.07	1.23	7.9
n04	OH1 → 3, 4 → 2, OH5 → 1		2.06	2.09	1.25	7.7
n05	OH1 → 3, 4 → 2, OH5 → 1		2.79	3.03	1.83	2.9
n06	OH1 → 3, 4 → 2, OH5 → 1		1.86	1.97	1.85	2.8
n07	OH1 → 3, 4 → 2, OH5 → 1		0.47	0.70	2.19	1.6
n08	OH1 → 3, 4 → 2, OH5 → 1		3.74	3.93	2.25	1.4
n09	OH1 → 3, 4 → 2, OH5 → 1		4.20	3.95	2.44	1.0
n10	OH1 → 3, 4 → 2, OH5 → 1		1.47	1.80	2.74	0.6
n11	OH1 → 3, 4 → 2, OH5 → 1		2.96	3.02	2.90	0.5

Torsion angles (°) are listed in Table . Only 11 local minima with the relative Gibbs free energy (ΔG) < 3 kcal mol–1 are listed.

Each H-bond type n → m stands for the H-bond between the H donor (e.g., the amide H atom for backbone) of the residue n and the H acceptor (e.g., the carbonyl O atom for backbone) of the residue m. In addition, OH1 and OH5 represent the hydroxyl H atom of the side chain of the Tyr1 residue and the carboxylic H atom of the Leu5 residue, respectively. 5OH represents the O atom of the carboxylic OH group of the Leu5 residue.

βII′34 stands for the type II′ β-bend at the Gly3–Phe4 sequence, which is stabilized by the 5 → 2 H-bond.

Relative electronic energies in kcal mol–1.

Relative enthalpies in kcal mol–1 at 25 °C.

Relative Gibbs free energies in kcal mol–1 at 25 °C and 1 atm.

The population of each conformer was calculated by its ΔG at 25 °C.

Figure 2

Preferred structures of the neutral Leu-EnK optimized at the M06-2X/6-31+G(d) level of theory in the gas phase. The relative Gibbs free energy (ΔG in kcal mol–1) at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) level of theory is shown in parentheses behind each conformation name. For clarity, all nonpolar hydrogen atoms are omitted. All H-bonds are represented by dotted lines with distances in Å.

In the gas phase, the most preferred conformer n01 with ΔE = 0.00 and ΔG = 0.00 kcal mol–1 (populated at 63.6%) adopted a folded structure with three H-bonds such as an OH1 → 3 type of the 16-membered (i.e., C16) side chain-to-backbone H-bond between O–H(Tyr1) and C=O(Gly3) with the distance of 1.80 Å; a 4 → 2 type of the C7 backbone-to-backbone H-bond between N–H(Phe4) and C=O(Gly2) with the distance of 2.09 Å; and an OH5 → 1 type of the C16 backbone-to-backbone H-bond between O–H(Leu5) and C=O(Tyr1) with the distance of 1.69 Å (Figure ). In particular, the formation of the second C7 4 → 2 H-bond seems to affect the torsion angle ψ3 of Gly3 moved to the less negative value of −73° (cf. the standard value of −120° reported in the literature[86]) and to prevent the formation of a type II′ β-bend (i.e., βII′34) at the Gly3–Phe4 sequence (therefore, it is a βII′34-like motif without a 5 → 2 H-bond). In addition, there was a close contact with a distance of about 3.5 Å between the phenyl ring of the side chain of Phe4 and the Hβ proton of the side chain of Leu5.

The second preferred conformer n02 with ΔG = 1.18 kcal mol–1 (populated at 8.6%) exhibited a folded structure with four H-bonds such as OH1 → 3, 4 → 2, 5 → 2, and OH5 → 1 types of H-bonds with the distances of 1.82, 2.23, 2.11, and 1.73 Å, respectively (Figure ). The structure of conformer n02 is quite similar to that of the most preferred conformer n01, except for the formation of a bifurcated H-bond with 4 → 2 and 5 → 2 types of H-bond and a βII′34 bend at the Gly3–Phe4 sequence (Table ). However, conformer n02 was 2.48 kcal mol–1 less stable in conformational energy (ΔE) than conformer n01.

The third-to-eleventh preferred conformer n03-to-n11 adopted a similar folded structure to that of the most preferred conformer n01 with three H-bonds such as OH1 → 3, 4 → 2, and OH5 → 1 types of H-bond with the distances of 1.80–1.84, 2.01–2.09, and 1.69–1.80 Å, respectively (Figure ). However, there were somewhat greater differences in torsion angles of the backbone and side chains of Phe4 and Leu5, which would result in the increase of ΔE values as 0.47 to 4.20 kcal mol–1 and ΔG values as 1.23 to 2.90 kcal mol–1 (populated at 7.9 to 0.5%, respectively) compared with that of the most preferred conformer n01 (Table ). In particular, the seventh conformer n07 exhibited a folded structure quite similar to conformer n06 but it had the second lowest energy of 0.47 kcal mol–1.

Hence, the neutral Leu-Enk dominantly adopted a folded structure stabilized by a side chain-to-backbone (OH1 → 3) H-bond and two backbone-to-backbone (4 → 2 and OH5 → 1) H-bonds with a βII′-like motif at the Gly3–Phe4 sequence and a close contact between the side chains of Phe4 and Leu5 in the gas phase.

Comparison with Earlier Works

There have been several computational studies performed to identify the lowest-energy conformation of the neutral Leu-Enk using the conformational search and energy minimization with the ECEPP3 force field[43,50] and DFT methods.[59,62] Using the initials of author’s last names, the conformers suggested by the ECEPP3 force field were denoted as MM1–MM4[43] and VG;[50] and those by DFT methods as J1–J3[59] and WH1–WH3[62] at B3LYP/6-31G(d) and EDF1/6-31+G(d) levels of theory, respectively. Conformer J1, J2, and J3 correspond to s10ndp, s4ndp, and s7ndp in ref (59), respectively; and conformer WH1, WH2, and WH3 correspond to qm-1-ke, qm-2-ebt, and qm-3-e in ref (62), respectively. In addition, the first nine preferred conformers predicted by the ECEPP3 force field in the present work were denoted as ECEPP3-1 to ECEPP3-9. All of these structures were reoptimized at the M06-2X/6-31+G(d) level of theory.

The torsion angles of optimized structures that were initially generated by ECEPP3 force field and DFT methods are listed in Table . The types of H-bond and thermodynamic properties of the corresponding structures calculated at the DSD-dTZ level of theory are shown in Table . The corresponding data calculated at the M062X-dTZ level of theory are shown in Table S4 in the Supporting Information. The corresponding absolute electronic energies, enthalpies, Gibbs free energies, and single-point energies are listed in Table S5 in the Supporting Information. The corresponding optimized structures are shown in Figures S1–S3 in the Supporting Information.

Table 3

Torsion Angles of Local Minima of the Neutral Leu-EnK Optimized at the M06-2X/6-31+G(d) Level of Theory from the Structures Obtained by the ECEPP Force Field and DFT Methods in the Gas Phasea

		Tyr1		Gly2		Gly3		Phe4			Leu5
conf.		ψ1	χ11	ϕ2	ψ2	ϕ3	ψ3	ϕ4	ψ4	χ41	ϕ5	ψ5	χ51
ECEPP3-1	present work	–38	165	113	0	88	–76	–83	–58	–179	–122	157	–179
ECEPP3-2	present work	–38	164	116	7	86	–84	–83	–41	–61	–144	148	–180
ECEPP3-3	present work	–36	161	112	5	83	–75	–90	–47	174	–120	–16	–176
ECEPP3-4	present work	133	–171	–84	63	70	–110	–62	–40	177	–171	–61	172
ECEPP3-5	present work	–25	173	58	–124	–109	1	–176	147	40	–57	150	–180
ECEPP3-6	present work	146	–171	–81	54	75	–97	–96	–13	64	–117	148	–179
ECEPP3-7	present work	–38	164	115	10	86	–88	–87	–32	65	–146	141	178
ECEPP3-8	present work	–21	–165	51	–137	–107	33	–169	137	23	–56	144	–178
ECEPP3-9	present work	122	180	–91	34	87	–48	–117	–6	66	–167	145	180
MM1b	ref (43)	132	176	–119	29	77	–95	–85	–10	62	–119	143	176
MM2b	ref (43)	–25	173	58	–124	–109	1	–176	147	40	–57	150	–180
MM3b	ref (43)	–140	–167	–87	64	86	7	–164	169	64	–151	–33	–175
MM4b	ref (43)	–169	58	–75	94	107	–15	–163	–174	–151	–138	162	–66
VGc	ref (50)	–31	–65	169	148	57	–137	–97	16	55	–123	15	–60
J1d	ref (59)	92	–59	63	14	103	–19	–84	84	–64	–81	58	–53
J2d	ref (59)	85	–64	63	13	103	–16	–72	102	–69	–145	162	–70
J3d	ref (59)	–26	–66	63	24	103	–23	–88	76	–65	–80	58	–51
WH1e	ref (62)	–39	170	74	2	109	–31	–165	127	–178	–161	156	–174
WH2e	ref (62)	–29	–178	57	–133	–108	35	–55	144	177	–150	156	–177
WH3e	ref (62)	137	170	–151	153	158	149	–163	169	–162	–54	148	179

Torsion angles (°) are defined in Figure .

Optimized from the ECEPP2 structures of Meirovitch and Meirovitch.

Optimized from the ECEPP3 structures of Vengadesan and Gautham using the MOLS technique.

Optimized from the structures of Jalkanen obtained at the B3LYP/6-31G(d) level of theory.

Optimized from the structures of Watson and Hirst obtained at the EDF1/6-31+G(d) level of theory.

Table 4

Type of H-Bond, β-Bend, and Thermodynamic Properties of Local Minima for the Neutral Leu-EnK Optimized from the Structures Obtained by the ECEPP Force Field and DFT Methods at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) Level of Theory in the Gas Phaseah

conf.b	H-bond typec	β-bendd	ΔEe	ΔHf	ΔGg
ECEPP3-1	OH1 → 4, 4 → 3, OH5 → 1		10.64	10.61	9.46
ECEPP3-2	OH1 → 4, 4 → 3, OH5 → 1		10.50	10.93	8.85
ECEPP3-3	OH1 → 4, 4 → 3, OH5 → 1		11.94	11.86	11.31
ECEPP3-4	OH1 → 4, 2 → 5, 3 → 1, 5 → 2	βII′34	10.32	10.72	11.25
ECEPP3-5	OH1 → 3, 4 → 1	βII′23	17.75	17.47	15.62
ECEPP3-6	OH1 → 4, 2 → 5, 3 → 1, 4 → 2, 5 → 2	βII′34	4.22	4.57	4.76
ECEPP3-7	OH1 → 4, 4 → 2, OH5 → 1		8.89	8.85	6.60
ECEPP3-8	4 → 1	βII′23	20.59	19.97	17.65
ECEPP3-9	OH1 → 3, 2 → 5, 4 → 3		6.90	6.92	5.88
MM1	OH1 → 3, 2 → 5, 4 → 2, 5 → 2	βII′34	7.24	7.44	5.11
MM2	OH1 → 3, 4 → 1	βII′23	17.75	17.47	15.62
MM3	1 → 4, 3 → 1		17.69	17.45	13.78
MM4	1 → 4, 3 → 1		17.08	16.72	14.30
VG	2 → OH5, 5 → 2	βII′34	11.50	11.32	7.42
J1	1 → 4, 4 → 1, 5 → 3, OH5 → 4	βI′23	13.99	14.18	12.16
J2	1 → 4, 4 → 1, 5 → 3	βI′23	13.70	13.62	10.73
J3	4 → 1, 5 → 3, OH5 → 4	βI′23	11.59	11.22	9.39
WH1	OH1 → 3, 4 → 1	βI′23	12.23	12.31	10.61
WH2	OH1 → 4, 4 → 1	βII′23	11.60	11.41	9.26
WH3	OH1 → 3		18.08	17.58	13.52

Torsion angles (°) are listed in Table .

The definition of conformations is noted in footnotes b–e of Table .

Each H-bond type n → m stands for the H-bond between the H donor (e.g., the amide H atom for backbone) of the residue n and the H acceptor (e.g., the carbonyl O atom for backbone) of the residue m. In addition, OH1 and OH5 represent the hydroxyl H atom of the side chain of the Tyr1 residue and the carboxylic H atom of the Leu5 residue, respectively. 5OH represents the O atom of the carboxylic OH group of the Leu5 residue.

βII′23 and βI′23 stand for the type II′ and I′ β-bend at the Gly2–Gly3 sequence, which are stabilized by the 4 → 1 H-bond. βII′34 stands for the type II′ β-bend at the Gly3–Phe4 sequence, which is stabilized by the 5 → 2 H-bond.

Relative electronic energies in kcal mol–1.

Relative enthalpies in kcal mol–1 at 25 °C.

Relative Gibbs free energies in kcal mol–1 at 25 °C and 1 atm.

The population of each conformer was calculated by its ΔG at 25 °C.

The most preferred conformer was ECEPP3-6, MM1, VG, J3, and WH2 for each group of computations with ΔG = 4.76, 5.11, 7.42, 9.39, and 9.26 kcal mol–1, respectively, which are shown in Figure . Conformer ECEPP3-6, MM1, and VG adopted different folded structures with different types of H-bond (i.e., five H-bonds such as OH1 → 4, 2 → 5, 3 → 1, 4 → 2, and 5 → 2 types of H-bond for ECEPP3-6; four H-bonds such as OH1 → 3, 2 → 5, 4 → 2, and 5 → 2 types of H-bond for MM1; and two H-bonds such as 2 → OH5 and 5 → 2 types of H-bond for VG), despite the same βII′34 bend at the Gly3–Phe4 sequence stabilized by the 5 → 2 H-bond. However, conformer J3 and WH2 exhibited a folded structure with βI′23 and βII′23 bends at the Gly2–Gly3 sequence, respectively, stabilized by the 4 → 1 H-bond. The former had three H-bonds such as 4 → 1, 5 → 3, and OH5 → 4 types of H-bond, whereas the latter had two H-bonds such as OH1 → 4 and 4 → 1 types of H-bond.

Figure 3

Structures of the neutral Leu-EnK obtained by (i) the search methods with the ECEPP force field (ECEPP3-6 in the present work, MM1 from ref (43), and VG from ref (50)) and (ii) the DFT calculations (J3 from ref (59) and WH2 from ref (62)), which were optimized at the M06-2X/6-31+G(d) level of theory in the gas phase. The relative Gibbs free energy (ΔG in kcal mol–1) at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) level of theory is shown in parentheses behind each conformation name. For clarity, all nonpolar hydrogen atoms are omitted. All H-bonds are represented by dotted lines with distances in Å.

In particular, all of these five structures possessed different types of H-bonds despite the formation of the common β-bend motifs, which are also quite different from that of the most preferred conformer n01, as described above (see Figures and 3). Hence, the most preferred conformer n01 is a new folded structure for the neutral Leu-Enk in the gas phase, not predicted in earlier works.

Zwitterionic Leu-Enk in Water

There were the 43 local minima for the zwitterionic Leu-Enk with the relative Gibbs free energy (ΔGw) < 5 kcal mol–1 at the DSD-PBEP86-D3BJ/def2-TZVP//SMD M06-2X/6-31+G(d) (abbreviated as DSD-dTZ/SMD hereafter) level of theory in water. The torsion angles of these 43 local minima optimized at the SMD M06-2X/6-31+G(d) level of theory are listed in Table S6 in the Supporting Information. The types of H-bond, β-bend, thermodynamic properties, and population of the corresponding structures calculated at the DSD-dTZ/SMD level of theory are shown in Table S7 in the Supporting Information, in which the thermodynamic properties and population at the M06-2X/def2-TZVP//SMD M06-2X/6-31+G(d) (abbreviated as M062X-dTZ/SMD hereafter) level of theory are also included. The corresponding absolute electronic energies, enthalpies, Gibbs free energies, and single-point energies are listed in Table S8 in the Supporting Information.

Conformational Preferences of Zwitterionic Leu-Enk

We found the 16 structures for the zwitterionic Leu-Enk with ΔGw < 3 kcal mol–1 at the DSD-dTZ/SMD level of theory in water, whose torsion angles optimized at the SMD M06-2X/6-31+G(d) level of theory in water are listed in Table . The types of H-bond and thermodynamic properties of the corresponding structures are shown in Table . The 3D structures of the first eight preferred conformers with ΔGw < 2 kcal mol–1 in water are depicted in Figure , whose Cartesian coordinates are also listed in the Supporting Information.

Table 5

Torsion Angles of 16 Local Minima for the Zwitterionic Leu-EnK Optimized at the SMD M06-2X/6-31+G(d) Level of Theory in Watera

	Tyr1		Gly2		Gly3		Phe4			Leu5
conf.	ψ1	χ11	ϕ2	ψ2	ϕ3	ψ3	ϕ4	ψ4	χ41	ϕ5	ψ5	χ51
zw01	159	63	55	–141	–68	–13	–76	–10	–50	69	23	–40
zw02	164	64	68	–122	–78	7	–152	170	52	–58	134	180
zw03	163	61	72	–118	–80	6	–155	167	51	–59	140	–61
zw04	145	175	63	–134	–63	–19	–84	–3	–61	72	15	–41
zw05	143	178	–80	93	–66	–14	–69	–12	68	–133	–179	56
zw06	142	–175	–79	104	–69	–10	–66	–13	59	–89	38	–57
zw07	–54	–69	169	144	57	–136	–68	–22	67	–82	–43	–174
zw08	141	179	–98	19	88	–46	–115	–12	63	–162	148	–178
zw09	129	176	–125	32	89	–70	–95	2	64	–161	170	76
zw10	142	180	–99	22	87	–45	–118	–11	63	–160	144	179
zw11	173	–167	165	–175	–69	–10	–69	–6	59	–151	1	–168
zw12	131	173	–110	17	85	–71	–68	–20	66	–162	144	–179
zw13	–180	–160	84	–69	–59	–37	–99	140	173	–125	–35	–61
zw14	131	174	–107	14	83	–70	–56	–39	175	–157	150	–177
zw15	142	–180	–101	18	86	–60	–79	–40	–63	–152	165	–64
zw16	149	–166	–150	172	–70	–9	–71	–11	68	–141	13	55

Torsion angles (°) are defined in Figure . Only 16 local minima with the relative Gibbs free energy (ΔGw) < 3 kcal mol–1 are listed.

Table 6

Type of H-bond, β-Bend, Thermodynamic Properties, and Population of 16 Local Minima for the Zwitterionic Leu-EnK Calculated at the DSD-PBEP86-D3BJ/def2-TZVP//SMD M06-2X/6-31+G(d) Level of Theory in Watera

conf.	H-bond typeb	β-bendc	ΔEwd	ΔHwe	ΔGwf	wg
zw01	1 → 5, 4 → 1, 5 → 2	βII′23, βI34	0.84	0.61	0.00	48.1
zw02	1 → 4, 1 → 5, 4 → 1	βII′23	1.15	1.39	0.69	15.0
zw03	1 → 4, 1 → 5, 4 → 1	βII′23	1.78	1.74	0.69	15.0
zw04	1 → 5, 4 → 1, 5 → 2	βII′23, βI34	3.66	3.50	1.54	3.6
zw05	OH1 → 5, 5 → 2	βI34	0.00	0.00	1.70	2.8
zw06	OH1 → 5, 5 → 2	βI34	1.16	1.46	1.77	2.4
zw07	1 → 5, 2 → 5, 5 → 2	βII′34	3.49	3.38	1.94	1.8
zw08	OH1 → 3, 2 → 5, 4 → 2		1.65	1.69	1.96	1.8
zw09	OH1 → 3, 2 → 5, 4 → 2		1.71	1.67	2.12	1.4
zw10	OH1 → 3, 2 → 5, 4 → 2		3.34	3.26	2.15	1.3
zw11	OH1 → 5, 5 → 2	βI34	1.25	1.44	2.27	1.1
zw12	OH1 → 3, 2 → 5, 4 → 2		1.75	2.02	2.27	1.0
zw13	OH1 → 5, 2 → 4, 3 → 1		1.55	1.52	2.50	0.7
zw14	OH1 → 3, 2 → 5, 4 → 2		3.05	3.20	2.77	0.5
zw15	OH1 → 3, 2 → 5, 4 → 2		3.35	3.70	2.97	0.3
zw16	OH1 → 5, 5 → 2	βI34	2.37	2.15	2.98	0.3

Torsion angles (°) are listed in Table . Only 16 local minima with the relative Gibbs free energy (ΔGw) < 3 kcal mol–1 are listed.

Each H-bond type n → m stands for the H-bond between the H donor (e.g., the amide H atom for backbone) of the residue n and the H acceptor (e.g., the carbonyl O atom for backbone) of the residue m. In addition, OH1 represents the hydroxyl H atom of the side chain of the Tyr1 residue.

βII′23 stands for the type II′ β-bend at the Gly2–Gly3 sequence, which are stabilized by the 4 → 1 H-bond. βI34 and βII′34 stand for the type I and II′ β-bend at the Gly3–Phe4 sequence, respectively, which are stabilized by the 5 → 2 H-bond.

Relative electronic energies in kcal mol–1.

Relative enthalpies in kcal mol–1 at 25 °C.

Relative Gibbs free energies in kcal mol–1 at 25 °C and 1 atm.

The population of each conformer was calculated by its ΔGw at 25 °C.

Figure 4

Preferred structures of the zwitterionic Leu-EnK optimized at the SMD M06-2X/6-31+G(d) level of theory in water. The relative Gibbs free energy (ΔG in kcal mol–1) at the DSD-PBEP86-D3BJ/def2-TZVP//SMD M06-2X/6-31+G(d) level of theory is shown in parentheses behind each conformation name. For clarity, all nonpolar hydrogen atoms are omitted. All H-bonds are represented by dotted lines with distances in Å.

In water, the most preferred conformer zw01 with ΔGw = 0.00 kcal mol–1 (populated at 48.1%) exhibited a folded structure with three H-bonds such as a 1 → 5 type of the C17 head-to-tail H-bond between zwitterionic N–H+(Tyr1) and C=O–(Leu5) species with a short distance of 1.72 Å; 4 → 1 and 5 → 2 types of the C10 backbone-to-backbone H-bond between N–H(Phe4) and C=O(Tyr1) and between N–H(Leu5) and C=O(Gly2) with the same distance of 2.04 Å, respectively. In particular, the two H-bonds of the latter induced the formation of double β-bends such as a βII′23 bend at the Gly2–Gly3 sequence and a βI34 bend at the Gly3–Phe4 sequence, respectively (Figure ).

The second and third preferred conformers zw02 and zw03 were iso-Gibbs free energy structures with ΔGw = 0.69 kcal mol–1 (populated at 15.0%). Both conformers adopted a folded structure with three H-bonds such as a 1 → 4 type of the C14 backbone-to-backbone H-bond between N–H+(Tyr1) and C=O(Phe4); a 1 → 5 type of the C17 head-to-tail H-bond between N–H+(Tyr1) and C=O–(Leu5); and a 4 → 1 type of the C10 backbone-to-backbone H-bond between N–H(Phe4) and C=O(Tyr1). The distances of 1 → 4, 1 → 5, and 4 → 1 types of H-bond were calculated as 1.82, 2.29, and 2.20 Å for conformer zw02, respectively; and 1.82, 2.33, and 2.25 Å for conformer zw03, respectively (Figure ). In particular, two conformer zw02 and zw03 possessed in common a βII′23 bend at the Gly2–Gly3 sequence stabilized by the 1 → 4 H-bond but different orientations of Phe4 and Leu5 residues.

The fourth preferred conformer zw04 adopted a folded structure quire similar to the most preferred conformer zw01 stabilized by three H-bonds such as 1 → 5, 4 → 1, and 5 → 2 types of H-bond with the distances of 1.69, 2.05, and 1.97 Å, respectively, but having the different orientation of the Tyr1 residue (Figure ). However, conformer zw04 was 1.54 kcal mol–1 less stable in ΔGw (populated at 3.6%) than conformer zw01, although both conformers had the same types of βII′23 and βI34 bends.

Conformer zw05 and zw06 were nearly iso-Gibbs free energy structures with ΔGw = 1.70 and 1.77 kcal mol–1, respectively (populated at 2.8 and 2.4%, respectively). Both conformers adopted a folded structure stabilized by OH1 → 5 and 5 → 2 types of H-bond, the latter of which induced the formation of a βI34 bend at the Gly3–Phe4 sequence (Figure ). The distances of OH1 → 5 and 5 → 2 types of H-bond were calculated as 1.71 and 2.03 Å for conformer zw05, respectively; 1.70 and 1.94 Å for conformer zw06, respectively. In particular, it should be noted that conformer zw05 was 0.84 kcal mol–1 more stable in ΔEw than the most preferred conformer zw01 (Table ).

The seventh preferred conformer zw07 with a βII′34 bend at the Gly3–Phe4 sequence was stabilized by three H-bonds such as 1 → 5, 2 → 5, and 5 → 2 types of backbone-to-backbone H-bond with the distances of 1.70, 1.90, and 2.07 Å, respectively, but it was 1.94 kcal mol–1 less stable in ΔGw (populated at 1.8%) than conformer zw01. The eighth preferred conformer zw08 was stabilized by three H-bonds such as OH1 → 3, 2 → 5, and 4 → 2 types of H-bond with the distances of 1.97, 1.85, and 2.22 Å, respectively. However, this conformer did not have any β-bend and it was 1.96 kcal mol–1 less stable in ΔGw (populated at 1.8%) than conformer zw01.

Hence, the zwitterionic Leu-Enk predominantly adopted a folded structure stabilized by three 1 → 5, 4 → 1, and 5 → 2 types of H-bond in water, of which the second and third H-bonds induced the formation of double β-bends such as a βII′23 bend at the Gly2–Gly3 sequence and a βI34 bend at the Gly3–Phe4 sequence, respectively. This most preferred conformer of the zwitterionic Leu-Enk in water was quite different from that of the neutral Leu-Enk in the gas, as described above. It was suggested that the dominant conformation of Leu-Enk is a folded structure with a type βI34 or βII′34 bend stabilized by the 5 → 2 H-bond in DMSO by NMR experiments.[13,14,16,18] Furthermore, the side chains of Phe4 and Leu5 residues were predominantly populated in the gauche– (g–) configuration in DMSO.[14,16] In addition, the distance restraints obtained from the combination of NOESY and ROESY experiments supported that Leu-Enk exists as an ensemble of single β-bend, double β-bend, and extended structures in water.[24] These conformational preferences deduced by NMR experiments in DMSO and water are reasonably consistent with those of the preferred conformer zw01, zw02, and zw03 (populated at 48.1, 15.0, and 15.0%, respectively) in water in the present work (see Table ).

Several studies have focused on the preferred conformation of the zwitterionic Leu-Enk in the crystalline state and in water. Leu-Enk adopted three different conformations, namely, single β-bend (CSD ID LENKPH11[7]), double β-bend (CSD ID GEWWAG[8]), and extended (CSD ID BIXNIF10[11] and FABJEX[12]) structures in the crystalline state. The three representative conformations (denoted by MD-AA1, MD-AA2, and MD-AA3 in the present work) were suggested from the MD simulations of the zwitterionic Leu-Enk in water.[53] All of these seven structures were reoptimized at the SMD M06-2X/6-31+G(d) level of theory in water, whose torsion angles are listed in Table . The types of H-bond and thermodynamic properties of the corresponding structures calculated at the DSD-dTZ/SMD level of theory are shown in Table . The corresponding data calculated at the M062X-dTZ/SMD level of theory are shown in Table S9 in the Supporting Information. The corresponding absolute electronic energies, enthalpies, Gibbs free energies, and single-point energies are listed in Table S10 in the Supporting Information. The corresponding optimized structures are shown in Figures S4 and S5 in the Supporting Information.

Table 7

Torsion Angles of Local Minima of the Zwitterionic Leu-EnK Optimized at the SMD M06-2X/6-31+G(d) Level of Theory from the Structures Obtained by X-ray Diffraction, MD Simulation, and Molecular Docking in Watera

		Tyr1		Gly2		Gly3		Phe4			Leu5
conf.		ψ1	χ11	ϕ2	ψ2	ϕ3	ψ3	ϕ4	ψ4	χ41	ϕ5	ψ5	χ51
LENKPH11b	ref (7)	114	–56	62	21	99	7	–139	128	–62	–141	–20	–67
GEWWAGb	ref (8)	124	170	–67	–22	–66	–15	–87	–10	–58	–60	141	–177
BIXNIF10b	ref (11)	160	62	170	–169	–178	–177	–151	151	–58	–58	159	–61
FABJEXb	ref (12)	164	60	–173	–177	–177	–179	–150	138	–58	–167	170	73
MD-AA1c	ref (53)	142	–177	–161	137	80	–111	–104	–5	62	–161	–9	–118
MD-AA2c	ref (53)	151	–170	–145	169	71	–179	–63	138	–64	–87	118	–174
MD-AA3c	ref (53)	142	–162	–90	–129	–86	–102	–57	–40	174	–133	–28	–68
Docking1d	ref (66)	178	–151	105	–14	–91	4	–85	85	–69	59	20	–118
Docking2e	ref (67)	139	175	63	–148	–159	175	–152	147	–54	–53	142	179
Docking2-ce	ref (67)	144	177	115	160	138	–167	–111	135	–62	–76	167	–84

Torsion angles (°) are defined in Figure .

Optimized from X-ray structures with CSD IDs.

Optimized from the structures obtained by MD simulations in water.

Optimized from the structure of Leu-Enk bound to the human δ-opioid receptor 7TM (PDB ID 4N6H, ref (64)) by molecular docking.

Optimized from the structure of Leu-Enk bound to the active δ-opioid receptor crystal structure with the potent opioid agonist peptide KGCHM07 (PDB ID 6PT2, ref (65)) by molecular docking. The structure obtained by optimization with a constraint of four torsion angles (ϕ2, ϕ3, ϕ4, and ϕ5) fixed at X-ray values is denoted by “Docking2-c”, whereas the structure obtained by the fully relaxed optimization is represented by “Docking2”.

Table 8

Type of H-bond, β-Bend, and Thermodynamic Properties of Local Minima for the Zwitterionic Leu-EnK Optimized from the Structures Obtained by X-ray Diffraction, MD Simulation, and Molecular Docking at the DSD-PBEP86-D3BJ/def2-TZVP//M06-2X/6-31+G(d) Level of Theory in Watera

conf.b	H-bond typec	β-bendd	ΔEwe	ΔHwf	ΔGwg
LENKPH11	1 → 4, 4 → 1	βI′23	5.75	5.22	4.55
GEWWAG	4 → 1, 5 → 2	βI23, βI34	4.61	4.25	4.57
BIXNIF10	ext		14.07	13.95	8.05
FABJEX	ext		14.67	14.25	8.75
MD-AA1	OH1 → 3, 2 → 5		4.58	5.11	5.91
MD-AA2	OH1 → 4		7.71	7.79	3.90
MD-AA3	OH1 → 4, 2 → 5, 3 → 5		3.47	4.08	3.60
Docking1	5 → 3		16.59	16.92	16.14
Docking2	ext		12.66	12.12	8.39
Docking2-c	ext		18.93	15.70	13.22

Torsion angles (°) are listed in Table .

The definition of conformations are noted in footnotes b–e of Table .

Each H-bond type n → m stands for the H-bond between the H donor (e.g., the amide H atom for backbone) of the residue n and the H acceptor (e.g., the carbonyl O atom for backbone) of the residue m. In addition, OH1 represents the hydroxyl H atom of the side chain of the Tyr1 residue.

βI23 and βI′23 stand for the type I and I′ β-bend at the Gly2–Gly3 sequence, which are stabilized by the 4 → 1 H-bond. βI34 stands for the type I β-bend at the Gly3–Phe4 sequence, which is stabilized by the 5 → 2 H-bond.

Relative electronic energies in kcal mol–1.

Relative enthalpies in kcal mol–1 at 25 °C.

Relative Gibbs free energies in kcal mol–1 at 25 °C and 1 atm.

Two conformers of LENKPH11 with a single βI′23 bend and GEWWAG with double βI23 and βI34 bends were nearly iso-Gibbs free energy structures with ΔGw = 4.55 and 4.57 kcal mol–1 in water, respectively (Table ). However, two conformers adopted different folded structures stabilized by 1 → 4 and 4 → 1 types of H-bond with the distances of 1.78 and 2.06 Å, respectively, for conformer LENKPH11; 4 → 1 and 5 → 2 types of H-bond with the distances of 2.17 and 2.10 Å, respectively, for conformer GEWWAG (Figure ). Two conformers BIXNIF10 and FABJEX adopted extended structures with ΔGw = 8.05 and 8.75 kcal mol–1 in water, respectively (Table and Figure S4 in the Supporting Information).

Figure 5

Structures of the zwitterionic Leu-EnK obtained by (i) the X-ray diffractions (LENKPH11 from ref (7) and GEWWAG from ref (8)) and (ii) the MD simulations (MD-AA3 from ref (53)), which were optimized at the SMD M06-2X/6-31+G(d) level of theory in water. The relative Gibbs free energy (ΔG in kcal mol–1) at the DSD-PBEP86-D3BJ/def2-TZVP//SMD M06-2X/6-31+G(d) level of theory is shown in parentheses behind each conformation name. For clarity, all nonpolar hydrogen atoms are omitted. All H-bonds are represented by dotted lines with distances in Å.

The conformer MD-AA3 adopted a folded structure with ΔGw = 3.60 kcal mol–1 in water (Table ) stabilized by three H-bonds such as OH1 → 4, 2 → 5, and 3 → 5 types of H-bond with the distances of 1.78, 1.96, and 2.41 Å, respectively, without any β-bends (Figure ). The conformer MD-AA2 exhibited a folded structure with ΔGw = 3.90 kcal mol–1 in water stabilized by a single OH1 → 4 H-bond with the distance of 1.76 Å without any β-bends (Figure S5 in the Supporting Information). The conformer MD-AA1 adopted a folded structure with ΔGw = 5.91 kcal mol–1 in water (Table ) stabilized by two H-bonds such as OH1 → 3 and 2 → 5 types of H-bond with the distances of 1.81 and 1.83 Å, respectively. In particular, it had a βII′34-like motif at the Gly3–Phe4 sequence with a longer distance of 2.80 Å for the 5 → 2 contact due to the more positive shift of torsion angle ϕ3 from the standard value of 60° (Figure S5 in the Supporting Information).

Hence, none of these seven structures optimized in water from X-ray diffractions and MD simulations is consistent with the most preferred conformer zw01 in water. The most preferred conformer zw01 with double β-bends such as βII′23 and βI34 bends is a new folded structure for the zwitterionic Leu-Enk in water, not predicted in earlier works.

To check whether the conformational preference of the zwitterionic Leu-Enk in water is consistent with the structure deduced from NMR experiments, we calculated the ensemble-averaged distances of d(CTyr1α···CPhe4α) and d(CGly2α···CLeu5α) from 43 conformers with ΔGw < 5 kcal mol–1 in Table S6 in the Supporting Information, which can be useful to monitor the formation of β-bend at Gly2–Gly3 and Gly3–Phe4 sequences, respectively. The calculated ensemble-averaged values were 5.74 and 6.13 Å for d(CTyr1α···CPhe4α) and d(CGly2α···CLeu5α), respectively. The latter is consistent with the value of 6.0 Å estimated from the simulated annealing using the 44 distance constraints derived from NOESY spectra.[36] MD simulations starting from this NMR model structure indicated that the zwitterionic Leu-Enk exists as an ensemble of the structures with β-bends at the Gly2 and Phe4 residues and a β-bend at the Gly3 residue in water.[36] This may imply the formation of a double β-bend at Gly2–Gly3 and Gly3–Phe4 sequences and a single β-bend at the Gly2–Gly3 sequence, of which the former is consistent with conformer zw01 and the latter with conformers zw02 and zw03 in the present work.

Structures of Leu-Enk Bound to δOR

There were two studies for the molecular docking of δOR complexed with Leu-Enk using the crystal structures. Sanfelice and Temussi[66] performed the molecular docking of δOR complexed with Leu-Enk using the crystal structure of PDB 4N6H. Recently, Sharma et al. reported the molecular docking of δOR complexed with the Leu-Enk derivative bearing the Tyr1-ψ[(Z)CF=CH]-Gly2 substitution using the crystal structure of PDB 6PT2.[67] These two structures of the zwitterionic Leu-Enk bound to δOR were denoted as Docking1 and Docking2 in the present work, respectively, which were reoptimized at the SMD M06-2X/6-31+G(d) level of theory in water. In addition, Docking2 was reoptimized with the constraint of torsion angles ϕ2, ϕ3, ϕ4, and ϕ5 fixed at the values of the corresponding structure of Leu-Enk bound to δOR at the same level of theory in water, which was denoted as Docking2-c. The optimized torsion angles of these three conformers are listed in Table . The types of H-bond and thermodynamic properties of the corresponding structures calculated at the DSD-dTZ/SMD level of theory in water are shown in Table . The corresponding data calculated at the M062X-dTZ/SMD level of theory are shown in Table S9 in the Supporting Information. The corresponding absolute electronic energies, enthalpies, Gibbs free energies, and single-point energies are listed in Table S10 in the Supporting Information. The three optimized structures are depicted in Figure .

Figure 6

Structures of the zwitterionic Leu-EnK from those bound to δOR optimized at the SMD M06-2X/6-31+G(d) level of theory in water. The relative Gibbs free energy (ΔG in kcal mol–1) at the DSD-PBEP86-D3BJ/def2-TZVP//SMD M06-2X/6-31+G(d) level of theory is shown in parentheses behind each conformation name. For clarity, all nonpolar hydrogen atoms are omitted. H-bond is represented by the dotted line with the distance in Å.

Conformer Docking1 adopted a partially folded structure with ΔGw = 16.14 kcal mol–1 in water (Table ) stabilized by a single 5 → 3 type of H-bond with the distance of 2.28 Å (Figure ). In particular, there was a close contact (about 3.4 Å) between two phenyl rings of Tyr1 and Phe4 residues for conformer Docking1, of which the aromatic ring of the Tyr1 residue was suggested to be essentially stacked with the aromatic ring of the Trp274 in the binding pocket of δOR.[66]

On the other hand, conformer Docking2 exhibited a nearly extended structure in water that did not possess any types of H-bond nor close contacts (Table and Figure ), but it was 3.93 and 7.75 kcal mol–1 more stable in ΔEw and ΔGw, respectively, than conformer Docking1. However, the orientations of the side chain of the Phe4 residue and the carboxylate group of the Leu5 of the fully optimized conformer Docking2 were quite different from those of conformer Docking2 bound to δOR.[67] So, we located another conformer Docking2-c by the constrained optimization with four torsion angles ϕ2–ϕ5 fixed at the values of conformer Docking2 bound to δOR. The values of ΔEw and ΔGw for conformer Docking2-c were 6.27 and 4.83 kcal mol–1 higher than those of the fully optimized conformer Docking2, respectively.

According to the refined model of the zwitterionic Leu-Enk bound to δOR (PDB ID 5C1M),[67] there were favorable interactions of the charged groups NH3+(Tyr1), C=O(Leu5), and C–O–(Leu5) of Leu-Enk with the side chains of Asp128, Lys108, and Arg291 residues of δOR, respectively, as well as a favorable CAryl···H interaction of the Phe4 residue of Leu-Enk with Trp284 of δOR instead of the π···π interaction. Hence, these favorable interactions would induce the conformation of the zwitterionic Leu-Enk (i.e., the folded structure of conformer zw01 with a βII′23 bend at the Gly2–Gly3 sequence and a βI34 bend at the Gly3–Phe4 sequence) isolated in water into the “bioactive conformation” like as conformer Docking2-c when binding to δOR.

Conclusions

The conformational preferences of Leu-Enk were explored by the conformational search and DFT calculations. First, the 738,234 initial structures of the neutral Leu-Enk were generated by a combination of low-energy conformers of each residue and optimized using the ECEPP3 force field in the gas phase, resulting in 2,633 feasible conformations with the relative energy <5 kcal mol–1. Then, these structures were reoptimized at the HF/3-21G(d) level of theory and followed by the optimizations at the M06-2X level of theory with 6-31G(d) and 6-31+G(d) basis functions. We finally identified the 139 structures with the relative energy <10 kcal mol–1 in the gas phase, from which the structures of the corresponding zwitterionic Leu-Enk were generated and reoptimized at the SMD M06-2X/6-31+G(d) level of theory in water. The conformational preferences of the neutral and zwitterionic Leu-Enk were analyzed using Gibbs free energies corrected by single-point energies calculated at the double-hybrid DSD-PBEP86-D3BJ/def2-TZVP level of theory in the gas phase and in water.

In the gas phase, the neutral Leu-Enk dominantly adopted a folded structure stabilized by a side chain-to-backbone (OH1 → 3) H-bond and two backbone-to-backbone (4 → 2 and OH5 → 1) H-bonds with a βII′-bend-like motif at the Gly3–Phe4 sequence and a close contact between the side chains of Phe4 and Leu5. On the other hand, the zwitterionic Leu-Enk exhibited a folded structure in water with a head-to-tail (1 → 5) H-bond between zwitterionic N–H+(Tyr1) and C=O–(Leu5) species and two backbone-to-backbone (4 → 1 and 5 → 2) H-bonds between N–H(Phe4) and C=O(Tyr1) and between N–H(Leu5) and C=O(Gly2). In particular, the two H-bonds of the latter induced the formation of double β-bends such as a βII′ bend at the Gly2–Gly3 sequence and a βI bend at the Gly3–Phe4 sequence, respectively. Hence, the preferred conformations of the neutral and zwitterionic Leu-Enk are new folded structures not predicted by earlier computational studies. However, the calculated ensemble-averaged distance between CGly2α and CLeu5α of the zwitterionic Leu-Enk in water is consistent with the value estimated from the simulated annealing using the distance constraints derived from NOESY spectra in water.

According to the refined model of the zwitterionic Leu-Enk bound to δOR, there were favorable interactions of the charged groups NH3+(Tyr1), C=O(Leu5), and C–O–(Leu5) of Leu-Enk with the side chains of Asp128, Lys108, and Arg291 residues of δOR, respectively, as well as a favorable CAryl···H interaction of the Phe4 residue of Leu-Enk with Trp284 of δOR. Hence, these favorable interactions would induce the structure of the zwitterionic Leu-Enk with double β-bends isolated in water into the “bioactive conformation” like an extended structure when binding to δOR.