Hydrogen Bond Dynamics of Histamine Monocation in Aqueous Solution: How Geometric Parameters Influence the Hydrogen Bond Strength[dagger]

Abstract

Chemometric statistical approaches involving multiple linear regression (MLR) and principal component analysis (PCA) were employed on a set of 42 distinct snapshot structures of the physiological histamine monocation in aqueous solution along the Car-Parrinello molecular dynamics trajectory, in order to obtain a better insight into the relationship between the geometry parameters of the system and the resulting ν^sub NH^ stretching frequencies. A simple 2D linear regression of ν^sub NH^ with N^sub amino^...O^sub water^ distances gave a very poor correlation (R^sup 2^ = 0.42), but both MLR and PCA with the inclusion of four directly bonded water molecules offered a notably predictive model that is even able to distinguish two classes of structures based on the Cl^sup -^ counterion position. Taking into account waters from the first, second and third solvation shells, sequentially diminished the overall predictive ability of the model, yet increased the number of useful predictors that, in the largest model with 51 solvent molecules, all correspond to bulk water, implying that both chemometric methods are consistent in suggesting that fundamental histamine N-H stretching vibrations are very complex in nature and strongly coupled to the fluctuating environment.

Headnote

Abstract. Chemometric statistical approaches involving multiple linear regression (MLR) and principal component analysis (PCA) were employed on a set of 42 distinct snapshot structures of the physiological histamine monocation in aqueous solution along the Car-Parrinello molecular dynamics trajectory, in order to obtain a better insight into the relationship between the geometry parameters of the system and the resulting ν^sub NH^ stretching frequencies. A simple 2D linear regression of ν^sub NH^ with N^sub amino^...O^sub water^ distances gave a very poor correlation (R^sup 2^ = 0.42), but both MLR and PCA with the inclusion of four directly bonded water molecules offered a notably predictive model that is even able to distinguish two classes of structures based on the Cl^sup -^ counterion position. Taking into account waters from the first, second and third solvation shells, sequentially diminished the overall predictive ability of the model, yet increased the number of useful predictors that, in the largest model with 51 solvent molecules, all correspond to bulk water, implying that both chemometric methods are consistent in suggesting that fundamental histamine N-H stretching vibrations are very complex in nature and strongly coupled to the fluctuating environment.

Keywords: hydrogen bond, hystamine monocation, multiple linear regresion, principal component analysis

INTRODUCTION

Histamine is a biogenic amine and is mostly known for its role as an inflammation mediator. Its functions as a neuromodulator and neurotransmitter in the brain are less well established, but its known scope of actions is expanding. Recently it has gained importance as a signaling molecule in the processes of sleep-wake cycle, appetite control, learning, memory and emotion. Moreover, its signaling paths seem to be involved in conditions such as depression, schizophrenia, Alzheimer's disease and epilepsy.1 Consequently, processes involving histamine synthesis, transport, metabolism and binding to macromolecules have a significant physiological relevance. In the brain, histamine is synthesized from histidine by a specific enzyme L-histidine decarboxylase, and is taken up into synaptic vesicles by the vesicular monoamine transporter-2.2 It is released into synaptic cleft upon depolarisation stimuli and like every other neurotransmitter it has to be removed from synaptic cleft in a relatively short time. The primary histamine metabolising enzyme in the brain seems to be an intracellular histamine-N-methyl-transferase, which produces N-methyl-histamine that is later converted to the corresponding aldehyde by monoamine oxidase (MAO). In other pathways, histamine is degraded by diamine oxidase (DAO) and to some extent directly by MAO.

