Léa Rougier

Alain Milon

Valérie Réat

Franck Jolibois

Modelling the influence of hydrogen bonds network on chemical shielding tensors description. GIAO -DFT study of WALP23 transmembrane -helix as a test case

Density Functional Theory (B3LYP/6-31G(d,p)) calculations of 15 N amide and 13 C carbonyl NMR chemical shielding tensors have been performed on WALP23 trans-membrane -helix peptide and compared to solid state NMR experiment performed on [ 13 C1-Ala13, 15 N-Leu14] specifically labelled peptide powder sample. Using either theoretical result obtained on the whole peptide or experimental data as reference, several simplest chemical models have been explored in order to reduce the computational cost while maintaining good theoretical accuracy. From this study, it appears that the hydrogen bond (N-H…O=C) network that exists in the -helix has a major influence on the chemical shielding tensor and more specifically on the carbonyl 13 C 22 eigenvalue. We show that a small truncated WALP_7 model is not adequate for 13 C1 NMR description. The application of an external electric field in order to model the hydr ogen bond network allows calculating chemical shielding tensors with accurate eigenvalues while the associated eigenvectors are largely modified. Finally, a 23 residues polyglycine peptide that includes the Alanine and Leucine residues for which NMR parameters must be calculated is proposed as the chemical model. This model is sufficient to mostly reproduce the calculation performed on WALP23 with major gain in computational time. Moreover, the application of an external electric field allows reaching the experimental accuracy for the determination of the eigenvalues.

Introduction

Accurate determination of "static" chemical shielding anisotropic (CSA) tensors components remains a challenge in the field of solid state NMR (SSNMR) applied to biological systems whereas these quantities can be directly related to structure, dynamic and/or orientation [2][3][4] . Not only isotropic 5,6 but also anisotropic chemical shifts 7,8 are useful experimental data to probe secondary structure and short or long range electrostatic interactions. Since few years, 13 C and 15 N chemical shielding tensors of labelled biomolecules have received increasing interest for structure and dynamic analysis purpose using SSNMR.

If the backbone only is considered, carbonyl 13 C1 and amide 15 N CSA tensors are not only sensitive to secondary structure but can reveal insight into the nature of the hydrogen bonding pattern that can exits locally. More particularly, the 13 C1 σ22 tensor eigenvalue that is almost collinear to the C=O bond vector strongly depends on the distance of the hydrogen bonding partner 9,10 thus giving valuable information concerning the backbone geometry 11 .

The same trend is also observed for 15 N CSA tensor that are sensitive to local geometry but also to the distance to the C=O hydrogen bond partner [12][13][14][15] . 13 C and 15 N CSA tensors are not only helpful for structure analysis, but can also be employed for the study of the dynamic and/or the 50 orientation of peptides and proteins. Concerning 15 N CSA tensor, the main SSNMR experiment used for this purpose is the PISEMA technique which uses 15 N chemical shift and 1 H- 15 N dipolar coupling to define the orientation of a peptide structure relative to an oriented lipid bilayer 16 . We 55 have recently demonstrated that 13 C CSA tensors can also be used to determine orientation angle and order parameter of a biomolecule immerged in a lipid bilayer 4 . In this study, dynamically averaged 13 C CSA tensors were measured and combined to static 13 C CSA tensors obtained from quantum 60 chemistry calculations in order to define the orientation and dynamical parameters of ergosterol relative to a 1,2dimyristoyl-sn-glycero-3-phosphocholine bilayer. A similar approach has also been proposed combining 13 C1, 15 N chemical shift anisotropies and 2 H quadrupolar couplings 65 dynamically averaged by large amplitude movements for the study of the dynamic and orientation of a transmembrane peptide 17,18 . These studies have highlighted the need of accurate static CSA tensors in order to understand the dynamic of a membrane compound. Furthermore, 70 determination of the high resolution structure of transmembrane peptides or proteins is also constrained by the knowledge of such static CSA tensors 16,19 . While static CSA tensor eigenvalues can be determined experimentally on 13 C and/or 15 N specifically labelled small dry peptides, the associated eigenvectors i.e. the tensor orientation relative to the molecular frame may be difficult to obtain using experimental approaches 12 . However, such quantities can be "easily" accessible through the use of theoretical approaches such as quantum chemical methods.

Since several decades, quantum chemical approaches have been employed for the calculations of NMR parameters with real success in the fields of chemistry, physical chemistry and biochemistry. While the whole chemical shielding tensor is generally computed, the large 10 majority of NMR theoretical studies concerned isotropic chemical shift values. The chemical shielding tensor eigenvalues are less analyzed and this tendency is even more pronounced for 15 N [20][21][22][23][24][25][26] and for 13 C1 carbonyl 24,27,28 CSA tensors compared to aliphatic 13 C ones [29][30][31][32][33][34][35][36][37] . Moreover, the overall agreement for the 15 N and 13 C1 CSA tensor calculations with experiment is not as great as for aliphatic Cα tensors. If the targeted information is the chemical shielding tensor eigenvectors, the theoretical orientation of CSA tensors show large variations from molecule to 20 molecule [22][23][24][25] .

