Interaction of Cysteine with Li+ and LiF in the Presence of (H2O) n (n = 0-6) Clusters.

The interaction between cysteine with Li+ and LiF in the microcosmic water environment was investigated to elucidate how ions interact with amino acids and the cation-anion correlation effect involved. The structures of Cys·Li+(H2O) n and Cys·LiF(H2O) n (n = 0-6) were characterized using ab initio calculations. Our studies show that the water preferentially interacts with Li+/LiF. In Cys·Li+(H2O)0-6, Li+ interacts with amino nitrogen, carbonyl oxygen, and hydrophobic sulfur of Cys to form a tridentate mode, whereas in Cys·LiF(H2O) n , Li+ and F- work in cooperation and interact with carbonyl oxygen and hydroxyl hydrogen of Cys to form a bidentate type. The neutral and zwitterionic forms are essentially isoenergetic when the water number reaches three in the presence of Li+, whereas this occurs at four water molecules in the presence of LiF. Further research revealed that the interaction between Li+/LiF and Cys was mainly electrostatic, followed by dispersion, and the weakest interaction occurs at the transition from the neutral form to zwitterionic form. Natural population analysis charge analyses show that for Cys·Li+(H2O) n , the positive charge is mostly concentrated on Li+ except for the system containing three water molecules. For Cys·LiF(H2O) n , the positive charge is centered on the LiF unit in the range n = 0-6, and at n = 5, electron transfer from Cys to water occurs. Our study shows that the contribution of anions in zwitterionic state stabilization should be addressed more generally along with cations.

Introduction

Amino acids are vital to life because they are the building blocks of peptides and proteins, which contain both a carboxyl (COOH) group and an amino (NH2) group. The neutral form of the amino acid is the most stable one in the gas phase. However, in the aqueous environment, the molecule may undergo an intramolecular proton transfer, leading to the formation of the zwitterionic form containing COO– and NH3+ groups.[1] Stepwise hydration has been widely studied by gradually adding water molecules to an isolated amino acid.[2−6] The quantity of water molecules required to stabilize the zwitterionic form is a subject of debate, which has been discussed in many investigations. For example, Kayi and co-workers investigated the low-energy conformers of the glycine isomer by computing the relative energies of many different structures in the presence of 1–10 water molecules and found that the glycine switch to the zwitterion being the more stable form occurs when there are eight or nine water molecules;[7] however, Gochhayat et al. observed that glycine required six water molecules to transform into the zwitterion.[8] Sun and co-workers studied the structural elucidation, stability, properties, and proton transfer processes of neutral and zwitterionic Gly(H2O) (n = 1–6) complexes and reported that glycine could be fully solvated by five discrete water molecules,[9] which is inconsistent with Aikens’s results that eight molecules cannot completely solvate glycine.[10]

Additionally, the conformational behavior and function of amino acids are often influenced by the presence of ions.[11−13] The interactions of amino acids and monovalent cations are critical in a variety of chemical and biological processes, such as the osmotic equilibrium of cells, the electrical excitability of nerves and muscles, and the active transport of glucides, which typically involve static electricity, hydrogen bonds (HBs), van der Waals forces, or dispersion interactions. Therefore, further investigations on how ions interact with biological molecules and how these interactions affect, adjust, or control the functions of biological molecules are essential for a better understanding of the role and effect of ions in biological systems.[14] Cation–amino acid clusters have been intensively studied in solution and in the gas phase over recent years to provide insights into ion–amino acid interactions.[14−25] However, the study of anionic species is substantially less intensive, despite the fact that proteins and enzymes usually are negatively charged.[26−29] The investigation[29] of the complexation of the halide ion to the gas-phase amino acids reveals that the zwitterionic and canonical minima are very close in energy for arginine·Br–, whereas the canonical form is significantly lower for arginine·Cl–. It is noteworthy that the cation–anion correlation effect profoundly affects the ion–amino acid interactions and needs to be taken into consideration.

Cysteine (Cys) is an important amino acid carrying an amino (NH2) group and a carboxylic acid (COOH) group, which can be the donor or acceptor of HBs. Besides, the thiol (SH) group of Cys can also donate and accept protons.[30] The interaction of Cys and alkali metal cation complexes has already been investigated both theoretically and experimentally.[13,14,22,31] In the gas phase, zwitterionic forms of Cys (Z-Cys) are reported to be neutralized during interaction with alkali and alkaline earth metal cations by Shankar et al.,[14] who determined that the structure of cysteine is not modified upon metal ion substitution, but the metal ion binding site changes noticeably. Metal ion coordination bond energies decrease as we move down the alkali metal group, and they appear to be inversely related to ion size and charge density.[13] Only the tridentate charge solvated structure was observed for Li+ and Na+, whereas for bigger ions, Rb+ and Cs+, a combination of conformers was discovered, including the bidentate carboxylic acid bound form.[31]

Although extensive studies on the complexes formed by the combination of various amino acids and alkali metal cations in the gas phase and solution have been conducted,[5,11−16,19−23,31−33] there are also a few studies on the interaction of anions and amino acids.[26−29] However, to the best of our knowledge, how anions affect the interaction between cations and amino acids and the cation–anion correlation effect have not been addressed in detail. Particularly, a systematic study on the interaction of Li+/LiF and distinct Cys coordination modes in the presence of zero–six water molecules has not been realized yet, and no research has focused on the cooperation of the metal cation and halide anion. In this study, we carried out studies on the interactions of Cys–Li+/LiF–water and mainly discuss the structures, energies, and properties of Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6), and we also investigated the number of water molecules needed to stabilize Z-Cys under the influence of Li+ or LiF, as well as the effect of the cation(Li+)–anion(F–) cooperative work on them.

Computational Details

The Gaussian 16 program package was used for all of the geometry optimization and frequency calculations of Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6), and the visualization of the optimized structure is completed using GaussView.[34]