Concerning the nature of molecular interactions between histamine and its larger biological systems, these are dominated by the histamine hydrogen bonding ability,3 just as it is the case when it is surrounded with solvent water molecules under physiological conditions. Therefore, the hydration of histamine is of considerable interest, particularly with regard to structures and vibrational spectra. Histamine is composed of an imidazole ring and an aminoethyl side-chain, and according to its p/Y, values (5.8 for the imidazole ring nitrogen atom and 9.4 for the aliphatic amino group),4 histamine is predominantly a monocation (96 %) protonated at the free-amino position.5 In cells, histamine recognition induces (at least) partial desolvation, implying that a portion of water molecules from hydration shells must be removed to make a way for the interact tion with proteins and enzymes. A very illustrative example is provided by the histamine receptors H1-H4, which are all activated through the hydrogen bonding to histamine,6 suggesting that it is highly desirable to properly understand its interactions with physiological water and biological macromolecules. The latter is in line with recent molecular dynamics simulations by Keseru and co-workers7 that emphasized the role of water-histamine interactions at the binding site of the human histamine H4 receptor, being in agreement with several previous reports.8

Recently, we modeled the nature of the hydration of physiological histamine monocation by using CarParrinello molecular dynamics simulations (CPMD)9 and a static cluster model.10 In the former study, we simulated the N-H stretching band envelope, which are degrees of freedom that predominantly influence the strength of the hydrogen bonding as probed by the vibrational spectroscopy. The approach involved the construction of one dimensional proton potentials on snapshot structures extracted from a CPMD trajectory, followed by solving the one-dimensional vibrational Schroedinger equation. By separating the amino N-H vibrations of the protonated ethylamino chain from those within the imidazole ring in the experimental spectra, we clearly demonstrated that the former moiety forms stronger hydrogen bonds with the surrounding water molecules,9,10 which should be significant for the biological function of histamine. It remains a challenge to identify how the fluctuating intra- and intermolecular geometric parameters influence the N-H stretching frequencies, which is the aim of the present work. The most important factor seems to be the proton donor-proton acceptor distance among histamine amino nitrogens and water oxygen atoms. Careful comparisons between the experimental donor-acceptor distances (taken from the frozen crystal structures) and the corresponding stretching frequencies for a series of related systems are known as Novak11 and Mikenda12 relations. Analysis of the CPMD trajectory offers an advantage that it is performed for several configurations of the same system and that all other geometric features could be considered. This approach was previously used in the case of a very strong hydrogen bonding in the crystals of sodium hydrogen bis(sulfate),13,14 where, remarkably, no appreciable correlation between the instantaneous geometric parameters and the pertinent OH stretching frequency could be found, demonstrating the extreme complexity of the short and strong hydrogen bond. In contrast, aqueous histamine forms relatively long and weak hydrogen bonds with water;910 hence part of the motivation for the present study is a comparison of the geometry-frequency correlation patterns between strong and weak hydrogen bonds.

COMPUTATIONAL DETAILS

