Hydrogen Bond Strengthens Acceptor Group: The Curious Case of The C-H...O=C Bond

1. Introduction

Hydrogen bonding is one of the most important forces governing molecular integrity in the physical and biological worlds. It describes the sharing of a hydrogen that is covalently bound to an electronegative Donor atom, with another electronegative Acceptor atom that has a lone pair of electrons: D−H···A.

In the biological world, the acceptor group is often composed of nitrogen or oxygen bound to a carbon. Upon H-bonding, the C−O or C−N bonds are weakened due to the pulling exerted on the electronegative atom by the hydrogen. In the current study, we identified the opposite effect, wherein the acceptor of an H-bond is surprisingly strengthened upon H-bonding. Consequently, by looking into the governing factors of H-bond formation, we aimed to examine the roots of this unexpected and intriguing observation.

Theoretical and experimental findings describe the formation of hydrogen bonds due to several energy components: electrostatic attraction, electron delocalization, charge transfer, dispersive interactions, cooperative effects, steric effects, and secondary interactions [1].

H-bonds in proteins exhibit variations in these energy components, tailored to fulfill the specific physiological role of the molecule [2]. The nature of H-bonds in proteins is further influenced by the molecular makeup of the donor and acceptor, the protein backbone, and the side chain [3,4].

One can classify protein and peptide H-bonds as canonical or noncanonical. Canonical or conventional H-bonds are formed by regular and predictable H-bonding interactions, like the i to $i+4$ H-bonds formed between the C=O and N−H groups in -helices [5]. Noncanonical H-bonds exhibit a significantly larger variation of donor and acceptor groups and a multiplicity of partners. For example, an over-coordinated system entails several donors bonded to a single acceptor, while multifurcation involves a single donor with multiple acceptors. Experimental and computational analyses have shown that some of these noncanonical bonds are highly important for protein structure [6,7].

An important class of noncanonical H-bonds in proteins are those formed between the a C hydrogen and a carbonyl oxygen [8,9,10,11,12,13,14,15]. X-ray and neutron crystallography have confirmed the identity of C−H···O=C hydrogen bonds within the protein environment in the presence of a solvent, other proteins, and ligands [16,17], suggesting an essential role in catalytic activity [18]. These important hydrogen bonding schemes necessitate detailed analyses to provide a quantitative understanding of their role in macromolecules.

The C−H···O=C hydrogen bond in a transmembrane helical dimerization interface, was first identified in glycophorin A [19]. The interface has a GxxxG motif, one of the most prevalent oligomerization factors in transmembrane helices [20,21]. In the GxxxG motif of glycophorin A, the C−Hs of glycines 79 and 83 form H-bonds with the carbonyls of isoleucine 76 and valine 80 in the opposing helix, respectively. Such bonds form the basis of the GxxxG dimerization motif [13].

Previous attempts to measure the energetics of the C−H···O=C H-bonds resulted in varying estimates. Using a mutagenesis cycle analysis, Bowie and coworkers found that such bonds are not stabilizing [22], while an empirical FTIR-based approach indicated a bond enthalpy of 0.88 kcal/mol [23].

Considering the prevalence of the C−H···O=C interaction and its contested contribution to protein stability, we decided to examine the energetics of C−H···O=C hydrogen bonds using a combination of experimental and computational methods. Experimentally, we employed vibrational spectroscopy, an exquisitely sensitive tool to measure H-bond strength. Computationally, we use Molecular Dynamic Simulations and DFT calculations that can yield detailed energetics and vibrational shifts that can be verified with experiments. Specifically, we targeted the frequency change of the C=O acceptor upon H-bond formation, made possible by isotopic labelling to resolve the specific carbonyl group [24,25]. This combined approach enabled us to both measure the strength of this interaction and explain the source of its curious spectroscopic behavior.

2. Results and Discussion