There are a lot of low-lying conformations in the systems we are investigating. The initial structures of the small clusters such as Cys·Li+(H2O)0–1 and Cys·LiF(H2O)0–1 were built by altering the locations of the ligand (Li+ and LiF) or the water molecules. The results of the smaller ones [Cys·Li+(H2O)0–1 and Cys·LiF(H2O)0–1] served as a guide to generate large clusters by placing water molecules in various positions. In addition, molecular dynamics (MD) simulations based on the semiempirical tight-binding method GFN-xTB were used to search for possible stable structures in each system’s conformation space by taking different snapshots from trajectories to ensure that all possible structures were identified.[35]

In order to consider the effect of dispersion in structures and estimations of the physicochemical properties, calculations have been performed using Grimme’s dispersion treatment with the original D3 damping function.[36] At the B3LYP-D3 (BJ) level, the geometries of Cys·Li+(H2O) and Cys·LiF(H2O) clusters were optimized by using the 6-311G(d, p) basis set, the optimized structures were checked by calculating the harmonic vibration frequency at this level to verify the true minimum, and then, the single point energy was calculated using the 6-311+G(2d, p) basis set. To evaluate the uncertainty of the B3LYP-D3(BJ) method, the single point energies of Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–3) clusters were further calculated at the CCSD(T)/may-cc-PVTZ level, which gave consistent results in terms of energy ordering of the isomer.

The lowest-energy clusters of Cys·Li+(H2O) and Cys·LiF(H2O) were analyzed using the noncovalent interaction method (NCl),[37] and the interaction energy calculations, energy decomposition (EDA),[38] atoms in molecules (AIM) analysis,[39] and electrostatic potential (ESP) analysis[40] were carried out using Multiwfn software[41] to gain insights into the interaction properties of Cys, the ligand (Li+ and LiF), and water. Additionally, natural population analysis (NPA)[42] and charge model 5 (CM5)[43] methods are selected to analyze the atomic charge of clusters with the lowest energy of Cys·Li+(H2O) and Cys·LiF(H2O), and the results obtained using the NPA method are discussed in detail.

Results and Discussion

Analysis of the Structure

Structures of l-Cysteine

Due to the existence of many possible intramolecular HBs and single bond rotators, cysteine has a large number of conformations. The neutral conformation of Cys (N-Cys) has been calculated by many groups. We reoptimized these conformation isomers at the B3LYP-D3 (BJ)/6-311G(d,p) level based on the results reported in the literature[2,14,44,45] and calculated their single point energies at the CCSD(T)/cc-PVTZ level. All configurations of N-Cys and the single point energy values at the two different levels of theory relative to that of the lowest-energy isomer are shown in Figure .

Figure 1

Low-energy isomers of N-Cys. C in gray, N in blue, O in red, H in white, and S in yellow are used for all the structures. All of the structures are real minima based on the frequency calculation. Relative energies (eV) including zero-point energy obtained using B3LYP-D3 (BJ)/6-311+G(2d,p) and CCSD(T)/cc-PVTZ (in parentheses) methods are also shown.

The energy sequencing of configurations obtained at the B3LYP-D3(BJ) level is completely consistent with the results determined by Bachrach[2] and his colleagues and is basically consistent with our results obtained at the CCSD(T) level, except for a few configurations with high energy. The lowest energy conformer determined by both B3LYP-D3(BJ) and CCSD(T) methods is N(I) (see Figure ), in which hydroxyl (OH) interacts with the amine lone pair and SH points toward carbonyl to form HBs.

The zwitterionic amino acids cannot exist stably in the gas phase, and the solvent therapy is required. In the present work, the geometry was calculated by introducing the SMD implicit solvent model, which considers both polar and nonpolar parts of the solvent effect. Based on the results of the geometry optimized at the B3LYP/6-311G(d,p) level and the single point energy at the CCSD(t)/cc-PVTZ level, six different configurations of the zwitterion using SMD were obtained (Figure ), with the lowest-energy Z(I) being 0.107 eV more stable than N(I) using the SMD model, which was identified by proton transfer from the COOH group in N(I) to the NH2 group. From Table S1, in the case of using the SMD model, all Z-Cys units are more stable than N-Cys, which also accords with the conclusion that zwitterionic amino acids are more stable than neutral molecules in the water environment.

Figure 2

Low-energy isomers of Z-Cys. C in gray, N in blue, O in red, H in white, and S in yellow are used for all the structures. All of the structures are real minima based on the frequency calculation. Relative energies (eV) including zero-point energy obtained using B3LYP-D3 (BJ)/6-311+G(2d,p) and CCSD(T)/cc-PVTZ (in parentheses) methods using the SMD model are also shown.

Structures of Cys·Li+(H2O) (n = 0–6)

The optimized low-energy structures of Cys·Li+(H2O) (n = 0–6) are shown in Figure (more conformations are available in Figure S1). For Cys·Li+, the most stable configuration is N0-A, in which the most preferred position for the interaction of Li+ with N-Cys is in a tridentate manner with amino nitrogen, carbonyl oxygen, and thiol sulfur atoms, consistent with the previous reports.[13,14] Moreover, the O–Li+ coordination distance is found to be shorter than that of N–Li+ and S–Li+, indicating that the interaction of carbonyl oxygen with Li+ is stronger than that of N–Li+ and S–Li+. For the second lowest-energy conformation N0-B, the direction of the H-S is the opposite to that in N0-A, away from the COOH group and only 0.033 eV higher than that of N0-A. The third lowest-energy conformer is the zwitterionic form, labeled as Z-0A, which is 0.311 eV higher than N-0A. In Z0-A, Li+ interacts with both oxygens of the COO– group in a bidentate manner.

Figure 3