Consequently, one of the remaining issues in chemical shielding tensor calculation is not the development of new ab initio or DFT methods but mainly the correct modelling of molecular systems that are used in theoretical calculation 25 compared to experimental real molecule under investigation. More precisely, the local and non-local environmental effects must be correctly modelled in order to obtain accurate magnitudes as well as orientations of CSA tensors 38 . For instance, it is known that 15 N CSA tensors are not only sensible on primary and secondary structures but also on long range electrostatic interactions up to 10 Å 25 . In this paper, we present a theoretical DFT investigation of the amide 15 N and carbonyl 13 C1 CSA tensors of a synthetic α-helical transmembrane peptide, namely WALP23 (see Figure 1 for schematic representation), which has been introduced and extensively studied in the group of . In a recent work, we have employed this peptide in order to develop a new strategy for the analysis of the orientation and dynamic of transmembrane peptide 18 . This study was based on the analysis of the 13 C and 15 N dynamically averaged CSA tensors measured by solid state NMR experiments. It required static 13 C1 and 15 N CSA tensors obtained by a quantum chemical approach. In the present work, from an experimental point of view, tensor amplitudes (or eigenvalues) have been measured on two different samples 50 based on a 13 C1, 15 N specifically labelled peptide: dry WALP peptide and lyophilised WALP peptide inserted in detergent micelles. Using various chemical models, we have built an inexpensive and accurate strategy for the calculation of such NMR parameters. One of the most 55 important results is that, on one hand, local hydrogen bond is necessary and sufficient for the accurate 15 N chemical shielding tensor calculation 40 . On the other hand, the whole hydrogen bond network has an important influence on the magnitude of the σ22 component of the 13 C1 tensor and must 60 be explicitly taken into account during the calculation. The use of an external electric field is discussed and shown to improve the prediction of experimental data.

Experimental Section