Glycophorin A is the first membrane protein to have its sequence determined [26] and is the archetypical GxxxG containing transmembrane domain dimer [13,19,27]. The GxxxG motif enables close positing of the two helices with minimal steric hindrance due to the small size of the glycine hydrogen side-chain that H-bonds to a carbonyl residue in the other helix. Figure 1 depicts the bond between the C−H of Gly83 and the carbonyl acceptor of Val80. Note that the acceptor is also involved in a canonical H-bond with the backbone amide donor hydrogen Val84. Hence, Val80’s carbonyl group participates in an over-coordinated H-bonding interaction.

We spectroscopically isolated the C−H···O=C bond in Glycophorin A by editing of Val80’s carbonyl group with ¹³C and ¹⁸O. Due to the fact that the amide I mode is composed mainly of the C=O stretch [28], such labeling baseline resolves the isotope-edited mode from the unlabelled amide groups by shifting it more than 60 cm₋1, enabling detailed site-specific analysis of the labeled peak separately from the unlabeled peak. [25]. Accordingly, as shown in Figure 2, the ¹³C=¹⁸O amide I mode of Val80 within glycophorin A is shifted to 1599.4 cm₋1, which is 62 cm₋1 from the main amide-I envelope of the protein. The entire spectra of both peptides, in which all other vibrational modes arising from the lipid and protein can be observed, are presented in Figure S1.

The peptide encompasses the transmembrane domain of wild-type glycophorin A, a strongly dimeric species [27,29,30] in which the C−H···O=C hydrogen bond exists [13,19]. Therefore, to determine the impact of the Gly83–Val80 inter-helical H-bond, we measured the spectrum obtained from a peptide containing the Gly79Leu monomerizing mutation which separates the helices and so the bond does not exist [27].

In the monomeric protein that does not contain the noncanonical H-bond, the isotope-labeled amide I mode resonates at 1596.8 cm₋1, which is lower than the H-bonded dimeric species by 1.7 cm₋1 (Figure 2). These results appear surprising since H-bond formation is expected to reduce the vibrational energy of the acceptor carbonyl and thus shift the isotope-edited carbonyl peak to lower frequencies when it is involved in a hydrogen bond.

We note that previous studies report a blue shifting of an H-bond donor. For example, using computational tools, Schlegel and colleagues have elucidated the electrostatic origin of this shift [31]. However, to the best of our knowledge, we are reporting for the first time the blue shift of the stretching frequency of an H-bond acceptor C=O in proteins.

To understand the root of the surprising frequency shift, we employed DFT-based computational tools that detail the characteristics of such H-bonds [32]. The size of a transmembrane helix system is beyond the current capabilities of detailed quantum calculations. Therefore, two mimetic systems that capture the specific H-bonding interactions were analyzed and compared (see supplementary Figure S4):

Monomer (single H-bond):

The canonical -helical H-bond between the carbonyl of valine 80 and the amide H of glycine 84 was mimicked by two N-methylacetamides.

Dimer (two H-bonds):

The inter-helical H-bonding system contained the canonical H-bond described above, with an additional non-canonical inter-helical hydrogen bond. The C−H of glycine 83 from the opposing helix is bonded to the same carbonyl of valine 80 that is also involved in the canonical hydrogen bond. This system was mimicked by two N-methylacetamides forming the canonical H-bond to which an acetylglycinemethylamide is hydrogen bonded.

To ensure that the calculations replicate the experimental system accurately, we sought to superimpose the coordinates of the atoms on the corresponding groups of the transmembrane protein. However, only the structure of the wild-type, dimeric glycophorin A, has been solved experimentally [19,33]. Moreover, the structures were elucidated in micelles or bicelles but not in a lipid bilayer, the native environment of the protein and the one in which the FTIR spectra were obtained. Therefore, we used molecular dynamics (MD) simulations in hydrated 1,2-dimyristoyl-sn-glycero-3-phosphocholine bilayers to determine the atom positions in the wild-type, dimeric glycophorin A and a monomeric glycophorin A peptide that contains the G79I mutations. In both instances, the experimentally determined structure in bicelles was used as a starting point [33]. The results of the MD simulation can be seen in supporting information, Figure S2, depicting structures of both monomer and dimer species and Figure S3 showing the root mean square deviation for each species.