Three lowest-energy conformers of Cys·Li+(H2O) (n = 0–6). The growth pattern of solvated clusters with n waters from the (n – 1) clusters is indicated with the nth water circled. C in gray, N in blue, O in red, H in white, and S in yellow are used for all the structures. All of the structures are real minima based on the frequency calculation. Relative energies (eV) including zero-point energy obtained using B3LYP-D3 (BJ)/6-311+G(2d,p) and CCSD(T)/may-cc-PVTZ (n = 0–3, in parentheses) methods are also shown. More isomers can be found in the Supporting Information.

For Cys·Li+(H2O)1, the most stable isomer is N1-A, in which the water molecule binds directly to Li+ to form tetracoordinate metal ions. The second and third lowest-energy structures are zwitterionic forms, designated as Z1-A and Z1-B, respectively. In Z1-A, Li+ sits between Z-Cys and the water molecule, while interacting with the COO– group and water molecule, where the configuration of Z-Cys is Z(I). Z1-B has a Z-Cys configuration of Z(IV), which is different from that of Z1-A, but the manner of interaction of water and Li+ with Z-Cys is identical to that in Z1-A. The relative energies of Z1-A and Z1-B with respect to that of N1-A are 0.243 and 0.306 eV, respectively.

For Cys·Li+(H2O)2, N2-A has the lowest energy, the newly added water molecule breaks the direct interaction between Li+ and the SH group and builds a bridge between them, and Li+ coordinates with amino nitrogen, carbonyl oxygen, and two water molecules. Z2-A and Z2-B are generated by adding an additional water molecule in a similar way to that for Z1-A and Z1-B respectively, in which the second water molecule interacts with Li+ and the COO– group via oxygen and hydrogen atoms, respectively. They are higher in energy by 0.013 and 0.077 eV than N2-A, respectively.

For Cys·Li+(H2O)3 complexes, the most stable isomer N3-A has one water bridging Li+ and the OH group and another water bridging Li+ and the SH group, and Li+ is tetracoordinated by binding to the N-Cys carbonyl and three waters. The next two low-lying energy conformers are both zwitterionic types, which are essentially degenerated with N3-A (just about 0.009 and 0.01 eV above N3-A). Z3-A is generated by attaching Z2-A to the third water molecule that attacks Li+ and receives a proton of the SH group. N3-B is derived from N2-A, in which the newly added water molecule interacts with the carboxyl group by forming an O···H–O bond.

Interestingly, the first three low-energy isomers are all zwitterionic forms and can be thought of as evolving from Z3-A as the water number is increased to four. In Z4-A, the fourth water molecule is used as a bridging molecule between the COO– group and the adjacent water. A cyclic interaction mode is formed between Li+ and Z-Cys through water molecules. In Z4-B, the fourth water molecule acts as a bridge between the SH group and two water molecules interacting with Li+. The newly added water molecules in Z4-C bridge the two water molecules interacting with Li+ and Z-Cys, respectively. Z-4B and Z-4C lie 0.008 and 0.044 eV above Z4-A, respectively.

For Cys·Li+(H2O)5, the most stable isomer Z5-A is yielded from Z4-A, in which the newly added water molecule serves as either a proton donor or a proton acceptor and interacts with the three water molecules interacting with Li+ to form three HBs. Z5-B is developed from Z4-B, and the fifth water molecule bridges the two water molecules that interact with SH and COO– groups, respectively, to form a cyclic HB network. Z5-C can be seen as derived from Z4-C, and the fifth water molecule interacts with the NH3+ group. Z5-B and Z5-C have energies that are 0.018 and 0.043 eV greater than that of Z5-A, respectively.

For Cys·Li+(H2O)6, the most stable isomer Z6-A can be considered as evolved from Z5-A, where the newly added water molecule bridges NH3+ and SH groups of Z-Cys, and generates a cubic type structure. In both Z6-B and Z6-C, the sixth water molecule, located in a second solvent shell, binds to water molecules interacting with Li+ or directly interacting with Z-Cys and links with the existing cyclic HB networks to form a cubic HB network structure. Z6-B and Z6-C are about 0.028 and 0.030 eV more energetic than Z6-A, respectively. In addition, from the structural evolution of Cys·Li+(H2O) (n = 0–6), we found that the maximum coordination number of Li+ is four.

Structures of Cys·LiF(H2O) (n = 0–6)

The three lowest-energy structures of Cys·LiF(H2O) (n = 0–6) are shown in Figure (more conformations are available in Figure S2 in the Supporting Information). For Cys·LiF complexes with ligand LiF instead of Li+, the first three low-energy structures are all neutrals. In N-0A′, the anion F– breaks the tridentate configuration that existed in N-0A of N-Cys·Li+ complexes and forms a new bidentate configuration between LiF and COOH groups of N-Cys, in which Li+ interacts with the nearest O atom and F– interacts with the H atom of the OH group. The Li-F bond length is 1.68 Å. The distance between Li+ and carbonyl oxygen is 1.87 Å and that between F– and the H atom of the OH group is 1.31 Å. N-0B′ and N-0C′ are generated from N(IV) and N(II), respectively, and both of them interact with LiF in a similar way to that in N-0A′. From Figure , for the zwitterionic form, the lowest energy of Z-Cys·LiF complexes is higher than that of N-Cys·LiF by 0.077 eV, and more optimized zwitterionic structures can be found in Figure S2.

Figure 4

Three lowest-energy conformers of Cys·LiF(H2O) (n = 0–6). The growth pattern of solvated clusters with n waters from the (n – 1) clusters is indicated with the nth water circled. C in gray, N in blue, O in red, H in white, and S in yellow are used for all the structures. All of the structures are real minima based on the frequency calculation. Relative energies (eV) including zero-point energy obtained using B3LYP-D3 (BJ)/6-311+G(2d,p) and CCSD(T)/may-cc-PVTZ (n = 0–3, in parentheses) methods are also shown. More isomers can be found in the Supporting Information.

Figure 6