Experimental procedures 65 13 C1, 15 N isotopically labelled WALP23 [ 13 C1-Ala13, 15 N-Leu14] peptide (acetyl-GW2(LA)8LW2A-NH2) was synthesized using Fmoc/tBu solid-phase synthesis 18 . The purity of the peptide was better than 95%. Two powder samples have been prepared for 13 C1 and 15 N static CSA eigenvalue evaluation. 70 6.5µmoles of WALP23 peptide were dissolved in a mixture of tert-butanol/water (1:1) and lyophilized. The resulting powder was then packed in 3.2mm MAS rotor. For the second dry sample, 0.8µmoles of WALP23 and 54.9 µmoles of dodecylphosphocholine (DPC) were mixed in 400µL of 2,2,2-75 trifluoroethanol (TFE). The organic solvent was evaporated under a nitrogen flow and further removed under vacuum overnight (ca. 1*10 -2 mbar). The sample was hydrated with milli-Q water and lyophilized to yield a fluffy powder. The powder was packed in 3.2mm MAS rotor. 80 Solid-state magic angle spinning (MAS) experiments were carried out on a Bruker Avance narrow bore spectrometer operating at 700.13 MHz for 1 H. A Bruker 3.2mm MAS tripletuned solenoid coil was used for 13 C1 and 15 N experiments. CP-MAS spectra were acquired using a 1 H excitation pulse length 85 of 2 µs and a CP spin-lock field-strength of 50 kHz. The CP contact time was 2 ms and 1.25 ms for 13 C and 15 N spectra, respectively. The 1 H radio-frequency field-strength for heteronuclear two-pulse phase-modulation (TPPM) decoupling was 125 kHz during acquisition. The repetition delay was 2s. 90 Spinning sideband 15 N and 13 C MAS spectra were recorded at 3 and 4 kHz, respectively. 13 C spectra were referenced to DSS and 15 N spectra to liquid ammonia, using the 1 H chemical shift of DSS in water as an external reference. The Static CSA eigenvalues 11, 22 and 33 (with |33-iso| ≥ |11-iso| ≥ |22-iso| 95 where iso, the isotropic chemical shift) [41][START_REF] Mehring | Principles of High Resolution NMR in solids 2nd[END_REF][START_REF] Spiess | [END_REF] were obtained by the analysis of spinning sidebands intensity. This analysis was performed using dmfit freeware (http://crmht-europe.cnrsorleans.fr/dmfit/). Eigenvalues determined for both samples are given in Table In order to evaluate the influence of the chemical environment on both CSA values, dry powder of "pure" peptide and dodecylphosphoglycol (DPC) lyophilized micelles that contains WALP peptides have been employed (see previous paragraphs). The second sample allowed us to remove possible 105 peptidepeptide interactions. DPC has been employed instead of standard DMPC phospholipids in order to avoid difficulties arising from the lipids carbonyl 13 C resonances.

Quantum calculation

All calculations (geometry optimization and NMR properties determination) have been performed using the Gaussian03 suite of programs [START_REF] Gaussian | [END_REF] . Because WALP23 experimental structure has not been determined by standard methods, perfect canonical helix (φ = -58° and ψ = -47°) has been used as starting geometry for the optimization procedure. The molecular system being large (367 atoms), the nature of the stationary point found by geometry optimization procedure has not been verified by Hessien matrix calculation as usually. Three type of molecular systems have been investigated, the entire WALP peptide (noted WALP23), a model peptide reduced to 7 amino acids (noted WALP_7) and a model peptide of 23 residues where all amino acids were replace by Glycine except the specifically labeled ones (noted polyGLY-AL). For the entire peptide (23 amino acids), geometry optimization have been performed using a 2 layers ONIOM approach 45, 46 . 20 The high and low layers were defined according to experimental 13 C and 15 N labelling of WALP23. Because we are interested in calculating the NMR parameters of Alanine 13 carbonyl atom (noted Ala13 13 C1) and of the Leucine 14 amide nitrogen atom (noted Leu14 15 N), the principal structural feature 25 that must be conserved is the intramolecular hydrogen bonds between one residue i and the residues i+3 and i-3 along the alpha helix. Consequently, the high layer of the ONIOM procedure was composed of residues 10 to 17, the cut being performed at the alpha carbon level. The low layer was constituted by the rest of the peptide. From a computational point of view, the high layer part was computed using the hybrid Density Functional Theory B3LYP method 47,48 associated with a Pople type double ξ basis set augmented by polarization functions on all atoms (namely 6-31G(d,p)) 49 while the lower part was treated at the semi empirical AM1 level [START_REF] Dewar | [END_REF][51][52] . During the optimization process, Cartesian coordinates associated to the nuclei of the lower part were frozen in order to prevent suspicious geometrical deformation. Therefore, only the geometry associated to the high layer part 40 (residues 10 to 17) was optimized. If the model peptide with 7 amino acids (namely WALP_7) is considered, geometry optimization (when performed) has been carried out at the B3LYP/6-31G(d,p) theoretical level, the Cartesian coordinates of the extremity residues being frozen for the same reason that previously exposed. For polyGLY-AL system, WALP23 optimized structure has been used in order to build this new molecular system. All removed residue side chains have been replaced by hydrogen atoms. The additional hydrogen atoms have been placed on the C positions and the C-H distances have been adjusted to 1.07 Å. No further geometry optimization has been performed for this peptide. NMR shielding tensors have been computed using the hybrid Density Functional Theory B3LYP method associated with a Pople type double ξ basis set augmented by polarization 55 functions on all atoms (namely 6-31G(d,p)). NMR chemical shielding tensors have been computed using the Gauge Including Atomic Orbital method (GIAO) for the numerous advantages it offers [53][54][55][56][57][58] . According to the Haeberlen, Mehring and Spiess [41][START_REF] Mehring | Principles of High Resolution NMR in solids 2nd[END_REF][START_REF] Spiess | [END_REF] convention, eigenvalues have been classified 60 according to the following inequalities: |σiso -σ33| ≥ |σiso -σ11| ≥ |σiso -σ22| where σiso is the isotropic chemical shielding.

Results and Discussion

The 15 N and 13 C CSA tensors have been determined experimentally, by standard spinning side band analysis 65 method. Two sample have been compared, pure peptide and peptide "diluted" in a DPC environment to assess possible intermolecular interaction effects. It appears that both samples give similar eigenvalues with a maximum difference of 1 and 0.5 ppm for 15 N and 13 C eigenvalues, 70 respectively (Table 1). Consequently, it can be argued that no specific external interaction significantly influences CSA tensors eigenvalues of 15 N amide and 13 C carbonyl atoms located on the peptide backbone and that comparison with theoretical values based on an isolated peptide is relevant.

75

The first part of this work was devoted to the determination of the Ala13 13 C1 and Leu14 15 N NMR parameters of WALP23. First calculations have been performed on the entire peptide after geometry optimization (see materials and methods for details). Results are presented in Table 1 80 and exhibit a relatively good agreement when compared to experimental values. Comparison between experimental and theoretical isotropic chemical shift is dependent of the choice and the accuracy of the chemical shift reference and do not necessarily reflect the quality of tensor calculations.

85

In order to avoid this problem, differences between chemical shielding eigenvalues and their associate isotropic chemical shielding have been calculated and will be the only data discussed thereafter. In the case of 13 C carbonyl CSA tensor, if absolute differences between theory and 90 experiment are considered, the errors for σii-σiso are relatively small (between 1 and 7 ppm). Moreover if the RMSD is considered, the overall agreement is rather good with a RMSD of 5.0 ppm for Ala13 13 C1. From an experimental point of view, -helix and -sheet secondary 95 structures exhibit significant differences if one consider 13 C carbonyl CSA tensors. In the case of poly-L-alanine peptide, a difference of 13 ppm has been experimentaly determined for σ22-σiso depending on the secondary structure (14.7 and 1.9 ppm for -helix and -sheet respectively) 59,60 . 100 This difference is less pronounced for the two other tensor components (5 and 8 ppm for σ11-σiso and σ33-σiso, respectively). Compare to these experiments, our calculations reproduce correctly the general trend of 13 C carbonyl CSA tensor of -helix secondary structure. In the 105 case of 15 N amide, CSA tensor the difference between secondary structures is less pronounced and no tendency can be clearly extracted from experimental 15 N amide CSA tensor components 61 . In our case, the errors for σii-σiso are the same than for 13 C atom (between 2 and 6 ppm) but with 110 a smaller RMSD (4.5 ppm for Leu14 15 N). It should be stated that our calculations were performed in vacuum and did not include any environment contribution such as others peptides or DPC molecules that may influence spectroscopic data. One way to include these contributions would have been to use a QM/MM approach with an extended MM part including the surrounding molecules. This kind of approach, commonly employed in the field of biology to calculate properties of one fragment of a protein or to analyse the reactivity at a catalytic site, needs the knowledge of a three dimensional structure extracted from crystallographic data or classical molecular dynamics simulations. However, because no structure already exists concerning our specific samples, the calculation starting from scratch would have 10 been too expensive for the present work. In conclusion, on can consider that theoretical eigenvalues are calculated with enough accuracy (4% and 3.5% errors relative to the 13 C and 15 N tensor experimental spans (tensor span = σ33-σ11)) and can then safely be used for the study of the dynamic and 15 orientation of a transmembrane peptides 17,18 .

Table 1 Experimental ssNMR and theoretical (see text for details) 13 C1 and 15 N chemical shielding tensors characteristic (in ppm) for WALP23 specifically labelled at position Ala 13 ( 13 C1) and Leu 14 ( 15 N). σiso is the isotropic chemical shielding. σii-σiso (i=1-3) have been calculated in order to avoid isotropic reference problem for theoretical and experimental values. |Max Error| (in ppm) is the maximum absolute value of error found for the eigenvalue indicated in the round brackets. RMSD (in ppm) is the root mean square difference calculated for each atom using the three σii-σiso values. 13 In an effort to reduce computational time and to find a simpler molecular model allowing the reproduction of all experimental NMR parameters of WALP23 we tried using a smaller peptide consisting of 7 amino acids (WALP_7, see materials and methods for details). NMR parameters calculations reveal a poorer agreement between theoretical and experimental values if 13 C is considered (See Table 1). More particularly, Ala13 13 C1 σ22-σiso and σ33-σiso values exhibit an error almost twice larger than the maximum error found for WALP23. The same tendency is also observed for the RMSD value (11.4 ppm) compared to the one calculated for WALP23 (5.0 ppm). On the other hand, the discrepancy between WALP_7 and WALP23 molecular model is considerably smaller in the case of 15 N atom with a RMSD value compared to experiment equal to 4.9 and 4.5 for WALP_7 and WALP23 respectively. It can be concluded that this simpler molecular model and more specifically the inclusion of the local hydrogen bonding pattern alone is not sufficient to reproduce the entire peptide and that some information has been lost by removing the 16 remaining amino acids.

In order to recover part or totality of this information, two theoretical physico-chemical perturbations have been applied on WALP_7. The first one has consisted in 45 simulating the rest of the peptide (16 amino acids) using a background charge distribution. For the second perturbation, a homogenous electric field has been applied along the helix axis.

Background charge perturbation. 50 One way to correct the WALP_7 molecular model was to simulate the part of the molecule that has been removed from the entire peptide by using a background charge distribution. In this context, atoms of this part of the peptide have been replaced by partial atomic charges located at the 55 nuclei Cartesian coordinates. Three type of atomic charge models have been employed: Mulliken type atomic charges, natural atomic charges extracted from Natural Population Analysis (noted NBO 62 ) and MK atomic charges which are charges fitted to reproduce the molecular electrostatic potential at a number of points selected around the molecule according to the Merz-Singh-Kollman scheme 63 . These three type of charges were extracted from a DFT B3LYP/6-31G(d,p) calculation performed on WALP23. It is well known that Mulliken type atomic charges are strongly dependent on the basis set and these values must be treated with care. However, the two other types of atomic charge are less sensitive to the choice of the basis set. Thus, results associated to these two type of atomic charges can be considered more reliable. 13 C and 15 N chemical shift tensors parameters evolution compared to experiment calculated for WALP_7, WALP23 and WALP_7 surrounded by charge distribution are presented in Figure 2. Concerning 15 N nucleus, charge distribution has little effect on discussed NMR parameters, leaving these values closed to WALP23 as already obtained for WALP_7 alone. When Ala13 13 C1 is considered, charge distribution displays a more significant influence on NMR parameters. While σ11-σiso and σ33-σiso are less perturbed by the presence of charge distribution, σ22-σiso is clearly but not sufficiently improved, the difference between theory and experiment ranging from 15.8 ppm (WALP_7 alone) to 10.4 ppm (WALP_7 with Mulliken charges) compare to 5.5 ppm for WALP23. However, if partial charges are calculated for atoms of Ala13 13 C1 and Leu14 15 N using Mulliken approximation (Table 2), sum of atomic partial charges in N-H bond (that reflects the charge repartition inside the covalent bond) does not vary when applying charge distribution (difference compare with WALP23 remains equal to 10%). On the other hand, while for WALP_7 this sum is 50% larger than for WALP23 in the case of C=O bond, the application of a charge distribution reduces this discrepancy to about 25%. Thus, it appears that the background charge distribution has a real influence on the polarization of the hydrogen bond and more specifically on its carbonyl part but is not sufficiently accurate to represent the extremities of the peptide in order to reproduce NMR parameters of the entire peptide. These observations concerning H-bond polarization associated to the analysis of the behaviour of chemical shielding tensor discussed above highlight the fact that good modelling of hydrogen bond type interactions and more specifically the consideration of the whole alpha helix hydrogen bond network has to be taken into account in order to determine accurate NMR chemical shielding tensors. For this purpose, we propose to apply a uniform electric field in order to reproduce the influence of the hydrogen bond network. From a physical point of view, this idea is similar to the one proposed by Torri H theoretical work on retinal molecule 64 . In this study, an electric field is applied along the principal axis of the molecule in order to 55 reproduce the effects on its electronic absorption maximum of the dipolar residues that surround the retinal inside the opsin protein. In should be stressed that in our approach the role of the external electric field is not to reproduce the external environment but rather to model the effect on 60 bonds polarization of the whole hydrogen bonds network that has been lost when shortening WALP23 into WALP_7. This approach is different from others studies that employ atomic charges or charge dipoles in order to reproduce specific effect arising from charges or polar groups in 65 proteins or to generate dipolar perturbations [65][66][67][68][69][70][71] .

Uniform electric field perturbation

In order to essentially try to reproduce the influence of the 75 alpha helix hydrogen bonds network, the electric field has been applied to WALP_7 along the inertia axis of the helix fragment, from N terminal to C terminal atoms. Electric field intensities have been set from 0 to 300 10 -4 atomic unit (noted a.u. with 1 a.u. = 514.10 9 V.m -1 ). To estimate the 80 influence of the electric field on the polarisation of the reduced hydrogen bonds network of WALP_7, the average value per amino acids of the dipole moment has been calculated and is given as a function of the electric field intensity in Figure 3. While the dipole moment value is 85 definitively too low when WALP_7 is considered (3.7 Debye / amino acids), the electric field significantly induces a regular increase of this value that can reach and goes beyond the one obtained for WALP23 (5.2 Debye / amino acids) for electric field intensities larger than 100 10 -4 a.u..

90

This first observation clearly demonstrates the effect of the electric field on the polarisation of the functionnal groups involved in the hydrogen bonds and thus as a possible way to model the hydrogen bond network existing in the entire peptide. However, if the electric field intensity become too 95 large (superior to 200 10 -4 for WALP_7), the dipole moment evolution as a function of electric field intensity presents a discontinuity. Before this electric field value, the dipole moment evolves linearly with a slope equals to 155.8

Debye/(electric field a.u.) for WALP_7. After this electric field intensity value, the slope becomes larger (282.1 Debye/(electric field a.u.)). This reflects the fact that with such intensity, charge separations become too large reflecting the unphysical nature of this calculation for such "high" intensity electric field. The evolutions of 15 N and 13 C1 chemical shielding tensor parameters (σii-σiso) compared to experiment are given in Figure 4 as a function of electric field intensity. If 13 C1 nucleus is considered, while, as previously stated, the main problem arises from σ22-σiso and σ33-σiso values for WALP_7, the application of the electric field clearly improves these values compared to experiment when electric field intensity is set between 120 10 -4 and 180 10 -4 20 a.u.. It must be noted that calculation using an electric field intensity set to 100 10 -4 a.u. gives three eigenvalues almost identical to the ones obtained for WALP23. In the case of 15 N nucleus, the effect of the electric field on chemical shielding tensor parameters seems less sensitive, NMR 25 calculations on WALP_7 being already in relatively good agreement with experiment. While the absolute values of σii-σiso are modified, their differences with experiment remain in an acceptable range of ±6 ppm. The best agreement over the 6 σii-σiso values compared to experiment 30 data has been found for an electric field of 140 10 -4 a.u. (see values in Table 1) and is associated to the smallest total RMSD calculated using both atoms of 3.8 ppm. However, when compared to WALP23, a large deviation of all the 15 N σii-σiso is observed when an electric field is applied.

35 Consequently, while it is possible to reproduce WALP23 13 C CSA tensor with an appropriate electric field applied on WALP_7, such perturbation has an inverse results on 15 N CSA censor. Fig. 4 Ala13 13 C1 (a) and Leu14 15 N (b) σii-σiso (i=1-3) difference between theory and experiment as a function of electric field applied to WALP_7 (D[σiiσiso] = (σii-σiso)theo -(σii-σiso)exp). Circle i=1, square, i=2, star i=3. In solid line is represented the total RMSD between theory and experience calculated using the 6 σii-σiso values as a function of electric field.

In order to highlight the effect of the electric field on the 13 C1 and 15 N chemical shielding tensors eigenvalues, we have analyzed the eigenvalues variation as a function of electric field intensity taking WALP_7 eigenvalues without electric field as references. It can be stressed that the most influenced values are σ22 and σ33 for 13 C1 and 15 N, respectively. This can be explained by the fact that the tensor eigenvectors associated to these eigenvalues are more or less collinear to the direction of the applied electric field. In the case of 15 N, this eigenvector is about 20° deviated from the N-H bond direction while for 13 C1, the eigenvector direction is almost collinear to C=O bond (Figure 5). As a consequence, the applied electric field directly polarizes the electronic density along C=O and N-H bonds, the effect of the entire H-bond network being partly recreated. Thus, the electric field effect on 13 C1 σ22 eigenvalue is more important than on 15 N σ33.

Considering that DFT methods and/or the atomic orbital basis set employed may have difficulties in correctly representing weak electrostatic interactions such as hydrogen bond and/or that geometry optimization has not been performed on the whole peptide, it can be argued that 65 the application of an electric field on WALP_7 has two consequences: modelling the hydrogen bond network of the entire peptide and correcting the theoretical method (method, basis set, optimization) for weak interactions. Furthermore, if the type of sample used for NMR 70 experimental measurement is considered (powder sample or lyophilized DPC micelles), dipoledipole intermolecular interactions have a small influence on the polarization of the hydrogen bond network and thus slightly modify NMR parameters compared to an isolated peptide in gas phase 75 which is the computational condition (See experimental section). In order to verify these assumptions, NMR calculations have been performed on WALP23 in the presence of an applied electric field. The same tendency, improvement of NMR parameters with the application of electric field, is observed for WALP23, thus confirming the correction of our theoretical model by the use of an electric field (See Supp. Info). In this case, because the hydrogen bond network is explicitly present in the molecular structure, the intensity necessary for a good description of chemical shielding eigenvalues compared to experiment is less important, the best value being 20 10 -4 a.u. compared to 140 10 -4 a.u. for WALP_7. With these electric field intensities, differences of σii-σiso compared to experiment 10 are less than 6 ppm, the RMSD for the 6 eigenvalues being equal to 3.1 ppm and 3.8 ppm for WALP23 and WALP_7, respectively. Dipolar moment evolution as a function of electric field has also been analysed in the case of WALP23 (see Figure 3). As previously observed with WALP_7, an increase of the dipolar moment per amino acid is observed with a slope higher than for WALP_7 (222.7 Debye/(electric.field a.u.)). Similarly, for a high intensity electric field (superior to 60.10 -4 a.u.), one observes a discontinuity that indicates an unrealistic polarisation of all 20 the bonds. With intensity set to 20.10 -4 a.u. (the one for which best NMR values are obtained compared to experiment), the dipole moment per amino acid is close to the one determined with an intensity set to 140.10 -4 a.u. for WALP_7 that gives the best NMR values compared to Influence of peptide planes number.

In order to obtain further insight on the influence of the hydrogen bond network extent, we have performed NMR calculations on WALP peptides of increasing size. Our chemical models were based on Ala13-Leu14 fragment (Me-CO-Ala13-Leu14-NH-Me). One peptide plane has been added at each side of this initial peptide, increasing its size by 2 peptide planes each time until we reach the whole WALP23. The WALP23 structure has been used and no 40 geometry optimization has been performed in order to evaluate only the hydrogen bond network effect on NMR parameters. If 15 N σii-σiso NMR parameters are considered (Figure 6a), it appears that a minimum of 7 peptide planes is necessary for accurate results compared to WALP23 data.

45

This is consistent with our previous results obtained for WALP_7. Thus, in the case of 15 N nucleus, the local hydrogen bond is necessary and sufficient to reproduce the spectroscopic observables and no long range interaction affects the calculation of these parameters 40 . In the case of 50 13 C1 σii-σiso NMR parameters, the situation is more complicated (6b). According to our ordering of chemical shielding eigenvalues relative to the isotropic value, one can see that the first hydrogen bond involving directly the carbonyl group is absolutely necessary for a correct 55 description of the tensor. From one to five peptide planes, σ11 and σ33 are inverted thus leading to a chemical shielding anisotropy (σ33-σiso) with an inverse sign (positive instead of negative). Moreover, while a number of peptide planes higher than seven, does not significantly improve σ11 and 60 σ33 eigenvalues, σ22 eigenvalue needs a more extended hydrogen bond network (at least two helical turns) in order to achieve a good computational accuracy. Thus, the cooperativity, that proceeds from interactions among residues and that has been highlighted previously 72 , has its 65 most important effect on 13 C1 σ22 eigenvalue that is associated to e22 eigenvectors collinear to the hydrogen bond network.

Simplified WALP model: Polyglycine. 70 Previous calculations have demonstrated that one of the most important structural feature that must be conserved in the chemical model is the hydrogen bonds network. Indeed, while the alpha helix structure is conserved in the WALP_7 model, chemical shift tensor are not well reproduced. The 75 hydrogen bond network not only influences the entire structure in a specific secondary structure but also provides bonds polarisation that affects chemical shift tensors principal components. Consequently, a simpler model of WALP23 has been proposed by removing all the amino 80 acids side chains in a first approximation. As expected, chemical shielding eigenvalues are correctly calculated for 15 N compared to experiment (see Table 1) while a poorer agreement is obtained for 13 C1 eigenvalues and more specifically for the σ22-σiso value. In this case, an error of 85 14.1 ppm relative to experiment is obtained, compared to 15.8 ppm for WALP_7.

In order to correct this chemical model, Alanine 13 and Leucine 14 side chains were kept, the rest of the polypeptide amino acids being replaced by glycine residues 90 (noted polyGLY-AL). In this case, chemical shielding eigenvalues are very close to the one obtained for in vacuo WALP23 (max error equal to 7.5 ppm compared to 6.8 ppm, respectively). The Ala13 13 C1 σ22-σiso value has been considerably improved compared to calculation on polyglycine (difference compared to experiment of 7.5 ppm instead of 13.4 ppm). Finally, on this last chemical model (polyglycine with Ala13 and Leu14), an external electric field has been applied using an intensity of 20 10 -4 a.u., this intensity being the best one obtained for WALP23 (see previous paragraph). With this external perturbation, theoretical calculations are also improved compared to in vacuo calculations. 15 N σii-σiso NMR parameters are in excellent agreement with experiment (RMSD ( 15 N) = 2.6 ppm). Paradoxically, it appears that compared to experiment, the simple polyglycine model gives better results for 15 N nuclei than WALP23. This highlights the fact that side chains may play an significant role for chemical shift tensor properties and that the theoretical WALP23 model may not correctly reproduces these side chains contributions. It is suspected 20 that side chain dynamics, i.e. rotation about C-C bonds, may play an role and influence the 15 N chemical shift tensor. This dynamic aspect is under investigation.

Tensor orientation analysis.

As stated in the introduction, quantum chemical calculations are essential for the determination of the chemical shielding orientation with respect to molecular frame. Thus, in order to check our methodological approach, the description of chemical shift tensor eigenvectors in the molecular frame was determined. For practical reasons, the molecular frame has been chosen relative to the Ala13-Leu14 peptide plane for all chemical models in use: x-axis is collinear with N-H, zaxis is normal to the peptide plane and y-axis is almost collinear to the N-C(i-1) vector. For 15 N, N is the angle between e11 and the N-C peptide bond, βN the angle between e22 and the z-axis and γN, the angle between e33 and the N-H bond where eii is the eigenvectors associated to the σii eigenvalue (see Figure 5 for eii representation). For 13 C1, C is the angle between e11 and the N-C peptide bond, βC the angle between e22 and the C1-O bond and γC, the angle between e33 and the z-axis. Because no experimental values for eigenvectors orientation exist for WALP23 polypeptide, the in vacuo calculations has been taken as references for each molecular system. The angles obtained for in vacuo are in agreement with the standard values used to describe 15 N and 13 C1 chemical shift tensors eigenvectors However, when external electric field is applied, several angle values are significantly modified as a function of the field intensity. For clarity reason, only in vacuo results and calculations performed using the optimal electric field intensity (140 10 -4 for WALP_7 and 20 10 -4 for WALP23 and polyGLY-AL) are presented in Table 3. Others results with several electric field intensities are given in supplementary materials, (Figure S1). It should be stressed that for the three molecular models, the 15 N e33 eigenvector orientation is almost not modified by the application of the electric field. As stated in previous paragraphs, an electric field with an intensity of 140 10 -4 a.u. is necessary to obtain good eigenvalues for WALP_7. Nevertheless, using such intensity for the applied electric field modifies the eigenvectors orientation by up to 16° for 13 C1 and 35° for 15 N. Thus, while this electric field is adequate for correcting chemical shielding eigenvalues, it is not suitable for 10 correcting chemical shielding eigenvectors. In WALP23 or polyGLY-AL model, the intensity of the electric field (20 10 -4 a.u.) applied for chemical shift tensor correction is less important. Consequently, the deviation of the eigenvector orientation from unperturbed chemical shielding tensors is rather small (max difference = 4°) and one can safely use such computational strategy in order to obtain accurate chemical shielding eigenvalues and eigenvectors with minimal computational effort.

Conclusions 20

In this paper, an extensive computational study of the NMR properties at the GIAO -DFT level has been proposed on WALP23 -helix transmembrane peptide. For this purpose, solid state NMR experimental data of a dry specifically labelled peptide ( 13 C1 Alanine 13 and 15 N Leucine 14) have served as references in order to calibrate an accurate and "inexpensive" computational strategy. While calculation on the whole peptide reveals good agreement compared to experiment, this molecular model is too large to envisage performing calculations on this system routinely.

In order to reduce computational time, a smaller peptide WALP_7 has been proposed that includes residues for which NMR parameters must be calculated and on helical step in order to include local hydrogen bond associated to the peptide bond chemical groups (amide and carbonyl groups). Our calculations clearly show that this molecular model is sufficient to reproduce 15 N amide chemical shift tensor while the 13 C1 carbonyl chemical shielding tensor is calculated with less accuracy. It is shown that this discrepancy is due to the bad reproduction of the hydrogen 40 bonds network created by the whole peptide backbone that adopts an -helix type secondary structure. Application of an external electric field is proposed in order to model the hydrogen bond network effects on local hydrogen bond. While this approach allows calculating chemical shielding tensors with accurate eigenvalues, the associated eigenvectors are largely modified. Moreover, at this stage of the study, no control on the choice of the electric field intensity is possible. It should be necessary to perform similar studies on peptides of various size in order to highlight possible correlation between the applied electric field intensity and the peptide size. Consequently, WALP_7 peptide is not a good alternative for reproducing NMR characteristics of WALP23. In order to explicitly include the hydrogen bond network, a 23 residues polyglycine that 55 includes the Alanine and Leucine residues for which NMR parameters must be calculated is proposed as the chemical model. It is shown that this molecular model associated (or not) to a low intensity external electric field allows calculating 15 N and 13 C1 CSA tensors (eigenvalues and 60 eigenvectors) with a quantitative accuracy compared to experiment. These calculations are carried out with a computational cost divided by almost one order of magnitude compared to the calculations on the whole peptide.

65

In the future, the computational strategy presented in this paper will allow performing NMR calculations on larger and more complicated peptide in terms of sequence. However, further studies concerning the influence of amino acid side chains (more particularly, polar and charged 70 residues) remain necessary in order to estimate to what extent these residues perturb chemical shielding tensors of 15 N amide and 13 C carbonyl atoms that are located on the peptide backbone. This strategy will also permit to carry out, within a reasonable time, calculations on a large 75 number of molecular structures that can be extracted from molecular dynamic simulations, for example, in order to introduce fast nuclear motions in the theoretical description of chemical shielding tensors.

Fig. 1

 1 Fig. 1 Schematic representation of WALP23 transmembrane peptide.

40 A

 40 . Killian (see review39

1

 1

  .

  100

Fig. 2

 2 Fig. 2 Ala13 13 C1 (left) and Leu14 15 N (right) σii-σiso (i=1-3) difference between theory and experiment as a function of punctual charge model used for charge distribution perturbation of WALP_7 (D[σii-σiso] = (σii-70

Fig. 3

 3 Fig.3 Dipolar moment (Debye) of WALP_7 (circles) and WALP23 (squares) calculated per amino acid as a function of the applied electric field intensity. The dashed line represents the usual theoretical dipole moment associated to the peptide group (3.45 D) 1

Fig. 5

 5 Fig.5 Orientation of eigenvectors of the 2 CSA tensors (Ala13 13 C1 and Leu14 15 N) and direction of the applied electric field.

Fig. 6

 6 Fig. 6 Leu14 15 N (a) and Ala13 13 C1 (b) evolution of σii-σiso (i=1-3) as a function of number of peptide planes (helix size). Circle i=1, square, i=2, star i=3. Horizontal lines are theoretical WALP23 references.

Table 2

 2 Mulliken partial charges of the atoms of Ala13 13 C1 Leu14 15 N as a function of chemical models. WALP_7 means in vacuo WALP_7. NBO, MK, Mulliken mean WALP_7 with background charge distribution. WALP23 means in vacuo whole peptide.

	Alanine 13 WALP_7	NBO	MK	Mulliken WALP23
	C1	0.617	0.617	0.621	0.617	0.61
	O	-0.532	-0.536	-0.545	-0.547	-0.55
	C1=O	0.085	0.081	0.076	0.070	0.057
	Leucine 14 WALP_7	NBO	MK	Mulliken WALP23
	N	-0.556	-0.556	-0.554	-0.556	-0.55
	H	0.307	0.308	0.305	0.312	0.33
	N-H	-0.249	-0.248	-0.249	-0.244	-0.225

Table 3

 3 Ala13 13 C1 and Leu1415 N eigenvector orientation relative to the molecular frame determined without and with the optimal applied electric field (see text for details) for WALP_7, WALP23 and polyGLY-AL.

		( 13 C1)	( 13 C1)	(1 13 C1)	( 15 N)	( 15 N)	( 15 N)
	WALP_7
	In vacuo	36.4	4.2	2.8	20.3	16.6	18.7
	With field (140 10 -4 a.u.)	47.7	18.3	11.9	53.0	52.5	18.1
	WALP23
	In vacuo	31.8	2.4	2.0	20.7	17.7	19.5
	With field (20 10 -4 a.u.)	36.2	5.3	4.5	21.9	19.0	19.6
	polyGLY-AL
	In vacuo	32.5	2.4	2.4	19.2	15.7	19.3
	With field (20 10 -4 a.u.)	35.3	4.5	4.1	17.5	13.3	19.2

  14,73,74 .

	50
		Particularly, for 15 N, the angle between e33 and the N-H
		bond is equal to 19° and the e11 is tilted by about 19-20°
		from the peptide plane. For 13 C1, C angle ranges between
		32° and 36°, e22 is almost collinear to the carbonyl C1-O
	55	bond and e33 is perpendicular to the peptide plane. If in
		vacuo non perturbed chemical models are considered,
		chemical shift tensors orientations in the molecular frame
		are almost identical for WALP23, WALP_7 and polyGLY-
		AL the maximum difference being equal to 5° for C
		compared to WALP23 value (see Table 3).

60

Journal Name, [year], [vol], 00-00 | 3

Journal Name, [year], [vol], 00-00 | 7

| Journal Name, [year], [vol], 00-00This journal is © The Royal Society of Chemistry[year]

Journal Name, [year], [vol], 00-00 | 9

Acknowledgements 80

Authors wish to acknowledge financial supports from "Conseil Régional Midi-Pyrénées" and PRES "Université de Toulouse" for Léa Rougier Ph.D. funding and the CALcul en MIdi-Pyrénées (CALMIP) computing center for generous allocations of computer time.

85

Notes and references