As part of our interest in examining the hydrogen bond dynamics of the histamine monocation in aqueous solution,910 the present work analyzes the snapshot structures along the CPMD trajectory with chemometric methods - a combination of the mathematical and statistical approaches designed to provide maximum chemical information of the chemical system by analyzing chemical data. For the CPMD simulation details, extraction of snapshot geometries, and the calculation of instantaneous proton potentials and thus derived anharmonic frequencies, the reader is referred to reference 9. Briefly, we collected 42 distinct geometries during 10 ps simulation, involving the histamine monocation, 93 water molecules and CP counterion, and since each snapshot structure consists of four N(immo' * *Ov,(lk.| moieties, we tailed a set of 168 anharmonic frequencies.

In order to study the relationship between geometric parameters and the corresponding νNH frequencies, we started the analysis by considering all internal coordinates of the studied system (termed predictors throughout the text), based on the fact that the main type of interactions in the simulation box is electrostatic in nature, which is essentially long-range in character, therefore the proton dynamics is controlled by surroundings that extend beyond the nearest atoms. However, this made a visual aspect of such multi-dimensional correlation a very difficult and challenging task. To overcome this problem, we reduced the data dimension by employing the Principal Component Analysis (PCA), being recommended as the most appropriate method for preprocessing multi-dimensional data into two dimensions.15 We were also interested in finding a relationship between a single dependent variable (i\n) and several independent variables (geometry parameters), for which we utilized the Multiple Linear Regression (MLR) as a calibration technique.15,16

As first, MLR and PCA analysis have been applied on the matrix of 42 x n dimensions, corresponding to 42 snapshot structures with n internal coordinates in columns (see Table 1). The elements of particular columns are given as a quotient of the structural parameters (bond distances, valence angles, dihedral angles). The total number of variables is obtained in the following way. For the system with N atoms, the number of internal degrees of freedom is 37V-6, yielding 888 unique internal coordinates that describe each snapshot structure with 298 atoms.9 Since one internal coordinate of each H-bonded moiety served as a "probe" coordinate along which the proton was displaced in a stepwise manner when calculating the corresponding proton potentials, we omitted it from our analysis. For this purpose, the periodic boundary conditions and replacement of water molecules were also considered. Our initial inspection suggested poor correlation between all geometry parameters and vnh at best, giving evidence that the investigated HB interactions are extremely complex as was the case with the strong hydrogen bond in sodium hydrogen bis(sulfate).14 However, in contrast to the virtually nonexistent correlation reported in the latter study, here we analyzed each of the individual subsets of geometric parameters separately (bond distances, bond angles, dihedral angles), and found that when examining only bond distances there appears to be a certain correlation with vnh values. Therefore, MLR and PCA have been applied on the reduced matrices of 42 x m dimensions, where m denotes the number of bond distances, yielding matrices M1-M4 with m = 30, 60, 129 and 171, respectively (Table 1). We further refined the process by sequentially increasing the number of considered water molecules, starting with four directly bonded molecules, and proceeding by employing the distance criteria of ¿/(Namino- 'Owater) < 4.5 Â, ¿/(Namino'''Owater) 6.4 Á, and ¿/(Namino'''Owater) <7.2 Á for the first, second and third solvation shells, corresponding to 15, 37 and 51 waters, respectively (Table 1, Figure 1). In addition, we categorized all datasets by assigning class numbers, firstly to the range of Cl -O distances (class attributes ClOl: ¿/(Cl -O) < 3.5 Â; C102: ¿/(Cl -O) > 3.5 Á), and secondly, to the N-H stretching frequencies (class attributes NH1: vnh < 2400 cm"1; NH2: 2700 < vnh > 2400 cm"1; NH3: 3000 < vnh > 2700 cm-1; NH4: vNh > 3000 cm"1). All calculations and preparation of plots were done with the Teach/Me software17 using Teach/Me Data Analysis application.

RESULTS AND DISCUSSION

Review of the Calculated Anharmonic Frequencies and the N-H Stretching Envelope

The calculated fundamental excitation frequencies span a wide range from 3229 cm-1 to 2136 cm-1, being a consequence of very diverse snapshot geometries used to calculate proton potentials of different shapes.9 This Gaussian-type frequency distribution represents the N-H stretching envelope (Figure 2b) having the maximum absorption between 2820 cm-1 and 2730 cm-1, and the center of gravity at 2799 cm-1. According to the correlation scheme of Novak,18 the latter frequency corresponds to the Namino-' Owater distance of about 2.75 Â, which is around 0.1 Á shorter than the calculated CPMD-averaged value of 2.85 Â.9 However, it should be noted that in contrast to the dynamic aqueous environment of the present system, Novak and co-workers established their scheme on a dataset of solid crystalline structures, hence this offset is likely to be of little significance. Another possible reason for a disagreement is that the applied DFT method introduced systematic error in the calculation of proton potentials that is reflected in shifted frequencies. Nevertheless, our computed results are in excellent agreement with the experimental spectra in terms of both the shape and position of the N-H stretching band.9

Multiple Linear Regression (MLR)

In the first stage, vnh frequencies were visualized in a simple 2D linear regression with only Namino' * -Owater distances (Figure 2c), which yielded a very poor correlation (R2 = 0.42), but still considerably better than found for the strong hydrogen bonding in crystalline sodium hydrogen bis(sulfate) (R2 ~ 0.2 at best). Since we were particularly interested in whether the hydrogen bond donor-acceptor distance is a better predictor of vnh than other geometry parameters, at first, we applied the MLR technique to all four Namino"'Owater moieties solvated by four directly bonded water molecules (Figure 2a). Since the mentioned fragments are located in different environments, for example relative to the position of the Cl counterion, the obtained results differed, which led us to apply the MLR analysis to larger matrices of different sizes (see later). One way to express the overall predictive accuracy of multiple regression models is the R2 value,19 a measure of how much variation a model explains. With short Namino' "Owater distances, the model is quite reasonable yielding R2 values above 0.75 (Table 2), and for the hydrogen bonds labeled HB1-HB3 (Figure 2a) it provides significant correlation between geometry parameters and vnh, while for the last one (labeled HB4) we found a weak correlation between data. Interestingly, for each of the four Namino' "Owater moieties one should consider different predictors to describe the corresponding vnh values in the best way. In other words, we did not find a single set of predictors that would work well with all four N-H groups. In addition, it turned out that the data in the full matrix were distinguished according to the CL- Owater distances (see Computational Details), which led us to apply the MLR technique on smaller matrices with 17 x m (ClOl) and 25 x m (C102) dimensions (Figure 3), which increased the R2 values in both cases. From the data we can conclude that the length of the Ca-Namino distance is the most influential predictor in both the protonated -NH* and ring-amino group, labeled as d7 and dl in Figure 3, respectively.

Secondly, we upgraded our model by considering the coordinates of 15 water molecules from the first hydrated shell all with ¿/(Namino' "Owater) < 4.5 Â (Table 2 and Figure 4). It follows that the water molecule labeled with Y is a useful predictor of only one Namino'' Owater moiety, while water labeled X is not listed as predictor in studied cases (Figure 4b). For all four H-bonded moieties the common predictors within the histamine monocation are d2 and d7, together with dl9 from the bulk. Further improvement was made by including more distant second hydration shell solvent molecules in the analysis (¿/(Namino-"Owater) < 6.4 Â) yielding a cluster of 37 waters. Again, CL- Owater distances turned out to be good predictors, together with the majority of other predictors that are all corresponding to bulk water molecules. From the histamine monocation, one should consider d7, listed in HB1, HB2 and HB3, and dl listed in HB2 and HB3 moieties. It is interesting that no predictors are listed in HB4. Finally, third hydration shell (¿/(Namino" 'Owater) < 7.2 Â) with 51 water molecules was also considered. In this case, CL- Owater distances are not good predictors, whereas the majority of predictors are geometry parameters of bulk water molecules. It is interesting that, the predictor d8 (Namino "* Owater distance) can be found in HB1, HB2 and HB4, and not in HB3 moiety. Moreover, the important predictors from histamine monocation are again d7 in HB1, and d6 and d3 in HB4 group.

Summarizing this section, we can conclude that the linear correlation between geometry parameters and individual frequencies exists. It is also worth stressing out, that by increasing the number of solvent water molecules, the features related to the complexity of the proton motion and its coupling to the environment become more apparent.

Principal Component Analysis (PCA)

PCA was performed in order to get an insight into the overall correlativity of all 168 vnh values with the matching geometry parameters. Similarly to the MLR, PCA was applied for each H-bonded moiety HB1-HB4 on matrices composed of 42 rows and sets of 30, 60, 129, and 171 variables in columns (M1-M4), the latter corresponding to the number of considered internal bond distances in each case. The data was additionally preprocessed by using mean centered data, in other words, by subtracting the overall mean value from individual elements.

PCA analyzes data by yielding eigenvectors (principal compoments, PCs) and the corresponding eigenvalues (variances) by using the covariance matrix of the whole dataset. The majority of the information is usually gathered in the first few PCs; in our case, the smallest matrix M1 already exhibits such behavior as already around 98 % variance is explained in the first two principal components in all four Niimmil* * * Owater moieties (Table 3), with the most influencing variables being distances d10, d11, d12, and d13 (Figure 5). We explored a number of pair-wise plots of principal components and compared the corresponding score vs. score displays, and additionally categorized all datasets by assigning class numbers according to a range of both Cl *Owilter distances (ClO1 and ClO2) and N-H stretching frequencies (NH1-NH4) (see Supplementary Materials). Interestingly, in M1, apart from the obvious Namino* Owater distances, one of the most influencing parameters is Cl*- Owater distance (d12), in agreement with the MLR analysis presented earlier, which is basically responsible that the results clusters in ClO1 (having 25 structures) and ClO2 (having 17 structures) classes, with the prevalence of each determined by d9 and d14 variables, respectively (Table 3).

By increasing the number of water molecules in the model one also increases the number of the most influencing geometry parameters while diminishing the significance of the two classes based on the counterion position. In the case of the largest model M4 (51 water molecules), almost all significant geometry parameters correspond to bulk water molecules. The only exception is provided by d8, which is related to the Niimmil* * * Owater moiety.

Further, we were interested in how changing the number of variables affects variances in the first ten PC axes, and found that it exerts pronounced influence (Table 3). By increasing the number of variables, the content of variances decreases, while the percentage of the first three variances in the resulting PCs changes only distinctly. This could be ascribed to the CL counterion, which has a more profound effect in the smaller models, and to the complexity of the proton motion and its coupling to its environment. On the other hand, when comparing percentages of the first three variances within similar models (e.g. Ml model with the first shell waters) they are almost identical (Table 3). It follows that with PCA method we obtained similar conclusion as from the MLR study - instantaneous fluctuating structure has a pronounced effect on the proton motion and its dynamic, which is very complex and differently coupled to the environment.

Along with other established computational methods, the presented approach will contribute towards understanding of the complexity of the hydrogen bonding dynamics in polar environments reflected in the vibrational spectra20"24 and its relevance for biocatalysis.25

CONCLUDING REMARKS

The aim of our study was to evaluate the correlation between the frequencies of the N-H stretching band and the corresponding geometry parameters of the snapshot structures of the histamine monocation in aqueous solution. A simple two-dimensional correlation scheme between vnh and the matching Namin0-* Owater distances appeared not to have significant predictive ability. On the other hand, more advanced techniques, such as multiple linear regression (MLR) and principal component analysis (PCA), applied on preprocessed data yielded reasonably accurate predictors. For example, around 98 % variance in the data is well explained by the first two principal components in all four Namin0- * - Owater moieties. Moreover, with PCA method applied on the cluster involving only four directly attached water molecules we were able to distinguish two classes of structures based on the counterion position. When the analysis was additionally tuned by assigning class numbers according to the value of the CL- Owater distances, the percentage of the variance in the first two principal components in all Namino - * 'Owater moieties decreases to around 90 %.

By increasing the number of water molecules in the model one also increases the number of the most influencing geometry parameters and diminishes the relevance of the mentioned two structural classes that depend on the position of the (X counterion. With the largest model investigated here (51 water molecules), all influential geometry parameters arise from bulk water, and only around 20 % of the variance is explained in the first two principal components of four Nammo - Owater groups. This might lead to a conclusion that a correlation between instantaneous fluctuating histamine structure and the N-H stretching frequencies is likely to be poor, and that the dynamic of the proton motion is very complex and strongly coupled to the fluctuating environment. The results obtained from the MLR analysis support this assumption.

Despite the fact that the structure-frequency correlation found for a set of dynamically sampled instantaneous snapshot structures within one system is only moderate, it represents a clear advancement in comparison with the previously investigated system of crystalline sodium hydrogen bis(sulfate),14 and demonstrates once again that hydrogen bonding is an extremely complex interaction. For the mentioned inorganic material, no appreciable linear or even non-linear correlation could be found, indicating the limitations of the applied one-dimensional quantum treatment. Remarkably, the present system features much weaker and longer hydrogen bonds, which is reflected in the observed correlations between the instantaneous geometric parameters and the corresponding N-H stretching frequencies.

Supplementary Materials. - Supporting informations to the paper, containing selected pair-wise plots of principal components PC1 to. PC2 together with the corresponding score to. score displays, are enclosed to the electronic version of the article. These data can be found on the website of Croatica Chemica Acta (http://public.carnet.hr/ccacaa).

Acknowledgements. G. P., .1. S. and .1. M. would like to thank the Slovenian Research Agency for financial support in the framework of the program group P1-0012 and within the corresponding research project contract No. .11-2014. R. V. gratefully acknowledges the European Commission for an individual FP7 Marie Curie Career Integration Grant (contract number PCIG12-GA-2012-334493). Part of this work was supported by COST Action CM1103.

Footnote

[dagger] to Dr. Mirjana Eckert-Maksic on the occasion of her 70th birthday.

References

REFERENCES

1. K. Yanai and M. Tashiro, Pharmacol. Ther. 113 (2007) 1-15.

2. A. Merickel and H. Edwards, Neuropharmacology 34 (1995) 1543-1547.

3. P. Gilli, L. Pretto, V. Bertolasi, and G. Gilli, Acc. Chem. Res. 42 (2009) 33-44.

4. D. G. Cooper, R. C. Young, G. J. Durant, and C. R. Ganelin, in: Comprehensive medicinal chemistry: the rational design, mechanistic study and therapeutic application of chemical compounds, C. Hansch (ed.), Pergamon Press, Oxford, U.K., 1990, 323-421.

5. W. L. Nelson, in: Foye's Principles of Medicinal Chemistry. T. L. Lemke, D. A. Williams, V. F. Roche, and S. W. Zito (eds.), 7th edition, Baltimore, Philadelphia: Lippincott Williams & Wilkins, 2013.

6. S. J. Hill, C. R. Ganellin, H. Timmerman, J. C. Schwartz, N. P. Shankley, J. M. Young, W. Schunack, R. Levi, and H. L. Haas, Pharmacol. Rev. 49 (1997) 253-278.

7. B. Jójárt, R. Kiss, B. Viskolcz, G. M. Keseru, J. Chem. Inf. Model. 48 (2008) 1199-1210, and references cited therein.

8. (a) J. Henin, B. Maigret, M. Tarek, C. Escrieut, D. Fourmy, C. Chipot, Biophys. J. 90 (2006) 1232; (b) L. Rivail, C. Chipot, B. Maigret, I. Bestel, S. Sicsic, M. Tarek, J. Mol. Struct. (Theochem) 817 (2007) 19-26.

9. J. Stare, J. Mavri, J. Grdadolnik, J. Zidar, Z. B. Maksic, and R.Vianello, J. Phys. Chem. B 115 (2011) 5999-6010.

10. R. Vianello and J. Mavri, New J. Chem. 36 (2012) 954-962.

11. R. Blinc, D. Hadzi, and A. Novak, Bericht. Buns. Phys. Chem. 64 (1960) 567.

12. W. Mikenda, J. Mol. Struct. 147 (1986) 1-15.

13. G. Pirc, J. Stare, and J. Mavri, J. Chem. Phys. 132 (2010) 224506.

14. G. Pirc, J. Mavri, M. Novic, and J. Stare, J. Phys. Chem. B 116 (2012) 7221-7231.

15. I. T. Jolliffe, Principal component analysis, 2nd edition, Springer: New York, 2002.

16. Y. Dodge and F. H. C. Marriott, International Statistical Institute. The Oxford dictionary of statistical terms, 6th edition, Oxford University Press: Oxford, UK, 2003.

17. H. Lohninger, Teach/Me Data Lab. 2.002. in: SDL-Software Development Lohninger, Springer: Berlin, 1999.

18. A. Lautié, F. Froment, and A. Novak, Spectrosc. Lett. 9 (1976), 289-299.

19. H. F. Kaiser, Educ. Psychol. Meas. 20 (1960) 141-151.

20. M. Z. Brela, M. J. Wójcik, M. Boczar, and R. Hashim, Chem. Phys. Lett. 558 (2013) 88-92.

21. N. Rekik and M. J. Wójcik, Chem. Phys. 369 (2010) 71-81.

22. M. Boczar, R. Kurczab, and M. J. Wójcik, Vibrat. Spectr. 52 (2010) 39-47.

23. N. Issaoui, N. Rekik, B. Oujia, and M. J. Wójcik, Int. J. Quant. Chem. 110 (2010) 2583-2602.

24. N. Doslic, M. Abdel-Latif, and O. Kuehn, Acta Chim. Slov. 58 (2011) 411-424.

25. S. C. L. Kamerlin, J. Mavri, and A. Warshel, FEBS Lett. 584 (2010) 2759-2766.