Lowest-energy structures of N-Cys·Li+/LiF(H2O) (n = 0–6) and Z-Cys·Li+/LiF (H2O) (n = 0–6) as well as their EZ-EN values. From left to right are N-Cys·Li+(H2O), Z-Cys·Li+(H2O), N-Cys·LiF(H2O), and Z-Cys·LiF(H2O) (n = 0–6), respectively.

For Cys·LiF(H2O)1, the Cys composition in the most stable N-1A′ is the N(VII) conformation, where the first water molecule interacts with Li+ via the O atom and one H atom points toward the SH group. N-1B′ is derived from N-0A′, in which LiF is sandwiched between the additional water molecule and Cys. Sorted by energy, the third ranked structure is zwitterionic form Z-1A′, in which Li+ and F– atoms interact with COO– and NH3+ groups of Z-Cys, respectively, to form a new bidentate configuration, and the added water molecule interacts with Li+ via the O atom and provides a proton for another oxygen of the COO– group. The Li-F bond lengths in the three isomers are 1.79, 1.77, and 1.77 Å, respectively. N-1B′ and Z-1A′ are higher in energy by 0.084 and 0.121 eV than N-1A′, respectively.

For Cys·LiF(H2O)2, the most stable isomer is N-2A′, where the second water molecule away from Cys interacts with LiF via O–Li+ and H–F– interactions. N-2B′ is derived from N-1A′ with an energy 0.023 eV higher than that of N-2A′, in which the second water molecule provides a proton to the F– atom and forms an HB with the first water molecule. The structure ranked third in terms of energy is zwitterionic type Z-2A′, in which Li+ and the two oxygen atoms of the COO– group are bridged by two water molecules. The relative energy difference between Z-2A′ and the most stable structure N-2A′ is 0.053 eV. The Li–F bond lengths of these three clusters are 1.87, 1.84, and 1.75 Å, respectively, and compared with N-1A′, the Li–F bond length in N-2A′ and N-2B′ increased by 0.08 and 0.05 Å, respectively, but that in Z-2A′ was shortened by 0.02 Å compared with Z-1A′.

For Cys·LiF(H2O)3, the most stable conformer N-3A′ is from N-2B′, where the third water molecule interacts with Li+ via O and one H atom points toward the SH group. N-3B′ is an outgrowth of N-2A′, and the third water molecule forms HBs by interacting with the F atom and with the adjacent water molecule bridging the SH group and Li atom. The third lowest energy conformation is zwitterionic type Z-3A′ which is derived from Z-2A′, and the third water molecule interacts with the COO– group and two adjacent water molecules via HB interactions. The energies of N-3B′ and Z-3A′ are 0.031 and 0.049 eV higher than that of N3-A′. The Li–F bond lengths of these three clusters are 1.89, 1.93, and 1.75 Å, respectively.

For Cys·LiF(H2O)4, the most stable isomer N-4A′ can be considered as obtained from N-3A′ with the fourth water molecule interacting with F– and forming a HB with the adjacent water. The second lowest cluster is zwitterionic Z-4A′, which is obtained from Z-3A′ and has a 0.003 eV higher energy than that of the neutral N-4A′, in which the fourth water molecule interacts with the SH group, F–, and the nearest water molecule through HBs. N-4B′ is also derived from N-3A′, in which the fourth water molecule acts as a donor to provide a proton to the SH group and F–, respectively. The energy gap between N-4B′ and the most stable cluster is 0.010 eV. The Li–F bond lengths of these three clusters are 1.90, 1.78, and 1.89 Å, respectively.

For Cys·LiF(H2O)5, the most stable isomer changes from neutral to the zwitterionic type, and the lowest-energy Z-5A′ can be considered as generated from Z-4A′, in which the fifth water molecule interacts with COO– and NH3+ groups. The second lowest-energy cluster Z-5B′ is similar to Z-5A′, except that the SH group in Z-5A′ interacts with the adjacent water, while the SH group is free in Z-5B′. In N-5A′, the fifth water molecule bridges Li+ and another water molecule. The energies of Z-5B′ and N-5A′ are higher than that of Z-5A′ by 0.048 and 0.068 eV respectively. The Li–F bond lengths of these three clusters are 1.77, 1.77, and 1.87 Å, respectively.

For Cys·LiF(H2O)6, the most stable isomer is Z-6A′, obtained from Z-5A′, in which the sixth water molecule interacts with the COO– group and the water molecule that interacts with the SH group as a HB donor, respectively, and works cooperatively with the adjacent cyclic HBs to form a cubic HB network structure. The second lowest-energy isomer N-6A′ is obtained from N-5A′, in which the sixth water molecule forms a cyclic HB structure by interacting with two water molecules that interact with F–, and meanwhile, it interacts with the nearest cyclic HB to form a cubic HB network structure. Z-6B′ is similar to Z-6A′, except that the sixth water molecule does not directly interact with the COO– group. The energies of N-6A′ and Z-6B′ are 0.008 and 0.024 eV higher than that of the lowest-lying isomer Z-6A′, respectively. The Li-F bond lengths of these three clusters are 1.77, 1.90, and 1.77 Å, respectively. During the successive addition of 1–6 water molecules, the bond lengths of Li-F with low-energy structures were not significantly extended, indicating that Li+ and F– always exist in the form of contact ion pairs.

Effect of Li+ or LiF and Water Coordination on the Stability of Zwitterionic Cysteine