We extended the utility of the MD simulations beyond creating mimetic systems for density functional theory (DFT) calculations, applying them also to geometric analyses of bond parameters. The geometry of a bond, and in particular its distance and angle influences its energetics and, consequently, the vibrational frequency. Previous research, has highlighted the variability of H-bond distances along helical peptides, discussing a trend of bond shortening at helix midpoints due to the cooperative effect of H-bonding along the length of the helix [34]. Similarly, Tan et al. have shown how the strength of an H-bond is dependent on the donor-acceptor angle, introducing the concept of an `antecedent angle’, which varies across protein secondary structures [35].

We calculated hydrogen bond distances and angles within our simulated systems. Specifically, we found the average O···N distance between Val80’s C=O and Val84’s N−H along the same helix to be $2.88\pm 0.06$ Å in the dimeric species of glycophorin A, compared to $2.97\pm 0.16$ Å in the monomeric variant. On the other hand, the C−O···N bond angle was $163.35\pm 2.7$ in the dimer and $160.81\pm 9.4$ in the monomer.

Hydrogen bonding interactions between donor and acceptor moieties can be both dipolar and electronic. An electronic interaction is possible when the donor N-H anti-bonding orbital overlaps the C=O oxygen nonbonding orbital. Electronic interactions are optimal when the C−O···N angle is near 120 [36], facilitating maximal overlap of these orbitals. But in protein systems, there are typically geometric constraints that restrict such bond angles to around 155 [35,37].

In the the dimeric model, with the additional inter-helical H-bond, a 3 wider C−O···N bond angle is observed relative to the monomeric model. This finding points to a reduced C=O to N−H canonical hydrogen bond strength. On the other hand, the canonical H-bond distance upon dimer formation was shortened by 0.09 Å , potentially increasing orbital overlap and thus increasing the canonical C=O···H−N interaction strength. Given the competing influences observed, MD simulations alone do not provide a definitive conclusion on the dominant effect on if the canonical hydrogen bonding strength should increase or decrease in the dimer system versus the monomer system. Such ambiguity is perhaps not unexpected due to the simplicity of the point charge model of the MD force field). Therefore, we proceeded to conduct quantum mechanical calculations for each system in which the C=O oscillator strength and other parameters can be estimated with greater accuracy, and the hydrogen bond strength directly inferred.

We calculated the frequencies of the amide I stretching modes (_C=0) of both the monomeric and dimeric systems. We chose several points of time from the MD simulation trajectory to make model systems for DFT calculations and expressed the results as a time average form. (see Experimental section) The average Val80 _C=0 value for the dimeric glycophorin A is $1656.1\pm 2.4$ cm₋1, whereas for the monomeric species (with the G79I mutation) the frequency is $1652.4\pm 2.3$ cm₋1. Notably, there is a 3.7 cm₋1 increase in the C=O stretching upon dimer formation which corresponds with the FTIR experimental _C=0 shift upon dimer formation. Hence, the surprising C=O bond strengthening is obtained both experimentally and computationally. The observed shift is independent of isotopic substitution. The ¹²C=¹⁶O variant of the construct produced a shift of 4.1 cm₋1 comparable to the ¹³C=¹⁸O analog that yielded an average shift of 4.2 cm₋1, arising from the stretching frequencies $1719.6\pm 2.3$ cm₋1 for dimer and $1715.4\pm 2.6$ cm₋1 for monomer.

In a typical H-bond, the acceptor (C=O in this instance) transfers electrons to the donor anti-bonding orbital, reducing the double bond character of the C=O, reflected as a decrease in _C=0. The H-bond of a C=O···H−D undergoes a hyperconjugative C+−O−H···D− type interaction (analogous to charge transfer) that in turn lowers the double bond character of C=O. [38]

To explore the electronic basis of the observed acceptor vibrational mode strengthening of the dimer system (upon over-coordination with the additional inter-helical hydrogen bond), we employed Natural Bonding Orbital (NBO) calculations to analyze the estimated electron distributions and atomic orbitals in a localized way. Local hybrid orbitals can be extracted using the NBO approach, which can correlate with classical qualitative bonding theories [39]. Alabugin et al. have shown how imperfect H-bonds can be analyzed by the NBO approach by showing the hyperconjugation effects [38].

According to valence bond theory, a C=O bond should be formed by the overlap of two sp²-orbitals from the C and O atoms to form a -bond, and two p-orbitals from the C and O atoms to form a -bond. Yet in reality, this picture changes depending on the electronegativity of each participating atom, as explained by Bent’s Rule [40]. However, for an ideal situation of a double bond, a -bond (like in the Val84 amide C=O) should project an sp² hybrid orbital to the other atom (O in this instance) present in the covalent bond, making the bond $33{\displaystyle \frac{1}{3}}$% s and $66{\displaystyle \frac{2}{3}}$% p. A deviation from these percentages will weaken the bond. In contrast, the -bond should be formed by two unhybridized 100% p-orbitals with no s-character.

Upon taking the hyperconjugative effect into account, the C=O bond should lose some of its double bond character due to the presence of a C⁺−O−H species and the localized orbitals will be far from an ideal sp² overlap that we would expect in the bond.

We ran NBO calculations on the energy-optimized DFT structures collected from MD-simulation trajectories and extracted the values of the s-character of the C and O atomic orbitals of the C=O bonding interactions. The results of the calculations, shown in Figure 3, indicate that both the monomeric and dimeric species experience a deviation in the s-character of the and -bonds from the ideal, $66{\displaystyle \frac{2}{3}}$% and 100% s and p-characters, respectively.

The dimeric form has a higher s-character in the -bond and a higher p-character in the -bond compared to the monomeric form. The average s-character of the atomic orbitals contributing to bond of C=O are shown in Figure 3, upper panels. It can be seen that the C and O hybrid orbitals posses a slight but definite higher s-character than the corresponding monomeric species. The monomer has a 29.5% and 37.3% s-character in its C and O hybrid orbitals, respectively. On the other hand, the dimeric species has a 30% and 38.8% of s-character in the respective orbitals. Therefore, upon dimerization and formation of the additional C−H···O bond, the C−O bond becomes stronger due to its higher s-character which can be a major contributing factor for the increased _C=0 stretching frequency. This kind of s-character enhancement has been previously encountered and has been accounted as the main factor for bond strengthening of the H-donor [38].

On the other hand, the -bonding C and O atomic orbitals in the dimer are 98.4% and 97.8% of p, respectively, whereas for the monomer they are 97.6% and 96.4% of p, respectively (Figure 3). So the C=O_dimer has more p-character in the sidewise overlapping orbitals than in the C=O_monomer, suggesting reduced hyperconjugation, which should theoretically enhance the oscillator strength of the C=O bond, supporting the observed frequency increase [38]. Finally, we note that indistinguishable results were obtained upon using a dispersion correction (D3_BJ) [41] method with a larger basis set (cc-pVTZ) [42].

These computational results suggest that the C=O bond strength increases upon dimerization, and the formation of the second H-bonding C−H group, due to the increased s-character of the bond and a simultaneous increase in the p-character of the bond which reduces hyperconjugation. Altogether, this leaves the C=O bond more double bonded in nature as observed in the later dispersion-corrected calculation.

We followed by analyzing how H-bonding impacts the charges of Val80’s carbonyl group. When the C=O···H−C interaction is absent, either in the monomeric system or in a separated dimeric assembly (Figure S6), a polar distribution is obtained: C^+0.45=O^-0.71 or C^+0.60=O^-0.65, respectively. In contrast, when the C=O···H−C is present (dimeric assembly) a markedly different charge distribution is obtained: C^-0.21=O^+0.09.

Lowering the C=O bond’s polarity should increase its covalent character. This factor, together with reduced hyperconjugation, potentially contributes to a stronger C=O bond.

It is worth mentioning that less hyperconjugation and a reduction in bond polarity do not imply that the C−H···O H-bond is energetically unfavorable. H-bonding involves both electronic and polar counterparts. Electronic interactions involve the transfer of electrons and the formation of charge-transfer species, whereas dipolar interactions are solely electrostatic in nature. An analysis of the energy terms associated with the canonical and noncanonical H-bonds in both monomer and dimer species reveals this complexity.

There is a small destabilization in the canonical H-bond of the dimer compared to the monomer, with the dimer having a canonical H bond energy value of $5.36\pm 0.12$ kcal/mol compared to that of the monomer at $5.71\pm 0.18$ kcal/mol. Hence, both bonds impact one another and are consequently not orthogonal in nature [43]. On the other hand, the noncanonical C−H···O H-bond in the dimer has a bond energy of $0.65\pm 0.10$ kcal/mol which indicates an overall stabilization of the dimer species due to the incorporation of two H-bonds. In other words, the additional H-bond from the glycine C−H strengthens the interaction between the two helices. Finally, the value of this additional noncanonical C−H···O H-bond is similar to previous measurements [23].

3. Conclusions

Taken together, we examined the C−H···O=C H-bond at the atomic level, using a synergistic combination of experiment and computation focusing on the carbonyl acceptor group. To our surprise, both lines of investigation show a blue shift of the acceptor stretching frequency indicative of a strengthening of the C=O bond, despite previous studies that have shown that proton acceptor groups often have reduced internal bonds upon increased hydrogen bonding. Molecular dynamic simulations followed by DFT calculations reveal that the C-H group cannot undergo effective hyperconjugation thereby increasing the C=O oscillator strength.

4. Materials and Methods

4.1. Sample Preparation

4.1.1. Isotopic Label Synthesis

Isotopic labeling was conducted as described previously [44]. Briefly, 4.52 mmol of 3,5-dimethylpyridine hydrobromide (Sigma-Aldrich; Rehovot, Israel) in 2 mL of anhydrous N,N-dimethylformamide (Sigma-Aldrich) was mixed with 2.24 mmol of N-(3-(dimethylamino)propyl)-N’-ethylcaerbodiimide hydrochloride (EDC·HCl) (Sigma-Aldrich) and with 11.3 mmol of H₂¹⁸O (Sigma-Aldrich) under N₂ atmosphere. 225 µmol of 1−¹³C−N−FMOC L-valine (Cambridge Isotope Laboratories, Inc; Andover, MA, USA), dissolved in 3 mL of anhydrous N,N-dimethylformamide, was added into the above-stated mixture. Subsequently, the reaction vessel was kept at room temperature and stirred overnight. After allowing the reaction to proceed overnight, another portion of 2.24 mmol EDC·HCl was added in the morning, followed by a third addition of 2.24 mmol EDC·HCl after approximately seven hours, and the reaction was allowed to continue overnight. Afterward, approximately 30 mL of ethyl acetate (Gadot-group; Netanya, Israel) was added. The mixture was transferred to a separatory funnel and washed three times with 0.1 M citric acid, and then once with brine. Ethyl acetate was then added to the combined citric acid and brine portions and separated. The combined portions of ethyl acetate, containing the labeled amino acid, were dried over anhydrous sodium sulfate (Da-sit Group; Milan, Italy) and filtered to remove any extra water. The ethyl acetate was removed by rotary evaporation, creating an azeotrope with methylene chloride (Gadot-group).

4.1.2. Peptide Synthesis and Purification

The labeled valine (see above), was incorporated into two different peptides corresponding to the transmembrane domain of glycophorin A [29]. Two peptides included the native sequence with Valines 80 labeled as well as an G79I [27] mutant with Valines 80 labeled in boldface (sequence starts at residue 70):

EPEITLIIFGVMAGVIGTILLISYGIRRL

EPEITLIIFIVMAGVIGTILLISYGIRRL

Both of the peptides were synthesized separately with N-(9-fluorenylmethoxycarbonyl) solid-phase chemistry. Each peptide sample was purified with high performance liquid chromatography on a 20 mL Jupiter 300 Å C4 5 µm high-performance liquid chromatography column (Phenomenex; Torrance, CA, USA). The column was pre-equilibrated with 80:8:12 (by volume) water:acetonitrile:isopropanol, where all solvents contained 0.1% trifluoroacetic acid (TFA) (Merck; Darmstadt, Germany). Approximately 2 mg of protein sample was dissolved in 2 mL of TFA and injected into the column. The solvent gradient was linearly altered with the VWR Hitachi Chromaster 5160 Pump (Tokyo, Japan) to remove all water composition while retaining the acetonitrile:isopropanol ratio at 40%:60% with 0.1% TFA. Peptide elution was monitored at 280 nm using the VWR Hitachi Chromaster 5410 UV detector.

4.1.3. Peptide Reconstitution

All experimental measurements were performed on peptides in lipid vesicles. We used organic solvent co-solubilization in order to reconstitute each peptide in a membrane bilayer. Approximately 1 mg of purified protein and 10 mg of 1,2-dimyristoyl-sn-glycero-3-phosphocholine (Avanti Polar Lipids; Alabaster, AL, USA) were dissolved in 1 mL of 1,1,1,3,3,3-hexafluoro-2-propanol (HFIP) (Merck, Darmstadt, Germany). The mixture was rotary evaporated at 37${ }^{\circ }$C until all the HFIP evaporated. One milliliter of water was added, and the mixture was rotated at 37${ }^{\circ }$C to spontaneously form vesicles. The sample was then sonicated to ensure uniformly sized vesicles and no aggregation.

4.2. FTIR Spectroscopy

For each of the two samples of peptides in a membrane vesicle, separate FTIR spectra were collected. First, 200 µL of sample was deposited on a germanium trapezoid ATR plate ($50\times 2\times 20$ mm) with a 45${ }^{\circ }$ face angle (Wilmad; Vineland, NJ, USA). Following removal of bulk solvent, the crystal was incorporated into a 25 reflection variable angle ATR unit (Specac; Orpington, UK), which reflects the incoming FTIR beam 25 times before its exit from the crystal. The ATR unit was incorporated within a Nicolet iS10 FTIR spectrometer, with a mercury cadmium telluride detector (Thermo Scientific; Waltham, MA, USA), cooled with liquid nitrogen. The FTIR spectrometer was purged with water- and CO₂-depleted air, and spectra were collected at room temperature.

For each sample, 1000 scans were sampled and averaged at a data spacing of 0.241 cm₋1 with two levels of zero filling, N-B strong apodization, and Mertz phase correction. For each of the two samples of peptides in a membrane vesicle, separate FTIR spectra were collected at room temperature. Prior to deposition of the sample on the germanium ATR plate, background spectra were collected with an empty germanium ATR plate and used to subtract background IR absorption.

4.3. Computational Details

4.3.1. Molecular Dynamic Simulations

PDB structure of the dimeric glycophorin A dimer (2KPF) was downloaded from RCSB PDB. The monomer was created by stripping out chain B from the dimeric structure and mutating G79 to I79 using Swiss-PDBeditor [45]. Both of the structures were then fitted into a pre-equilibrated 1,2-dimyristoyl-sn-glycero-3-phosphocholine bilayer obtained from supporting data of the study by Poger & Mark [46] to satisfy and mimic the experimental conditions. The fitting was done in such a way that the helical bundle(s) remain perpendicular to the bilayer plane. The bilayer contains 228 lipid molecules with 4166 water molecules. After initial structural alignment, lipid molecules that coalesced with protein coordinates were removed using InflateGro methodology developed by Tieleman and co-workers [47]. In this method an automatic running algorithm deletes molecules within 2 Å range. Position restraints (100,00 kcal mol^-1 Å^-2) on heavy atoms of protein molecules were imposed to make sure that the protein molecules do not change positions during the energy minimization process. Energy minimization was done with steepest descent minimization algorithm with a tolerance of 500 kj mol^-1 nm^-1.

Molecular Dynamic (MD) Simulations of the monomeric and dimeric protein were done for 200 ns using Gromacs version 2022.3 [48,49,50,51,52] using an extended version of GROMOS96 53A6 force field [53]. The length and angles of H-atoms were restrained with the LINCS algorithm allowing an integration time step of 2 fs [54]. Atomic coordinates were saved every 1000 ps. Reference temperatures were kept at 323 K and solvent, lipids, and proteins were separately coupled to a Nosé-Hoover temperature bath [54,55] with a coupling constant value = 0.5 ps. Pressure coupling was done with a Parrinello-Rahman barostat with = 2 ps [56,57]. A 1.2 nm distance was set as cut-off for van der Waals interactions. At every 10 fs the neighbor list was updated. 4^th order Particle Mesh Ewald (PME) long-range electrostatics was used to calculate electrostatic parameters [58].

The final simulation box had 225 lipid molecules. Hydration with FLEXSPC model [59] water molecules and addition of Na⁺ and Cl^- as counter ions led to a number of atoms of 25598 and 37365 for the monomer and dimer protein systems, respectively.

4.3.2. Quantum Mechanical Calculations

The protein H-bonding system was modeled with smaller molecular mimics. Our focus is the Val80 residue which is the proton acceptor in the C−H···O bond while simultaneously participating in a canonical H-bond with its H-bond- donating Val84 counterpart from the same helix. We represent this canonical H-bond between Val80 and Val84 with two N-methylacetamide molecules. The noncanonical H-bond-donating glycine moiety from the neighboring helix was modeled with an acetylglycinemethylamide. Schematic representations of the chemical formulae are shown in Figure S4.

Besides atoms in the capping CH₃ groups, all other atoms of the three molecules are derived from the original protein. Coordinates were chosen from the structure of glycophorin A and its monomer from different frames of MD simulations trajectory based on minimal RMSD values. The C atoms of the capping CH₃ groups were positioned where the C and C atoms were located in the PDB structures. Atom coordinates were selected from chains A and B of the structure, and the terminal C atoms of the N-methylacetamides and glycine were fixed during minimization. The α-carbon of the glycine mimicking molecule was similarly fixed.

While determining H-bonding energy, we required breaking the C−H···O H-bond. Breaking was accomplished by moving the glycine residue 30 Å apart from the original position, shown schematically in Figure S6. The ¹³C=¹⁸O system of Val80 was achieved by using ¹³C and ¹⁸O isotopes in the carbonyl group of the Val80-mimicking N-methylacetamide residue.

The computation was done using a Q-chem [60] set up at the B3LYP level of theory [61,62,63,64] using the “aug-cc-pvdz” basis set [65]. Such large basis sets reduce the chance of basis set superposition error. Moreover, the basis set superposition error was corrected using the DFT-C method available in the Q-chem package [66]. Furthermore, no difference was obtained upon repeating the computational analyses at the same level of theory (B3LYP), using a dispersion correction (D3_BJ) [41] method and with a larger basis set (cc-pVTZ) [42] to further minimize the effects of basis set superposition error.

The dielectric constant was fixed at 4 to mimic the membrane environment. Section 1 of supporting information shows the detailed methodology of the computational procedure. The structures of the model molecular assembly before and after energy minimization can be seen in Figures S4 and S6 in supporting information.

Q-chem package provides a Natural Bonding Orbital (NBO) calculation tool for measuring natural bonding orbital contributions involved in the chemical bondings. NBO outputs were extracted from the Q-chem output as percent orbital contributions. NBO output also involves Natural Population Analysis (NPA) for each construct which was also extracted.

5. Patents

This section is not mandatory, but may be added if there are patents resulting from the work reported in this manuscript.