We examined the number of water molecules needed to stabilize Z-Cys in the presence of Li+ or LiF by adding explicit water molecules one by one. Figure shows the energy difference between the lowest-energy structures of the neutral complexes Cys·L(H2O)0–6 (L = Li+/LiF) and the corresponding zwitterionic form isomers (detailed structural information is presented in Figure ). For Cys·Li+(H2O)0–6, Figure shows that the value of EZ-EN decreases obviously when the first two water molecules are involved, indicating that the stability of Z-Cys increased, but two water molecules are not sufficient to stabilize Z-Cys as N-Cys is energetically more favorable. When the water number is increased to three, the lowest zwitterionic tautomer and the corresponding neutral are of nearly equal energy (with an energy difference of <0.01 eV) considering the inherent error in density functional theory calculations, indicating that three water molecules can basically convert the structure from neutral to zwitterionic. For Cys·LiF(H2O)0–6, the energy of neutral structures approaches that of the zwitterionic isomer when the water number reaches four, indicating that the zwitterionic structure is more difficult to form than Cys·Li+(H2O)0–6 because of the effect of anion F–. According to the Grotthuss type[46] mechanism, the formation of Z-Cys with NH3+ and COO– groups is due to the diffusion of protons through HBs in water, so we speculate that anion F– hinders the proton mobility.

Figure 5

Variation of the energy difference between the lowest-energy zwitterionic and neutral structures with the number of water molecules for Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6) complexes at the B3LYP-D3 (BJ)/6-311+G(2d,p) level of theory.

The thermochemical study was also carried out and yielded thermal corrections. Relative theoretical 0 K enthalpies and 298 K free energies were calculated at the B3LYP-D3(BJ)/6-311+G(2d,p)//B3LYP-D3(BJ)/6-311G(d,p) levels of theory, where the vibrational frequencies calculated were scaled by 0.9888.[47] As shown in Tables and 2, with regard to the free energy of typical low lying isomers at 298 K of Cys·Li+(H2O), the neutral conformer is the lowest-energy isomer for n = 0–2. When the water number approaches three, zwitterionic isomer Z-3B is more favorable than neutral N-3A by ∼0.7 kJ mol–1, indicating that the structure of the Cys·Li+(H2O) complex can coexist as neutral and zwitterionic isomers. However, for Cys·LiF(H2O)0–6, the system containing zero–four water molecules prefers to form neutral type conformers, and the structure of the complex begins to be characterized by zwitterions until the water number is boosted to five, which is consistent with the above results of Ez -EN analysis. It seems that with the cooperation of the cation and anion, more water molecules are required to stabilize zwitterions than in systems containing only metal ions. Therefore, the contribution of anions in zwitterionic state stabilization should be addressed more generally along with cations.

Table 1

Relative Theoretical 0 K Enthalpies and 298 K Free Energies (kJ/mol) of Cys·Li+(H2O) (n = 0–6) Using the B3LYP-D3(BJ)/6-311+G(2d,p)//B3LYP-D3(BJ)/6-311G(d,p) Method

complex	isomers	0 K	298 K
Cys·Li+	N-0A	0.0	0.0
	N-0B	3.0	2.7
	Z-0A	29.8	28.1
Cys·Li+ (H2O)1	N-1A	0.0	0.0
	Z-1A	23.5	20.8
	Z-1B	28.5	25.0
Cys·Li+ (H2O)2	N-2A	0.0	0.0
	Z-2A	2.9	1.0
	Z-2B	7.8	4.7
Cys·Li+ (H2O)3	N-3A	0.0	0.0
	Z-3A	5.9	9.1
	Z-3B	0.3	–0.7
Cys·Li+ (H2O)4	Z-4A	0.0	0.0
	Z-4B	0.3	0.2
	Z-4C	3.5	2.9
Cys·Li+ (H2O)5	Z-5A	0.0	0.0
	Z-5B	0.1	–0.8
	Z-5C	–1.7	–5.3
Cys·Li+ (H2O)6	Z-6A	0.0	0.0
	Z-6B	1.9	0.9
	Z-6C	2.9	–2.5

Table 2

Relative Theoretical 0 K Enthalpies and 298 K Free Energies (kJ/mol) of Cys·LiF(H2O) (n = 0–6) Using the B3LYP-D3(BJ)/6-311+G(2d,p)//B3LYP-D3(BJ)/6-311G(d,p) Method

complex	isomers	0 K	298 K
Cys·ys	N-0A′	0.0	0.0
	N-0B′	2.0	1.9
	N-0C′	3.3	2.7
Cys·LiF (H2O)1	N-1A′	0.0	0.0
	N-1B′	11.9	8.7
	Z-1A′	16.2	15.8
Cys·LiF (H2O)2	N-2A′	0.0	0.0
	N-2B′	8.8	6.5
	Z-2A′	13.2	14.3
Cys·LiF (H2O)3	N-3A′	0.0	0.0
	Z-3A′	10.6	14.2
	N-3C′	3.9	4.6
Cys·LiF (H2O)4	N-4A′	0.0	0.0
	Z-4A′	8.2	12.4
	N-4B′	0.5	0.6
Cys·LiF (H2O)5	Z-5A′	0.0	0.0
	Z-5B′	6.0	6.1
	N-5B′	2.0	0.6
Cys·LiF (H2O)6	Z-6A′	0.0	0.0
	Z-6B′	0.3	–0.8
	N-6A′	–1.7	–1.2

Interaction Discussion

Reduced Density Gradient Analyses

In order to examine the interactions between Cys, ions and water in Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6) complexes, we performed a series of weak interaction analyses based on their most stable isomers.

The noncovalent interaction index based upon the reduced density gradient (RDG) analysis is a powerful tool to explore noncovalent interaction intuitively and supplies more evidence of noncovalent interaction. The location and type of interatomic interaction can be recognized by observing the isosurface diagram, in which the deeper blue region indicates the stronger electrostatic attractive interaction, the red area represents the steric effect, and the green one represents the van der Waals effect. The corresponding scatter diagrams that represent the electron density multiplied by the sign of the second Hessian eigenvalue [symbol(λ2)ρ] are also created in conjunction with the RDG diagram, from which the isosurface is expressed quantitatively.

In this work, RDG analysis is realized using Multiwfn program, and the display is achieved using the visualization program VMD. The RDG analysis results of the Cys·Li+(H2O)0–6 and Cys·LiF(H2O)0–6 complexes are shown in Figures and 8. For Cys·Li+ (H2O)0–3, the isosurface between Li+ and carbonyl (CO) and that between the NH2 group and water molecules are shown in blue, and the corresponding spikes of these positions are located between −0.03 and −0.02 au on the scatter diagram, indicating that there are strong electrostatic attractions. The isosurface between the SH group and Li+ and that between the SH group and water molecules are shown in green, which corresponds to the spikes near −0.01 au on the scatter diagram, showing van der Waals interaction. For Cys·Li+(H2O)4–6, the regions between water molecules and COO– or NH3+ groups of Z-Cys are shown in blue color, which corresponds to a spike near −0.035 au on the scatter diagram, indicating the strong hydrogen bonding.

Figure 7

RDG surface plot and scatter plot of the lowest-energy isomers of Cys·Li+(H2O) (n = 0–6). The surfaces are colored on a blue–green–red scale according to values of sign(λ2)ρ, ranging from −0.035 to 0.02 a u. Blue indicates strong attractive interactions, and red indicates steric clash.

Figure 8

RDG surface plot and scatter plot of the lowest-energy isomers of Cys·LiF(H2O) (n = 0–6) isomers. The surfaces are colored on a blue–green–red scale according to values of sign(λ2)ρ, ranging from −0.035 to 0.02 a u. Blue indicates strong attractive interactions, and red indicates steric clash.

From Figure , for Cys·LiF, it can be clearly seen that the interaction between the OH group and F– is obviously a strong HB based on the blue color of the isosurface between them. As for Cys·LiF(H2O)0–4, the abscissa of a spike here exceeds the maximum negative value of the scatter plot of −0.05, indicating that there was a very strong electronic attractive interaction between F– and Cys carboxyl. In addition, the isosurface between Li+ and CO of N-Cys and that between Li+ and water molecules are shown in light blue, and the corresponding spikes of these positions are located between −0.03 and −0.02 au on the scatter diagram, which can also be considered as electronic interaction. The isosurface between water molecules and Cys is shown in green, corresponding to the spike between −0.02 and −0.01 au on the scatter diagram, which is weaker than the interaction between Li+ and Cys or F– and Cys. For Cys·LiF(H2O)5–6, the stability of the Z-Cys configuration is more advantageous. We determine that the interaction between F– and the NH3+ group shows strong hydrogen bonding based on the blue region of RDG and the corresponding spike being located at −0.04 of the scatter diagram.

Interaction Energy between Cys and Li+/LiF

To get a better understanding of the interaction between Cys and Li+ or LiF in water clusters, the interaction energy between Li+ and Cys and that between LiF and Cys in the lowest-energy complexes were calculated accurately through the counter-poise method[48] with correction for the basis set superposition error. The interaction energy between A and B in the presence of C, ΔEC(A–B), was calculated as ΔE(A–BC) – ΔE(A–C).[49] In this case, A represents Li+ or LiF, B represents Cys, and C represents H2O clusters.

From Table , for Cys·Li+(H2O) (n = 0–6), ΔEC(A–B) decreases linearly with the increase of water molecules from zero to three, by 44.85, 55.32, and 71.54 kJ/mol, respectively, indicating that the interaction strength decreased when going from n = 0 to n = 3. When the water number is increased to four, the most stable conformer changed to the zwitterionic form, and the ΔEC(A–B) no longer changes appreciably as the number of water molecules is increased.

Table 3

Corrected Interaction Energies (kJ/mol) for Li+ and Cys in the Presence of H2O Clusters in Cys·Li+(H2O) (n = 0–6) Using the B3LYP-D3(BJ)/6-311+G(2d,p) Methoda

	ΔE(A–BC)	ΔE(A–C)	ΔEC(A–B)
N-0A			–333.34
N-1A	–438.32	–149.83	–288.49
N-2A	–508.65	–275.47	–233.17
N-3A	–551.70	–390.07	–161.63
Z-4A	–584.04	–412.38	–171.67
Z-5A	–573.50	–406.60	–166.90
Z-6A	–547.35	–369.20	–178.15

A: Li+, B: Cys, and C: H2O clusters.

When the number of water molecules is five and six, the value difference is 4.77 and −11.25 kJ/mol, respectively, compared with that at n = 4. In addition, based on the results of EDA using the EDA-FF method (Table S2), it was found that the interaction between Li+ and Cys and that between Li+ and water are mainly electrostatic, supplemented by a small amount of dispersion. From Table , for Cys·LiF(H2O) (n = 0–6), ΔEC(A–B) increased from −206.44 to −326.39 kJ/mol when the first water molecule is added, and then, ΔEC(A–B) decreased with the increase of water molecules from one to four, by 77.44, 71.3, and 15.9 kJ/mol respectively, indicating that the interaction between LiF with Cys gradually weakened with the increase of the water number. The interaction energy of ΔEC(A–B) increases when the water number reaches five, indicating that the interaction between LiF with Cys increases when the structure shifts from neutral to zwitterionic. From the results of EDA using the EDA-FF method (Table S3), the interaction between LiF with Cys molecules is dominated by electrostatic force, followed by dispersion, which is similar to that of Cys·Li+(H2O) (n = 0–6).

Table 4

Corrected Interaction Energies (kJ/mol) for LiF and Cys in the Presence of H2O Clusters in Cys·LiF(H2O) (n = 0–6) Using the B3LYP-D3(BJ)/6-311+G(2d,p) Methoda

	ΔE(A–BC)	ΔE(A–C)	ΔEC(A–B)
N-0A′			–206.44
N-1A′	–397.69	–71.30	–326.39
N-2A′	–422.37	–173.43	–248.95
N-3A′	–440.58	–262.92	–177.65
N-4A′	–486.85	–325.10	–161.75
Z-5A′	–412.33	–224.89	–187.44
Z-6A′	–437.98	–244.72	–193.26

A: LiF, B: Cys, and C: H2O clusters.

Topological Analyses

In this paper, based on the AIM theory, the topological analysis of electron density (ρ) of the lowest-energy configurations of Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6) complexes were carried out by using Multiwfn software. Tables and 6 give the electron density (ρBCP) and Laplacian of electron density (∇ρBCP2) values at the bond critical point (BCP) between Li+/LiF and Cys (CO, OH, N, and S), which are closely related to the interaction strength. The 2D diagrams of ρBCP and ∇ρBCP2 are shown in Figures S3–S6, where the blue point represents the critical point of the bond between two atoms. The dotted line and solid line region in the 2D diagram of ∇ρBCP2 represent the regions where the ∇ρBCP2 is negative and positive, respectively.

Table 5

Electron Density (ρBCP) and Laplacian of Electron Density (∇ρBCP2) Values at the BCP between Li+ and Cys (CO, N, and S) in Cys·Li+ (H2O) (n = 0–6) Clusters

isomers	BCP	ρBCP	∇ρBCP2
N-0A	Li+···CO	0.0249	0.2119
	Li+···NH2	0.0240	0.1700
	Li+···SH	0.0107	0.0806
N-1A	Li+···CO	0.0208	0.1765
	Li+···NH2	0.0215	0.1521
	Li+···SH	0.0082	0.0655
N-2A	Li+···CO	0.0212	0.1807
	Li+···NH2	0.0207	0.1445
N-3A	Li+···CO	0.0214	0.1867
Z-4A	Li+···CO	0.0266	0.2365
Z-5A	Li+···CO	0.0284	0.2569
Z-6A	Li+···CO	0.0338	0.2823

Table 6

Electron Density (ρBCP) and Laplacian of Electron Density (∇ρBCP2) Values at the BCP between LiF and Cys (CO, OH, and N) in Cys·LiF(H2O) (n = 0–6) Clusters

isomers	BCP	ρBCP	∇ρBCP2
N-0A′	Li···CO	0.0334	0.2857
	F···OH	0.0526	0.5597
N-1A′	Li···CO	0.0311	0.2576
	F···OH	0.0700	0.4396
N-2A′	Li···CO	0.0242	0.1990
	F···OH	0.0638	0.5986
N-3A′	Li···CO	0.0230	0.1915
	F···OH	0.0476	0.5052
N-4A′	Li···CO	0.0206	0.1747
	F···OH	0.0411	0.4158
Z-5A′	F···NH3+	0.0389	0.3755
Z-6A′	F···NH3+	0.0390	0.3751

From Table , for N-0A, the ρ values at the BCP positions between Li+ and CO and Li+ and N are 0.0249 and 0.0240, respectively, which are greater than that of Li+···SH (0.0107), indicating that the interactions of Li+···CO and Li+···NH2 are stronger than that of Li+···SH. After the first water molecule was added, the ρBCP of these three places decreased. It is noted that ρBCP between Li+···CO increased when more water molecules are involved, whereas the direct interaction of Li+···SH and Li+···NH2 continued to weaken until disappeared. From Figures S3 and S5, for Cys·Li+(H2O) (n = 0–6), we can see that the BCP of the interaction positions of Li+ and Cys are located in the solid line area and ∇ρBCP2 values are in the range of 0.0655–0.2823 au. It is expected that the connection between Li+ and Cys can be considered as a closed-shell interaction, indicating that the interaction is dominated by electrostatic attraction, which agrees well with the results of RDG and interaction energy analyses.

From Table and Figures S4 and S6, for Cys·LiF(H2O) (n = 0–6), the BCPs between LiF and Cys are located in the solid line area of 2D Laplacian of the electron density map and the values of ∇ρBCP2 fall between 0.1915 and 0.5986 au, indicating that the interaction between LiF and Cys is also closed-shell. For Cys·LiF(H2O)0–4, we found that ρBCP of F–···OH is larger than that of Li+···CO, indicating that F–···OH interaction is stronger than that of Li+···CO. This is in line with the above RDG and interaction energy analyses. For Cys·LiF(H2O)5–6, Table shows that ρ values at the BCP position between F– and NH3+ are 0.0389 and 0.0390 and that ∇ρBCP2 values are 0.3755 and 0.3751, respectively, indicating that the interaction of F– with NH3+ is electrostatic attraction and noncovalent. However, no BCP was observed between Li+ and CO, indicating no direct interaction between Li+ and CO.

Charge Calculation

In this paper, NPA and CM5 methods are employed to analyze the atomic charge for the lowest energy of Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6) complexes, and the results obtained using the NPA method are discussed in detail. The charge values of each part are shown in Tables and 8, and the graphic displays are shown in Figures S7 and S8. The changes of charge values with the addition of water molecules in Cys·Li+(H2O)0–6 and Cys·LiF(H2O)0–6 complexes are shown in Figure a,b. The calculated results of CM5 can be found in Tables S6 and S7 for comparison.

Table 7

Calculated NPA Charge for the Lowest-Energy Complexes of Cys·Li+(H2O)0–6 Using the B3LYP-D3(BJ)/6-311+G(2d,p) Method

isomers	Cys	Li+	H2O
N-0A	0.0872	0.9128	0.0000
N-1A	0.0962	0.8766	0.0271
N-2A	0.0849	0.8872	0.0279
N-3A	–0.0923	0.4966	0.5957
Z-4A	0.0643	0.8863	0.0494
Z-5A	0.0697	0.8843	0.0461
Z-6A	0.0276	0.8747	0.0977

Table 8

Calculated NPA Charge for the Lowest-Energy Complexes of Cys·LiF(H2O)0–6 Using the B3LYP-D3(BJ)/6-311+G(2d,p) Method

isomers	Cys	LiF	H2O
N-0A′	–0.1343	0.1342
N-1A′	–0.1977	0.2229	–0.0252
N-2A′	–0.1396	0.1477	–0.0081
N-3A′	–0.0722	0.0858	–0.0137
N-4A′	–0.0444	0.0671	–0.0227
Z-5A′	–0.0048	0.0848	–0.0800
Z-6A′	–0.0095	0.0838	–0.0743

Figure 9

Variation of the NPA charge with the number of water molecules at the B3LYP-D3 (BJ)/6-311+G(2d,p) level of theory. (a) Cys, Li+, and H2O clusters of Cys·Li+(H2O) (n = 0–6) and (b) Cys, LiF, and H2O clusters of Cys·LiF(H2O) (n = 0–6).

For Cys·Li+(H2O)0–3, it can be seen from Table and Figure A, with the interaction between Cys and Li+ in N-0A, the electrons of Cys transfer to Li+, which reduces the positive charge of Li+ to 0.9128. In N-1A, the electrons of water molecules move to Li+ via the interaction between water molecules and Li+, lowering the positive charge of Li+ to 0.8766, which weakens the interaction between Li+ and Cys. As the water molecule number is increased to three, the positive charge of Li+ reduced significantly but that on water increased noticeably, which resulted in the weakest interaction between Li+ and Cys. When the water number reached four, the positive charge on Li+ increased again by about 0.3898 compared to that of the system containing three water molecules, and then, the charge distribution changed slightly in the range n = 4–6. The overall positive charges are mainly concentrated on Li+ except that of N-3A containing three water molecules, which is the transition from the neutral form to zwitterionic form.

For Cys·LiF(H2O)0–6, the positive charge is always localized on the LiF unit from n = 0 to n = 6, whereas the negative charge is mainly localized on the Cys until n = 4, and then, electron transfer from Cys to water starts and negative charge is mainly concentrated on water for n = 5 and 6.

Electrostatic Potential Map

According to Figure , the positive potential surface of Li+ in N-0A has obvious penetration into the van der Waals surface of atoms (N, S, and O) with negative ESP in N-Cys. It shows that there is a strong noncovalent interaction dominated by electrostatic interaction, which is in nice agreement with the above analysis in the present work. When the first water molecule is added, we found that the negative potential surface of oxygen of water and the positive potential surface of Li+ overlap. With the addition of more water molecules, for Cys·Li+(H2O)2–3, the potential surface of water molecules overlaps with the regions of COOH and SH groups of N-Cys that have the opposite potentials. For Cys·Li+(H2O)4–6, the potential surface of water molecules overlaps with the areas with the opposite potential of COO–, SH, and NH3+ groups in Z-Cys.

Figure 10

ESP map on the molecular surface of Cys·Li+(H2O) (n = 0–6) isomers with the isodensity surface value of 0.001 a u. The positive ESP is colored in red, and the negative ESP is colored in blue.

For Cys·LiF(H2O) (n = 0–6), from Figure , in N-0A′, the penetration between the F– and OH groups is very deep, indicating that there is a strong electrostatic interaction between them, which is obviously stronger than that between Li+ and CO groups. For Cys·LiF(H2O)1–4, at the interaction regions between water molecules and LiF, the van der Waals surfaces with the opposite potential penetrate each other. For Cys·LiF(H2O)5–6, only the COO– group in Cys shows negative potential, while all other positions have positive potential. The penetration degree between F– and NH3+ groups is deeper than that of any other interaction regions in the complex, which indicates that there is a strong electrostatic interaction between them.

Figure 11

ESP map on the molecular surface of Cys·LiF(H2O) (n = 0–6) isomers with the isodensity surface value of 0.001 a u. The positive ESP is colored in red, and the negative ESP is colored in blue.

Conclusions

The interaction of cysteine with Li+/LiF as well as the solvation of various complexes by one–six water molecules was investigated in the gas phase at the B3LYP-D3(BJ)/6-311+G(2d,p) and CCSD(T)/may-cc-PVTZ levels of theory, and the following important results have emerged from the detailed calculations.

For the Cys·Li+ complex, Li+ interacts with amino nitrogen, carbonyl oxygen, and thiol sulfur atoms of Cys in a tridentate manner; however, with the addition of F–, a new bidentate configuration is formed between LiF and Cys. For both Cys·Li+ and Cys·LiF, with the increase of the number of water molecules, the lowest-energy configuration has undergone a transition from a cyclic motif to a cubic type structure. The neutral and zwitterionic forms are essentially isoenergetic with three water molecules for Cys·Li+(H2O) (n = 0–6), whereas stabilizing Z-Cys in Cys·LiF(H2O) requires four water molecules.

Through RDG, AIM, d ESP, and EDA analyses, and the calculation of interaction energy, we found that the intramolecular interaction types in Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6) complexes are dominated by electrostatic interaction, followed by dispersion. For Cys·Li+(H2O), the interaction between Li+ and Cys gradually decreases from n = 0 to n = 3, increases slightly at n = 4, and then changes a little up to n = 6. For the system containing LiF, the interaction between LiF and Cys grew dramatically after the addition of the first water molecule, then gradually declined with the addition of the water molecules n = 1–4, and surged again at n = 5–6. For both Cys·Li+(H2O) and Cys·LiF(H2O) (n = 0–6), the weakest interaction between the ligand (Li+/LiF) and Cys occurs at the transition from the neutral form to zwitterionic form (N-3A and N-4A′). Additionally, NPA charge analyses show that for Cys·Li+(H2O)0–6, the positive charge is mainly localized on Li+ except the system with three water, which is the transition from the neutral form to zwitterionic form. For Cys·LiF(H2O)0–6, the positive charge is always concentrated on the LiF unit in the range n = 0–6, with the electron transfer from Cys to water starting at n = 5. The present results provide valuable insights into the synergistic effect of anions and cations with amino acids in a water environment, which might contribute to a further understanding of the mechanism with which ions interact with proteins and ligands in a saline solution of the living system.