CARBO_R Force Field for

Subscriber access provided by UNIV OF YORK

B: Biophysical Chemistry and Biomolecules CARBO/CARBO_R

Extension of the GROMOS 56a6 Force Field for Charged, Protonated and Esterified Uronates Karina Panczyk, Karolina Gaweda, Mateusz Drach, and Wojciech Plazinski J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.7b11548 • Publication Date (Web): 20 Mar 2018 Downloaded from http://pubs.acs.org on March 20, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Extension of the GROMOS 56a6CARBO/CARBO_R Force Field for Charged, Protonated and Esterified Uronates Karina Panczyk1, Karolina Gaweda1, Mateusz Drach2, Wojciech Plazinski1,* 1

Jerzy Haber Institute of Catalysis and Surface Chemistry, Polish Academy of Sciences, ul. Niezapominajek 8, 30-239 Cracow, Poland

2

Department of Theoretical Chemistry, Faculty of Chemistry, M. Curie-Sklodowska University, pl. M. Curie-Sklodowskiej 3, 20-031 Lublin, Poland

Abstract An extension of the GROMOS 56a6CARBO/CARBO_R force field for hexopyranose-based carbohydrates is presented. The additional parameters describe the conformational properties of uronate residues. The three distinct chemical states of the carboxyl group are considered: deprotonated

(negatively charged),

protonated

(neutral) and esterified (neutral).

The

parameterization procedure was based on quantum-chemical calculations and the resulting parameters were tested in the context of: (i) flexibility of the pyranose rings under different pH conditions; (ii) conformation of the glycosidic linkage of the (1→4)-type for uronates with different chemical states of carboxyl moieties; (iii) conformation of the exocyclic (i.e. carboxylate and lactol) moieties; (iv) structure of the Ca2+-linked chain-chain complexes of uronates. The presently proposed parameters in combination with the 56a6CARBO/CARBO_R set can be used to describe the naturally occurring polyuronates, composed either of homogeneous (e.g. glucuronans) or heterogeneous (e.g. pectins, alginates) segments. The results of simulations relying on the new set of parameters indicate that the conformation of glycosidic linkage is nearly unaffected by the chemical state of the carboxyl group, in contrary to the ring conformational equilibria. The calculations for the poly(α-D-galacturonate)-Ca2+ and poly(α-L-guluronate)-Ca2+ complexes show that both parallel and anitiparallel arrangements of uronate chains are possible but differ in several structural aspects.

*

Corresponding author. Tel.: +48815375685; fax: +48815375685 E-mail addresses: [email protected]; [email protected]

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1. Introduction Polyuronates are linear polymers of uronic acids that contain (1→4)-linkage (see Fig. 1).1 As their molecules bear the carboxyl groups, they are polyelectrolytes able to change the ionization state upon different pH conditions. The three most important natural derivatives of polyuronates are pectins, alginates and glucuronans. They have diverse biological functions in plants, related mainly to the preservation of the structure, texture and flexibility.2–4 Furthermore, they are able to form hydrophilic gels,5–7 typically in the presence of divalent metal cations. Pectins4–9 exhibit a complex sequence involving mainly homogenous segments of (1 → 4)linked α-D-galacturonate (α-D-GalU) residues but may also comprise neutral sugar molecules (e.g. rhamnose, galactose and arabinose) and include branched segments. Furthermore, the carboxyl groups of pectins may be methoxy-esterified to different extents which strongly influences their physico-chemical properties.10,11 Alginates3,5–7,12 are linear copolymers of β-D-mannuronate (β-DManU) and α-L-guluronate (α-L-GulU) residues in (1 → 4)-linkage that predominantly involve homogeneous segments of either types of residues and regular repeats of their dimer. Finally, glucuronans13,14 are linear, homogeneous polymers of (1 → 4)-linked β-D-glucuronate (β-D-GlcU) residues. The most relevant processes involving polyuronates occur in the aqueous environment, therefore the experimental determination of structural and conformational properties of these compounds is significantly hindered. In view of this fact, the methods of molecular modeling are a promising alternative offering insight into desired information at the molecular scale.15–21 However, the accuracy of such methods relies primarily on the quality of the molecular models that represent real physical systems. Although modeling the structure/conformation of uronates as well as their interactions with metal ions has involved a wide spectrum of theoretical methods from rigid docking22 to QM-MM-based simulations,18 the large dimensions of the typical, polyuronate-containing system make the molecular dynamics (MD) simulations performed within the classical force fields the most promising choice. The non-GROMOS carbohydrate-dedicated force fields applied to simulate the uronate molecules include: AMBER14,23 GLYCAM0624 (with the extensions proposed by Dentini et al.25 or Hackbusch et al.26), CHARMM,27 OPLS-AA,2 CVFF28 and (CHARMM-based) CSFF.29 The main field of their application is connected mainly with investigating the structural and conformational properties of naturally-occurring mono- (e.g. iduronic acid and its derivatives) and polyuronates (e.g. alginates or pectins) as well as their interactions with divalent metal ions.4,21,30– 36


Page 2 of 37

Page 3 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The most recent version of the carbohydrate-dedicated GROMOS force field is the 56a6CARBO set37 which in combination with the correction named 56a6CARBO_R38 describes the nonfunctionalized or O-alkylated pyranoses as well as di-, oligo- and polymers composed of (1→1)-, (1→2)-, (1→3)-, (1→4)- and (1→6)-linked residues. So far, no other parameters extending the applicability of this force field in the context of the functionalized carbohydrates have been proposed. Although the GROMOS force field exhibits a large degree of transferability,38–42 the direct adaptation of parameters for various functional groups from the existing sets may be problematic, especially in reference to the ring-inversion properties which are extremely sensitive to alterations in the substituent-substituent and ring-substituent interactions.38 The approach relying on the direct adaptation of the GROMOS-inherent parameters for carboxylate residues has been tested in the past34,35 in the case of the GROMOS 45a4 force field41 but without extensive, quantitative validation against the ab initio or experimental data. In the present article we introduce the parameters that extend the functionality of the GROMOS 56a6CARBO and 56a6CARBO_R sets for the general case of uronates, i.e. aldohexopyranoses substituted by the carboxyl group at position C5 (see Fig. 1 for atom numbering). The derived parameters take into account the possible different chemical states of the carboxyl moiety which can be deprotonated, protonated or esterified. Although the parameterization was focused mainly on the most relevant uronates (e.g. α-D-GalU, β-D-GlcU, αL-Gul)

the final set can be applied to any other uronate residue. Furthermore, because of

adaptation of the already existing parameters, the current set can be used to study the longer uronate chains as well as the chains of non-uniform composition, containing either different uronate residues (e.g. alginates) or both uronate and non-uronate residues (e.g. pectins). It is worth mentioning that the direct combination of the existing GROMOS parameters (the 53a6 set to describe the interactions involving carboxyl moieties and the 56a6CARBO/CARBO_R set to describe the rest of carbohydrate molecule) may not be a viable procedure when considering several different conformational descriptors. More precisely, a preliminary study has shown that such combination leads to inaccurate energy profiles associated with the rotation of a carboxyl group (see Fig. S1 in the Supporting Information). Furthermore, the ring-inversion free energies calculated for esterified and protonated uronates are often in disagreement with the available structural data. The newly-developed and validated set of parameters is used in order to investigate several issues related to the structure/conformation of uronate mono- and oligomers as well as to the interactions of natural uronate polymers (alginates and pectins) with calcium ions.



Fig. 1. (A) Atom numbering used in the article for the three different chemical states of carboxyl group of uronates. (B) Structural formulae of the uronate residues comprising the natural polyuronates (alginates, pectins and glucuronans) and the monomer of α-L-iduronate. Aliphatic hydrogen atoms (treated implicitly in the GROMOS force field) are not shown.


Page 4 of 37

Page 5 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 2. Comparison of energy profiles calculated at the classical (the presently proposed parameters) and QM (DFT/B3LYP/6-311++G**) levels of accuracy. Panels A−C present the energy profiles for the rotation of the anionic carboxyl groups of β-D-glucuronate, α-D-galacturonate and α-L-guluronate, respectively. Panels D and E present the profiles calculated for the protonated carboxyl group, whereas panels F and G for the esterified carboxyl group of α-D-galacturonate, respectively. Panel H shows the profiles for the protonated carboxyl group of β-D-glucuronate. Quadruplets of atoms that define considered dihedral angle are given in each of the panels. The QM- and force field-derived values are given in red and green, respectively. ACS Paragon Plus Environment


Tab. 1. Atom types and atomic partial charges for the carboxyl group in uronates*. See Fig. 1 for atom numbering. Atom

Atom type

Partial atomic charge [e-]

-COO- group C6

C

+0.284

O61

OM

-0.642

O62

OM

-0.642 -COOH group

C6

C

+0.874

O61

OM

-0.642

O62

OA

-0.642

H62

H

+0.410 -COOCH3 group

C6

C

+0.874

O61

O

-0.642

O62

OE

-0.464

C7

CH3

+0.232

* charge group is created by all atoms belonging to any of the above moieties.


Page 6 of 37

Page 7 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Tab. 2. Bond stretching, bond-angle bending, and improper-dihedral deformation parameters for the carboxyl group in uronates.* See Fig. 1 for atom numbering. Bond type

kb [106 kJ/mol/nm4]a

b0 [nm]

C5-C6

7.15

0.153

–COO- group C6-O61, C6-O62

13.4

0.125

–COOH group C6-O61

16.6

0.123

C6-O62

10.2

0.136

O62-H6

15.7

0.1

–COOCH3 group C6-O61

16.6

0.123

C6-O62

10.2

0.136

O62-C7

8.18

0.143

Bond-angle type

kθ [kJ/mol]b

θ0 [deg]

–COO- group C5-C6-O61, C5-C6-O62

635.0

117.0

O61-C6-O62

770.0

126.0

–COOH group C5-C6-O61

685

121

C5-C6-O62

610

115

O61-C6-O62

730

124

C6-O62-H6

450

109.5

–COOCH3 group C5-C6-O61

685

121

C5-C6-O62

610

115

O61-C6-O62

730

124

C6-O62-C7

380

109.5

Improper dihedral type

kζ [kJ/mol/deg2]c

ζ0 [deg]

C6-O5-C4-C5

0.051

35.264

C6-O61-O62-C5

0.102

0

* functional forms of the corresponding potentials are given in ref.37 by eqs. (1)-(3).



Page 8 of 37

Tab. 3. Torsional dihedral parameters for the carboxyl group in uronates*. See Fig. 1 for atom numbering. Torsional dihedral

kφ [kJ/mol]b

φ0 [deg]

m

120

2

type

–COO- group C4-C5-C6-O61

3.90

–COOH group C4-C5-C6-O61

3.90

120

2

C4-C5-C6-O61

3.50

−60

1

C5-C6-O62-H62

16.7

180

2

–COOCH3 group C4-C5-C6-O61

3.90

120

2

C4-C5-C6-O61

9.50

−60

1

C5-C6-O62-C7

16.7

180

2

* The functional form of the corresponding potential is given by eq. (1).

2. Methods 2.1 Parameterization strategy The parameters describing the interactions that do not involve carboxyl moiety were transferred from the GROMOS 56a6CARBO force field37,38 without any changes. This includes the 56a6CARBO_R correction38 which is relevant for non-monomeric residues as well as for O1-alkylated residues (see Fig. 1 for atom numbering). The parameters present in the original force field and involving the hydroxymethyl moieties were removed and replaced by the new parameters, related to the carboxyl moiety. For the clarity, we do not describe in details every attempt made in order to develop the parameter set but limit ourselves to the final procedure that has lead to successful results. The parameterization strategy relied on the following consecutive steps: (i) assigning atomic charges for atoms present in the carboxyl (–COO-), protonated carboxyl (–COOH) of esterified carboxyl (–COOCH3) groups; (ii) assigning the non-bonded (Lennard-Jones) interaction parameters for atoms present in –COO-, –COOH and –COOCH3; (iii) assigning the bonded (bond stretching, bond angle bending, dihedral angle motion, improper angle bending) interaction parameters involving –COO-, –COOH and –COOCH3; (iv) validation procedure performed by comparison of the energy profiles with the ab initio-originating data and by the analysis of the MD simulation results. Partial atomic charges were assigned on the basis of quantum-mechanical (QM) calculations at the DFT/B3LYP/6-311++G** level of theory43–47 and the subsequent RESP48 fitting of the resulting electrostatic potentials. The charge calculations involved α-D-galacturonate and β-DACS Paragon Plus Environment

Page 9 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


glucuronate with carboxylate groups in three different chemical states (see above). The RESP procedure involved the constraints on the allowed atomic charges which were introduced in order to: (i) maintain the charges on non-carboxyl atoms unchanged with respect to the original set; (ii) keep the charge groups neutral or of integer charge. The final atomic charges are reported in Tab. 1. The comparison with the charges inherent to the 53a6 set39 reveals that the differences in the case of deprotonated carboxyl groups are rather minor (-0.635 vs. -0.642 e and +0.270 vs. -0.284 e for the O6 and C6 atoms, respectively). However, both protonation and methylation of the O62 atom lead to a significant charge flux concerning mainly the C6 atom and large divergences from the default parameters present in the 53a6 set.39 For instance, the charges on C6 characteristic of the 53a6 set vary between +0.253 e (methylated carboxyl) and +0.33 e (protonated carboxyl). The current, uronate-oriented calculations predict much higher values (+0.874 e), closer to those inherent to the 53a6OXY set.49 This divergence has a large impact on the predicted energies associated with the rotation of the carboxyl group around the -C5-C6- bond (see Fig. S1 in the Supporting Information). The Lenard-Jones (LJ) parameters were unchanged with respect to the GROMOS 53a6/56a6CARBO/56a6CARBO_R sets of parameters. The LJ constants for atoms present in the newlyparameterized groups were taken from the 53a6 set and involve OM, OA, OE, CH3 and C atom types. Note that the original GROMOS 56a6CARBO force field contains the ‘exceptional’ LJ pairs which may involve the 1-4 and 1-5 neighbors. In the case of atoms not belonging to the carboxyl group, these type of interactions were kept unchanged with respect to the original set.37,38 The interactions involving C6 were also kept unmodified (the Hassel–Ottar repulsion terms between the C6-O1 and C6-O3 atom pairs as well as substituent-to-ring repulsion terms between the C6-C1 and C6-C3 atom pairs) with respect to 56a6CARBO. Note that this involves the incorporation of the additional atom type in the topological patterns N3 and N4 defined in ref.37 (see Tab. 5 therein). The potentials for bond stretching, bond-angle bending and improper dihedral deformation terms were retrieved directly from the GROMOS 53a6/56a6CARBO sets. The functional forms of these potentials are given elsewhere37,39–41 and the corresponding parameters [force constants (kb, kθ and kζ, respectively), optimal values of bond lengths and angles (b0, θ0 and ζ0, respectively)] are shown in Tab. 2. These parameters were obtained on the basis of the experimental spectroscopic (force constants) and X-ray diffraction (optimal distance and angle values) data for small molecules. The optimal bond-lengths and bond-angle values agree well with the results of our QM calculations (geometry optimization followed by the RESP charge fitting; see above). The torsional potentials, associated with the rotation of the particular elements of –COO-, –COOH and –COOCH3 were adjusted on the basis of the QM-derived values obtained due to relaxed scans of energy in the function of changing values of dihedral angles defined by the following quadruplets ACS Paragon Plus Environment


Page 10 of 37

of atoms: C4-C5-C6-O61 (anionic β-D-guluronate, α-D-galacturonate and α-L-guluronate), C5-C6O62-H62 (protonated α-D-galacturonate and β-D-guluronate) and C5-C6-O62-C7 (esterified α-Dgalacturonate). The functional form of the torsional potential is given by the following equation: ,

(1)

where ϕ is the dihedral angle value, m the multiplicity of the term, φ0 the associated phase shift, and kϕ the corresponding force constant. The values of m, φ0 and kϕ were adjusted on the basis of the QM data as shown in Fig. 2; the final parameters are given in Tab. 3. It is worth noting that a given dihedral angle may be involved in more than one torsional potential energy term with different multiplicities and/or phase shifts. Interestingly, and contrary to the default GROMOS parameters for carboxyl group (see e.g. ref.),34 the best adjustment of the QM data related to the rotation of carboxyl group is achieved when not assuming the single, 6-fold multiplicity but either the single, 2-fold potential or the combination of the two (1-fold and 2-fold) potentials. The performance of the direct combination of the charges and torsional potentials taken from the 53a6 set with those of 56a6CARBO/CARBO_R was tested in the context of energy profiles associated with the rotation of the carboxyl groups and their agreement with the QM data (analogously to the data shown in Fig. 2). The results are shown in Fig. S1 (Supporting Information). The QM calculations were performed in Gaussian0950 at the DFT/B3LYP/6-311++G** level of theory. The starting structures were drawn in Avogadro 1.1.151 and preoptimized at the UFF level of theory;52 they reflected the ring conformation characteristic of the series (i.e. 4C1 for D-saccharides

and 1C4 for L-guluronate). The orientation of the carboxyl group was adjusted

according to the most favorable rotamer found during the energy scans. The orientation of remaining groups (exocylic hydroxyl moieties) was adjusted in the positions that correspond to the lack of intramolecular hydrogen bonding (according to the findings described in ref.53).

2.2 Simulation details The compounds under consideration included: (i) monomers of all uronates of the

D

series; (ii)

monomers of O1-methylated α-L-iduronate and O1-methylated β-D-glucuronate in their deprotonated forms; (iii) homooctamers of (1→4)-linked residues of β-D-glucuronate (β-D-GlcU), β-D-mannuronate (β-D-ManU), α-L-guluronate (α-L-GulU) and α-D-galacturonate (α-D-GalU); (iv) parallel- and antiparallel-oriented molecular complexes composed of the two oligomers (14 units) of α-L-guluronate (in the anionic state) associated with the contribution of seven calcium ions; (v) parallel- and antiparallel-oriented molecular complexes composed of the two oligomers (14 units) of α-D-galacturonate (in the anionic state) associated with contribution of seven calcium ACS Paragon Plus Environment

Page 11 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ions. For compounds listed in points (i) and (iii) both their anionic and protonated (the latter denoted by subscript ‘H’) forms were taken into account; additionally, for α-D-GalU-containing compounds, the esterified forms were considered (denoted by subscript ‘Me’). Fig. 1 shows the chemical structures of the studied compounds or their monomeric units. The initial structures of complexes mentioned in points (iv) and (v) were based on the results from our previous studies15 and on ref.22 All MD simulations were carried out with the GROMACS 4.5.5 package54 in combination with PLUMED 1.355 (metadynamics) or GROMACS 5.056 with PLUMED 2.1.557 (unbiased MD). The parameters used for simulations are described in the previous section. Molecular systems under consideration consisted of one carbohydrate molecule and approximately 9001600 simple point charge (SPC) water molecules58 placed in a cubic (30 × 30 × 30 Å3 or 36 × 36 × 36 Å3) simulation box. The systems of smaller or larger dimensions corresponded to simulations of monomers or homooctamers, respectively. Additionally, we studied the complexes composed of the two uronate chains and calcium ions. The initial geometry of such complexes was taken from our previous studies15 and the number of water molecules (about 19160) as well as the edges of simulation cell (97 × 97 × 97 Å3) were increased to account for the larger dimensions of related systems. The unbiased simulations were carried out under periodic boundary conditions and in the isothermal-isobaric ensemble (unbiased MD) or at fixed volume (metadynamics simulations). The temperature was maintained close to its reference value (298 K) by applying the V-rescale thermostat,59 whereas for the constant pressure (1 bar, semi-isotropic coordinate scaling) the Parrinello-Rahman barostat60 was used with a relaxation time of 0.4 ps. The equations of motion were integrated with a time step of 2 fs using the leap-frog scheme.61 The solute bond lengths were constrained by application of the LINCS procedure62,63 with a relative geometric tolerance of 10-4. The full rigidity of the water molecules was enforced by application of the SETTLE procedure.58 The translational center-of-mass motion was removed every timestep separately for the solute and the solvent. The non-bonded interactions were calculated using a twin-range scheme,64 with short- and long-range cutoff distances set to 0.8 and 1.4 nm, respectively, and an update frequency of 5 timesteps for the short-range pairlist and intermediate-range interactions. The reaction-field correction65,66 was applied to account for the mean effect of the electrostatic interactions beyond the long-range cut off distance, using a relative dielectric permittivity of 61 as appropriate for the SPC water model.67 All the systems were preoptimized by a 0.5-1.5 ns constant-pressure MD equilibration at 1 bar and 298 K, ensuring an effective solvent density appropriate for these conditions in the subsequent production simulations. After equilibration, all unbiased simulations were carried out ACS Paragon Plus Environment


for 200 ns and the trajectory was saved every 10 ps. Apart from the regular GROMACS trajectory, the values of selected parameters were saved every 0.2 ps. These parameters included: the γ dihedral angle, describing the rotation of carboxyl group (defined by the O61-C6-C5-C4 quadruplet of atoms) and the φ dihedral angle, describing the rotation of lactol group (defined by the O5-C1-O1-H1 quadruplet of atoms). The distributions of dihedral angles were calculated on the basis of unbiased MD simulations initiated from the structures representing the energetically favourable chair conformer. To study the free energy profiles (FEPs) as a function of the Cremer-Pople θ puckering coordinate68 the parallel tempering69 scheme was combined70 with the well-tempered metadynamics approach71 as implemented in PLUMED 1.3. The well-tempered metadynamics relied on Gaussian local functions of widths 2.86 deg, an initial deposition rate of 0.01 kJ/mol/ps and a temperature parameter ∆T (see eq. (2) in ref.71) of 1788 K. The parallel-tempering69 relied on 16 metadynamics simulations carried out in parallel at different temperatures ranging from 298.0 to 363.2 K in steps of about 4.3 K, along with replica-exchange attempts performed at 2 ps intervals. For the FEP analyses along θ, the intervals (0-60 deg) and (120-180 deg) were ascribed to regular chair (4C1) and inverted chair (1C4) conformations, respectively. The free-energy difference associated with the chair-chair (4C1→1C4) distortion was defined as the difference between the values at the minima of FEP in the two corresponding regions and denoted as ∆Fi. The convergence of the final ∆Fi values was checked by the hand-written bash script and the PLUMED in-build sum_hills tool. The convergence plots are given in Fig. S2 (Supporting Information). The conformation of glycosidic linkages was studied in terms of 2D free energy maps (FEMs) calculated by using the combination of parallel tempering and well-tempered metadynamics, analogously to the procedure described above. The φ and ψ dihedral angles defining glycosidic linkages are created by the O5-C1-O1-C'4 and C1-O1-C'4-C'3 quadruplets of atoms, respectively. The simulations relied on Gaussian local functions of widths 18 deg and an initial deposition rate of 0.01 kJ/mol/ps. The 3JH,H coupling constants for vicinal protons were calculated by using the AltonaHaasnoot relation72–74 and the unbiased MD trajectories initiated from both the 1C4 and 4C1 conformers (Me-O1-α-L-IdoU and α-L-IdoU) or only from the 4C1 one (Me-O1-β-D-GlcU). The average 3JH,H value for Me-O1-α-L-IdoU and α-L-IdoU was calculated as weighted average of both 4

C1 and 1C4 contributions with the weight assigned on the basis of the separately calculated ∆Fi

value.


Page 12 of 37

Page 13 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3. Results and Discussion 3.1 Flexibility of the pyranose rings The chair-chair ring-inversion free energies (∆Fi) are given in Tab. 4 while some of their values calculated for monomers are graphically illustrated in Fig. 3. Regarding the ring-inversion properties of the deprotonated uronate monomers, the sign of the calculated ∆Fi values shows that the dominant conformation for 13 (out of 16) compounds is a regular chair. The remaining three (α-D-AllU, α-D-IdoU and α-D-AltU) exhibit the inverted chair geometry of the ring. The close-to-zero values of ∆Fi (|∆Fi| < 5 kJ/mol) for α-D-IdoU, α-D-AllU, α-D-ManU and α-D-GulU confirm the non-negligible populations of both 1C4 and 4C1 conformers. The protonation of the carboxyl moiety leads to a notable increase of the ring-inversion energies; this increase is the largest for those α-anomers which exhibited the lowest values of ∆Fi in the form of deprotonated compounds. This effect results in positive ∆Fi values for all neutral uronates and significant reduction of the population of non-4C1 conformers for the whole group of uronates as the lowest ∆Fi value equals now +8.4 kJ/mol. Therefore, our results indicate that pH-dependent protonation state may have an influence on the ring-inversion properties of uronates but this effect is limited to the α-anomers of compounds which exhibit most flexible rings. However, this group contains several among the most biologically-relevant and naturally abundant uronates (e.g. iduronate, mannuronate and guluronate).

Fig. 3. The ring-inversion free energies (∆Fi) calculated within the current set of parameters for both deprotonated and protonated uronates, compared to the corresponding values for unfunctionalized aldohexopyranoses (data taken from ref.38). ACS Paragon Plus Environment


Page 14 of 37

Tab. 4. The ring-inversion (4C1→1C4) free energies (∆Fi) in kJ/mol for the monomeric uronates and the uronate residues in a homooctameric chain, calculated using the presently proposed extension to the GROMOS 56a6CARBO force-field. The calculated values rely on the Cremer-Pople coordinates68 for the state assignment. The part of the data is illustrated graphically in Fig. 3. Monomers

α-anomers -

(-COO )

Residue in a chain

β-anomers -

(-COO )

Compound

α-anomers

β-anomers

(-COOH)

(-COOH)

Compound

∆Fi

∆Fi

D-AllU

-2.7

18.0

17.5

20.2

α-D-GalU

18.1

D-AltU

-13.5

8.1

8.9

14.5

α-D-GalUH

13.2

D-GlcU

9.6

25.5

18.7

22.3

β-D-GlcU

7.6

D-ManU

0.5

16.6

8.4

19.5

β-D-GlcUH

1.2

D-GulU

5.7

25.3

17.9

24.1

β-D-ManU

1.4

D-IdoU

-4.9

14.5

8.6

16.1

β-D-ManUH

-1.5

D-GalU

16.0

32.3

20.3

28.3

α-L-GulU

-11.4

D-TalU

14.1

21.9

12.7

19.9

α-L-GulUH

-22.7

O1-Me-L-IdoU

-4.1

-

-

-

-

O1-Me-D-GlcU

-

17.1

-

-

Interestingly, the population of the non-chair conformers (i.e. boat and skew-boat ring geometries) is always much smaller in comparison to both the 4C1 and 1C4 ones and such effect is independent of the protonation state. This means that such ring geometries do not contribute much into the overall conformational space of uronate monomers and the ring conformational equilibria in some very flexible uronates (e.g. IdoU and its derivatives)75–77 can be described by taking into account only the 4C1 and 1C4 shapes. The similar findings have been reported in refs.78,79 for unfunctionalized hexopyranoses and hexopyranose residues in a chain. When comparing the ∆Fi values calculated for uronate monomers with the analogous quantities obtained for corresponding unfunctionalized aldohexopyranose monomers (data taken from ref.38) one obtains a clear correlation between these two sets of data, especially in the case of protonated uronates. This correlation is slightly less for deprotonated uronates which (mostly in the case of the α-anomers) exhibit the lower ∆Fi values. The existence of such correlation between chair-chair inversion energies can be expected on rational grounds, mainly due to the additivity of the interactions between ring substituents. This effect has been confirmed experimentally and became a foundation for proposing semi-empirical schemes capable of predicting the ∆Fi values


Page 15 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


for various pyranoses.80–82 A similar additivity effect, resulting in a correlation between ∆Fi values, has been confirmed to occur for various groups of pyranoses.83 Additionally, we have calculated the ring-inversion properties for residues in a chain for selected types of homooctamers. The results are given in Tab. 4. Although in the majority of the cases the favorable ring shape is the same as for the corresponding monomer, some of residues in a chain are very prone to ring inversions (e.g. mannuronates and protonated glucuronate). This is also the case of some of unfunctionalized hexopyranoses in a chain (see Tab. 3 in ref.38) and seems to be an inherent feature of the GROMOS 56a6CARBO_R force field but is not necessarily in disagreement with the real conformational behavior of the polyuronates in aqueous solutions (note that there is no experimental data regarding the ring shapes in polyuronate molecules in solutions). Interestingly, the influence of the protonation state seems to be opposite to that observed in the case of monomeric uronates (Fig. 3) but due to the fact of limited number of considered oligomers, this may be ascribed to the decorrelation of the ring-inversion properties of residues in a chain relative to monomers. Similarly to the results described elsewhere38,53, also here we have not been able to distinguish any single interaction pattern responsible for this effect and its magnitude.

Fig. 4. Experimental (yellow bars) and calculated vicinal coupling constants for O1-methylated α-Liduronate and β-D-glucuronate (red bars) as well as for unfunctionalized α-L-iduronate (blue bars). The experimental data were taken from ref.75 whereas the Haasnoot-Altona equation72 was used to obtain the theoretical values.

To our knowledge, no quantitative data exist on the ring-inversion properties of any uronate monomers except of O1-methylated β-D-GlcU and α-L-IdoU.75 Ring three-bond H-H vicinal coupling constant (3JH,H) measurements revealed significant differences between the ring geometries of these two compounds. The larger values of 3JH,H for β-D-GlcU speak for the strongly dominating 4C1 conformation whereas the lower ones for α-L-IdoU suggest the ACS Paragon Plus Environment


Page 16 of 37

significant contribution of non-4C1 shapes. In ref.75 the latter finding has been interpreted as the existence of the equilibrium mixture of chair and skew-boat conformers. The results of our simulations offer an alternative interpretation based rather on the dynamic equilibrium between 1

C4 and 4C1 conformers with the negligible contributions of any boat or skew-boat conformers.

Fig. 4 shows the experimental data compared with the results of our MD simulations. The ∆Fi value calculated for Me-O1-α-L-IdoU is equal to -4.1 kJ/mol which suggests the domination of the 1

C4 conformation with a non-negligible (~16 %) contribution of the 4C1 conformer. The chair-

boat/skew-boat inversion is associated with much higher free energy (> 10 kJ/mol) and the corresponding contribution is neglected. The agreement with the experimental data is almost ideal in the case of β-D-GlcU and satisfactory in the case of α-L-IdoU. In the later case, the agreement is of comparable quality in comparison with the predictions of the GLYCAM force field24 according to Fig. 2 in ref.75 According to the predictions of the GROMOS 56a6CARBO_R force field, methylation at O1 greatly influences the ring-inversion properties of pyranose monomers in comparison to the unfunctionalized ones. More precisely, the Me-O1 substitution results in a systematic shift (~7.5 kJ/mol) of ring-inversion free energies (4C1 → 1C4) which is positive for the α and negative for βanomers of

D-hexopyranoses

(the trend is opposite for

L-hexopyranoses).

However, recent

studies53 indicate that the magnitude of such shift may be smaller than previously estimated. Due to this observation, we additionally calculated the 3JH,H constants for unfunctionalized α-L-IdoU which exhibits the ∆Fi value equal to ~5 kJ/mol. The results given in Fig. 4 are also in a good agreement with the experimental data (although the direction of deviations is opposite in comparison to the results obtained for Me-O1-α-L-IdoU). On the basis of the above experimental 3JH,H values related to α-L-IdoU and the two-state model one can calculate the theoretical proportion of the 4C1 and 1C4 conformers that matches the experimental results best. The minimal deviation between the theoretically- (via. Haasnoot-Altona relation72) and experimentally-derived 3JH,H values can be achieved by assuming the ∆Fi value equal to about +1.1 kJ/mol which corresponds to the 39%:61% (1C4: 4C1) equilibrium mixture of chair conformers. This is in a fair agreement with the predictions of the current set of parameters (∆Fi = -4.1 kJ/mol for O1-methylated and +4.9 kJ/mol for unfunctionalized α-L-IdoU) in spite that the signs of the calculated ring-inversion energies may be opposite. However, the deviations from the experimental data are relatively minor (up to 8 Hz for the whole set of data) if accepting any ∆Fi value from the range (-4.3 kJ/mol; +5 kJ/mol). This speaks for non-negligible populations of both regular and inverted chair conformers and significant conformational heterogeneity of the αL-IdoU

molecule, in contrary to the β-D-GlcU one. ACS Paragon Plus Environment

Page 17 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3.2 Orientation of the exocyclic groups

Fig. 5. The probability distribution plots related to the value of the φ dihedral angle that describes the conformation of the lactol group (defined by the O5-C1-O1-H1 quadruplet of atoms). All the plots correspond to the unbiased MD simulations in the explicit water solvent.

There is no direct experimental data related to the orientation of the lactol groups in uronate saccharides. However, the conformational distributions can be checked for agreement with the trends expected considering the exo-anomeric effect.84 Regarding the D-pyranose molecule in the regular chair conformation, due to this effect, the most populated rotamers are expected to be gfor the β-anomers (i.e. φ ≈ -60 deg) and g+ for the α-anomers (i.e. φ ≈ +60 deg), whereas the t rotamer (i.e. φ ≈ ±180 deg) should be essentially absent. These expectations are reflected very well by the results of the MD simulations shown in Fig. 5 when considering GlcU and GalU monomers in different protonation states. The agreement with the analogous data collected for unfunctionalized aldohexopyranoses (see e.g. Fig. 11 in ref.37) indicates that the influence of the carboxyl moiety on the orientation of the lactol group is very limited; the same is true in the case of the influence of the orientation of the hydroxyl group at C4 (which results in minor differences between results obtained for GlcU and GalU). This is additionally confirmed by the fact that also the protonation state is not relevant for the conformation of the lactol moiety. The deprotonated carboxyl groups of both GlcU and GalU exhibit a very high degree of conformational flexibility (Fig. 6). The anomery of the saccharide is irrelevant for the conformation of the -COO- group, in contrary to the orientation of the substituent at C4 which, when oriented axially, reduces rotational freedom of -COO-. However, in both possible cases, changing the value of the γ angle is associated with overcoming very low free energy barriers: ~5 kJ/mol (for GalU) and 2 kJ/mol (for GlcU). In the case of GalU, the least favorable conformations (γ = ~60 deg and γ = ~-120 deg) are associated with the close proximity of any of the carboxyl oxygen atoms and O4 (~ 0.26 nm) when they are in the syn position. On the contrary the most ACS Paragon Plus Environment


favorable conformations (γ = ~145 deg and γ = ~-35 deg) correspond to perpendicular orientation of the carboxyl group plane with respect to the C5-H covalent bond. In the case of GlcU, the most favorable orientations (γ = ~30 deg and γ = ~-150 deg) are represented by the perpendicular orientation of the carboxyl group plane with respect to the C5-O5 covalent bond whereas the least favorable ones by the perpendicular orientation of the -COO- plane with respect to the C5-C4 covalent bond.

Fig. 6. The probability distribution plots related to the value of the γ dihedral angle that describes the conformation of the carboxyl group (defined by the C4-C5-O6-O61 quadruplet of atoms). All the plots correspond to the unbiased MD simulations in the explicit water solvent.

The protonation of the carboxyl group is associated with a significant reduction of the available conformational space. For each of the compounds only two orientations exhibit the nonnegligible populations, i.e. the γ = ~70 deg and γ = ~130 deg values which correspond to primary rotamers for GlcU and GalU, respectively, whereas γ = ~-140 deg and γ = ~-40 deg which represent the secondary conformers. Again, the anomery is not relevant, in contrary to the orientation of the hydroxyl group at C4. The γ values characterizing the most populated rotamers are very close to those observed for deprotonated moieties of the corresponding compounds. However, due to the loss of internal symmetry of the carboxyl group (i.e. the attachment of proton to O62), the resulting asymmetry in distribution of rotamers is observed. Interestingly, the favorable orientations of the protonated carboxyl group are those that permit the formation of hydrogen bonding with hydroxyl group at C4 (which can be either axially- or equatorially-oriented as reflected by different positions of dominating peaks in Fig. 6). Although this effect is not necessarily imposed by the presence of such interactions,85,86 the behavior opposite to that typical for unfunctionalized hexopyranoses and hydroxymethyl group is worth emphasizing. Although the X-ray diffraction (XRD) structures cannot be treated in terms of a direct evidence for the conformational behavior of saccharide molecules in solutions, they still may ACS Paragon Plus Environment

Page 18 of 37

Page 19 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


provide valuable qualitative information about the conformational preferences of interest which includes the insight into the most favorable geometries or conformationally-restricted regions. Fig. 7 contains the XRD data extracted from the PDB database and represented by the value of the γ dihedral angle. The comparison with the results of MD simulations (Fig. 6) reveals a good agreement for a series of naturally abundant D-uronates, i.e. α-GalU, β-GlcU, β-ManU, as well as for α-L-GulU. For compounds with an equatorially oriented hydroxyl group at C4 (with respect to the ring in the 4C1 conformation, i.e. β-GlcU and β-Man) the available conformational space covers all the possible values of the γ dihedral angle. Upon changing the mutual orientation of the mentioned hydroxyl group and carboxyl group (i.e. for α-L-GulU and α-D-GalU), a much clearer preference for conformations around γ = -60/+120 deg is visible whereas some of the conformations are excluded (i.e. around γ = 0/±180 deg). This remains in agreement with the course of the probability distribution plots shown in Fig. 6. To our knowledge, there is no structural data on acidic uronates containing the explicitly-defined location of carboxylate-bound proton, to perform analogous comparison with the MD results. Note that the use of the default, carboxylate-related GROMOS parameters results in predicting a notably different rotational behavior of the carboxyl moiety; see Fig. 3 in ref.34 for the quantitative results. In particular, the free energy barriers separating particular rotamers are much higher than those predicted in the present study. Moreover, in contrary to the results presented in Fig. 6, the compounds with a possible syn orientation of both O4 and O61/O62 (i.e. GalU and GulU) exhibit the higher rotational flexibility of the carboxyl group in comparison to e.g. GlcU and ManU.34

Fig. 7. The values of the γ dihedral angle, describing the orientation of the carboxyl group in uronates, according to the XRD data deposited in the PDB database.



3.3 Orientation of the glycosidic linkages The conformation of longer, oligo- or polymeric saccharide chains is controlled mainly by the orientation of glycosidic linkages between residues. According to the results described elsewhere,19 in the case of uronate residues, the conformation of glycosidic linkage characteristic of dimers can be slightly different than that of longer chains. Therefore, in the present study we focused on the longer, octameric chains which better reflect the conformational behaviour of natural polyuronate molecules. Fig. 8 shows the calculated free energy maps relying on the glycosidic φ and ψ dihedral angles.

Fig. 8. Free-energy maps in the space of the glycosidic φ and ψ dihedral angles for nine different octamers composed on uronate residues exhibiting different protonation/esterification states, calculated using the presently proposed extension of the GROMOS 56a6CARBO force field in the presence of explicit water. The contours are plotted every 20 kJ/mol.


Page 20 of 37

Page 21 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The results indicate that the chemical state of the carboxyl group is not relevant in the context of the conformation of longer saccharide chains. Both the location of the main minima of the free energy and the dimension of the available conformational phase space remain nearly unchanged, independently if the carboxyl group is protonated, deprotonated or esterified. Moreover, the favourable conformations are very close to those characteristic of the corresponding unfunctionalized hexopyranoses37 which suggests that the type of exocyclic substituent is not crucial for the dynamic equilibrium related to the chain shape. Likeliness of free energy maps calculated for ManU and GlcU shows than the primary factor controlling the orientation of glycosidic linkage is the disposition of glycosidic oxygen atoms with respect to the ring of the given residue (i.e. axial vs. equatorial). The favourable conformation of non-monomeric saccharide connected by the given type of glycosidic linkage always corresponds to the same region of φ vs. ψ phase space. This is in agreement with the results of other studies, performed for unfunctionalized oligohexopyranoses78,87 and disaccharides.86 Furthermore, this effect is not disturbed much upon the inversion of the ring from 4C1 to 1C4 as indicated by extremely similar maps obtained for α-L-GulU (diaxially linked 1C4 rings) and β-D-GlcU (diequatorially linked 4C1 rings). This observation is also common for unfunctionalized aldohexoses.87 Fig. S3 (Supporting Information) contains the glycosidic φ and ψ dihedral angle values extracted from the XRD data measured for pectins and mannuronate-rich alginates. The comparison with the free-energy maps shown in Fig. 8 reveals an excellent agreement with respect to the location of the most favourable conformations for both β-D-ManU and α-D-GalU following the trends expected on the basis of the stereoelectronic effects governing the orientation of glycosidic linkages.88,89

3.3 Chain-chain complexes In order to obtain the correct Ca2+-carboxylate binding mode (i.e. the so-called 'tight' binding) and stable calcium-uronate complexes, the introduction of the Project’s correction90 in the default LJ parameters for the CA2+ and OM atom types is necessary. Otherwise, only the 'loose' binding, contradictory to the experimental observations5,6,10 can be observed.35 The analogous correction is required also in the case of other biomolecular force fields.90 The ability for binding calcium ions and creating the stable complexes composed of two polyuronate chains and calcium ions was validated in the context of stability of complexes created by the two chains of α-L-GulU and calcium ions. The geometry of such complexes is close to the so-called ‘egg-box’ structural model5,6,10,32,91 and is described in details elsewhere.15–17 The simulations performed with use of the current version of parameters indicate that: (i) the Ca2+ACS Paragon Plus Environment


guluronate complex is stable during the 200 ns-long MD run; (ii) all the basic structural parameters (i.e. coordination of the Ca2+ ions, conformation of the glycosidic linkages, helical twist of the whole complex, etc.) are very close to those reported in refs.15,16 and based on the application of the GROMOS 45a4 force field. This is understandable in view of the close similarities between the 45a4 and 56a6CARBO sets related to the treatment of the glycosidic linkages; these similarities are described in ref.37 The overall structure is different from the classical egg-box model91,92 and has been confirmed by the more accurate, QM-MM/MD simulations.18 As the present paper contains one of the first descriptions of the MD-based structure of the pectin-Ca2+ complexes,32 let us summarize briefly the most important structural parameters found in our investigation. At the same time, let us note that a more detailed study devoted to the conformational properties of pectins, including structural heterogeneity effects (i.e. esterification and/or presence of non-α-D-GalU units) is postponed to our future work. In the case of both parallel and antiparallel arrangements of the two Ca2+-linked α-D-GalU chains, we observed the interaction pattern similar to that observed in the Ca2+-guluronate complexes.15 The calcium ions are coordinated by eight oxygen atoms; the four of them belong to the two carboxyl moieties of the two paired chain whereas the remaining four to the water molecules. No other atoms contribute to the direct contact with the calcium ions. The antiparallel arrangement of the α-D-GalU chains results in the close proximity of the C2and C3-linked hydroxyl groups belonging to two different chains. This pattern involves both Ca2+binding residues and the remaining ones. Some of these proximities are associated with interchain hydrogen bonding. On the contrary, the parallel arrangement does not favor any direct contacts between chains and the whole structure is more loose than in the case of the antiparallel arrangement. In addition to the results described in ref.,15 let us summarize some further findings originating from the current simulations. Some of them are graphically shown in Figs. 9 and 10, whereas Fig. 11 presents the simulation snapshots illustrating typical structural properties of the investigated complexes. Regarding α-L-GulU-containing complexes, although both the parallel and antiparallel binding modes appear to be stable over the timescale of hundreds of nanoseconds, there exist several differences between them. Firstly, the conformational behavior considered in the context of single oligomeric chain may change upon complexation by the second chain and calcium ions; the degree of such change is dependent on the type of chain-chain pairing. The antiparallel arrangements of both α-L-GulU and α-D-GalU result in larger changes of the following conformational descriptors: (i) the magnitude of the helical twist of single chain (χ, defined by the C6-C2-C’2-C’6 quadruplets of atoms); (ii) the rotation of the carboxyl moieties ACS Paragon Plus Environment

Page 22 of 37

Page 23 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


engaged in calcium binding. Although the average twist exhibited by the α-L-GulU chains is roughly the same, independently of their molecular environment (χ ≈ -30 deg), the flexibility of the chain is significantly reduced upon anitiparallel binding. On the contrary, parallel arrangement does not significantly influence the helical twist. The same observation is also valid for the systems containing α-D-GalU chains. The only difference are the average twists which differ depending if either single chain (χ ≈ -30 deg) or parallel (χ ≈ 0 deg) and antiparallel (χ ≈ 0 deg) complexes are considered. Interestingly, the changes in the favorable χ value are not associated with the corresponding divergences in the glycosidic dihedral angle values which is in agreement with the results reported in ref.

19

. Note that the complexation of both α-L-GulU and α-D-GalU

chains may result in a small reduction of their flexibilities expressed in terms of the φ and ψ values. Furthermore, the antiparallel binding mode results in slightly modified pattern of rotational preferences of those carboxyl groups which contribute to the Ca2+-mediated chain-chain binding. Namely, the preferred conformation is shifted toward higher values of the γ angle (which results in the favorable eclipsed conformation with respect to the C4 and O6 atoms). Additionally, in the case of α-D-GalU, the rotational freedom of the carboxyl moiety is slightly reduced. The diverse results obtained for galactouronate- and guluronate-containing complexes in the context of their helical twist may contribute to the mechanistic interpretation of different gelling behaviors observed experimentally for pectins and alginates. A more linear form of chain-chain complexes typical for galactouronates may favor the association of further chains (also linked by the contribution of the Ca2+ ion) and stabilize the junction zone. The larger helical twist observed in the case of guluronates results in a relatively smaller area of possible contacts between two (or more) complexes.16 Nevertheless, due to a more heterogeneous structure of naturally-occurring pectins, a more detailed study is required in this matter. The differences between parallel and antiparallel arrangements of chains can be ascribed to the molecular details of the interactions between chains. Independently of the considered compounds, the antiparallelly oriented chains exhibit much closer contact with each other in comparison to the parallel arrangements. This is expressed by the value of the α angle which describes the orientation between the planes of the two chains in the manner depicted symbolically in Figs. 9 and 10 (the regular α angle is defined by the C2-Ca2+-C2 triplet of atoms where the C2 atoms belong to the two residues connected by the same Ca2+ ion). Strictly distinct distributions of this angle value speaks for much tighter steric arrangement observed in antiparallel complexes which is correlated with a more intensive network of chain-chain hydrogen bonds and close distance between those edges of the chains which are not involved in the Ca2+-mediated chainchain binding. These stronger chain-chain interactions result in alterations of the conformational



degrees of freedom, discussed above. On the contrary, the parallel arrangement results in a more loose structure of the whole complex, demonstrated by the values of conformation-related parameters close to those characteristic of single chains. Apart from using the distributions of the α angle values (Figs. 9 and 10), the same effect can be illustrated by calculating the radial distribution functions (RDFs) for the distance between the two residues bound by the same calcium ion. The results of such calculations are graphically presented in Fig. S4 (Supporting Information). Finally, let us note that a more loose structure of parallelly-arranged complexes is reflected by the location of the bound calcium ions. In the case of antiparallel chain-chain orientation, all bound Ca2+ ions are located in one straight row, i.e. the angle between any consecutive three of them is close to 180 deg. One the contrary, in the case of parallel complexes, the conformational rearrangements, occurring in a relatively large timescale (i.e. tens to hundreds of ns) can be seen, which correspond to local deformations of the complex and shift some of the calcium ions out of the ideal straight line. Such behavior was found to be independent of the initial structure and occurred multiple times during the single 200 ns-long simulation; thus, we concluded that it is an inherent feature of the parallel arrangement described in terms of the current force field parameters.


Page 24 of 37

Page 25 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 9. The probability distribution plots related to the values of the dihedral or regular angles describing: (i) the conformation of the glycosidic linkage (φ and ψ); (ii) the helical twist of the single chain either involved or not in formation of the complex (χ); (iii) the conformation of the carboxyl group (γ). Additionally, the regular α angle (defined by the C2-Ca2+-C2 triplet where the C2 atoms belong to the two residues connected by Ca2+ ion) describing the mutual orientation of the two chains that create the complex, is considered. The schematic plot shows the differences between both possible chain-chain arrangements. All the plots correspond to the unbiased MD simulations of the α-L-GulU-Ca2+ complexes in the explicit water solvent.



Fig. 10. The probability distribution plots related to the values of the dihedral or regular angles related to: (i) the glycosidic linkage (φ and ψ); (ii) the helical twist of the uronate chain (χ); (iii) the conformation of the carboxyl group (γ); (iv) the mutual orientation of the two chains creating the complex (α). All the plots correspond to the unbiased MD simulations of the α-D-GalU-Ca2+ complexes in the explicit water solvent. Other details as in Fig. 9.


Page 26 of 37

Page 27 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 11. The compilation of snapshots extracted from the unbiased MD simulations performed within the current set of force field parameters. Both view along the complex axis (right) and the central fragments of complexes (left) are shown. Calcium ions are represented as yellow balls. Water molecules (present during simulations) are not shown for clarity.



Conclusions The article presents an extension of the GROMOS 56a6CARBO/CARBO_R force field aimed at description of the uronate mono-, di-, oligo- and polysaccharides. The derived set of parameters takes into account the different chemical states of carboxylate groups (anionic vs. protonated vs. esterified). In combination with the previously proposed sets, the present extension can be used to study the conformational features of all naturally occurring uronate molecules including those composed of heterogeneous uronate residues (e.g. alginates) or incorporating non-uronate components (e.g. pectins). The parameterization procedure was based both on the quantummechanical calculations (in the context of carboxyl groups) and incorporating the already existing parameters for non-carboxyl parts of uronates. The validation tests proved the agreement with the external data in the context of glycosidic linkage conformation, rotameric freedom of carboxyl group, ring-inversion properties and the stability of the complexes formed by the Ca2+ ions and polyuronate chains. The results obtained for α-L-iduronate suggest that the conformational equilibria of the ring shapes is dominated by the chair conformers and the two-state model indicates that the composition of the equilibrium mixture is close to 39%:61% (1C4:4C1). The protonation state of the carboxyl moiety may influence the flexibility of the pyranose ring; the anionic monomers are usually more prone to the ring inversion in comparison to their protonated analogues and corresponding unfunctionalized aldohexopyranoses. The anionic carboxyl group is very prone to rotation, with the corresponding free energy barriers varying from ~2 to 5 kJ/mol. The main factor controlling its conformational preferences is the disposition of the hydroxyl group at the C4 atom. The protonation of the carboxyl group results in a significant reduction of the rotational freedom. Conversely, the chemical state of the carboxyl group has no influence on the conformation of the glycosidic linkages in uronate oligomers. Regarding the Ca2+-uronate complexes, both α-L-GulU and α-D-GalU are able to create stable structures according to either parallel or antiparallel arrangements of the chains. The pattern of interactions with calcium ions is the same as that reported in ref.,15 i.e. only the calcium-carboxylate contacts contribute to the tight ion binding. The antiparallel arrangements of both α-L-GulU and α-D-GalU chains result in more compact structures with closer contacts between chains; on the contrary, parallel arrangements are associated with more loose structures and less intensive chain-chain contacts. Differences between guluronate- and galacturonate complexes suggest the more facile contacts between the already paired chains in the latter case.


Page 28 of 37

Page 29 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Supporting Information The supporting information consists of plots showing the performance of the parameters adopted from the 53a6 set in comparison with the QM data, convergence plots related to the ring-inversion free energies, RDFs showing different chain-chain contact patterns for both parallel and antiparallel arrangements of chain-chain complexes and the PDB-extracted data showing distribution of glycosidic dihedral angles for various uronates. This material is available free of charge via the Internet at http://pubs.acs.org.

Acknowledgments The authors acknowledge the financial support of the Polish National Science Centre (contract financed in 2016–2020 under Project No. 2015/18/E/ST4/00234 SONATA BIS). The corresponding author acknowledges the partial support from Polish National Science Centre (Project No. 2015/17/B/NZ9/03589 OPUS 9).

References (1)

Ovodov, Y. S. Structural Chemistry of Plant Glycuronoglycans. Pure Appl. Chem. 2009, 42, 351–369.

(2)

Ramesh, H. P. F.; Tharanathan, R. N. Carbohydrates-the Renewable Raw Materials of High Biotechnological Value. Crit. Rev. Biotechnol. 2003, 23, 149–173.

(3)

Smidsrod, O.; Draget, K. I. Chemistry and Physical Properties of Alginates. Carbohydr. Eur. 1996, 14, 6–13.

(4)

Ridley, B. L.; O’Neill, M. A.; Mohnen, D. Pectins: Structure, Biosynthesis, and Oligogalacturonide-Related Signaling. Phytochemistry 2001, 57, 929–967.

(5)

Rees, D. A.; Welsh, E. J. Secondary and Tertiary Structure of Polysaccharides in Solutions and Gels. Angew. Chem. Int. Ed. Engl. 1977, 16, 214–224.

(6)

Morris, E. R. Molecular Interactions in Polysaccharide Gelation. Br. Polym. J. 1986, 18, 14– 21.

(7)

Rees, D. A. Polysaccharide Shapes and Their Interactions - Some Recent Advances. Pure Appl. Chem. 1981, 53, 1–14.

(8)

Willats, W. G.; McCartney, L.; Mackie, W.; Knox, J. P. Pectin: Cell Biology and Prospects for Functional Analysis. Plant Mol. Biol. 2001, 47, 9–27.

(9)

Pérez, S.; Mazeau, K.; Hervé du Penhoat, C. The Three-Dimensional Structures of the Pectic Polysaccharides. Plant Physiol. Biochem. 2000, 38, 37–55.



(10) Rees, D. A. Polysaccharide Conformation in Solutions and Gels - Recent Results on Pectins. Carbohydr. Polym. 1982, 2, 254–263. (11) Willats, W. G.; Orfila, C.; Limberg, G.; Buchholt, H. C.; van Alebeek, G. J.; Voragen, A. G.; Marcus, S. E.; Christensen, T. M.; Mikkelsen, J. D.; Murray, B. S., et al. Modulation of the Degree and Pattern of Methyl-Esterification of Pectic Homogalacturonan in Plant Cell Walls. Implications for Pectin Methyl Esterase Action, Matrix Properties, and Cell Adhesion. J. Biol. Chem. 2001, 276, 19404–19413. (12) Draget, K. I.; Smidsrød, O.; Skjåk-Bræk, G. Alginates from Algae. In Alginates from Algae in: Polysaccharides and Polyamides in the Food Industry. Properties, Production, and Patents.; WILEY-VCH Verlag GmbH & Co. KGaA: Weinheim, 2005. (13) Heyraud, A.; Dantas, L.; Courtois, J.; Courtois, B.; Helbert, W.; Chanzy, H. Crystallographic Data on Bacterial (1 → 4)-β-D-Glucuronan. Carbohydr. Res. 1994, 258, 275–279. (14) Braccini, I.; Heyraud, A.; Pérez, S. Three-Dimensional Features of the Bacterial Polysaccharide (1 → 4)-β-D-Glucuronan: A Molecular Modeling Study. Biopolymers 1998, 45, 165–175. (15) Plazinski, W. Molecular Basis of Calcium Binding by Polyguluronate Chains. Revising the Egg-Box Model. J. Comput. Chem. 2011, 32, 2988–2995. (16) Plazinski, W.; Rudzinski, W. Molecular Modeling of Ca2+-Oligo(α-L-Guluronate) Complexes: Toward the Understanding of the Junction Zone Structure in Calcium Alginate Gels. Struct. Chem. 2012, 23, 1409–1415. (17) Plazinski, W.; Drach, M. The Dynamics of the Calcium-Induced Chain-Chain Association in the Polyuronate Systems. J. Comput. Chem. 2012, 33, 1709–1715. (18) Plazinski, W.; Drach, M. Calcium-α-L-Guluronate Complexes: Ca2+ Binding Modes from DFT-MD Simulations. J. Phys. Chem. B 2013, 117, 12105–12112. (19) Plazinski, W. Conformational Properties of Acidic Oligo- and Disaccharides and Their Ability to Bind Calcium: A Molecular Modeling Study. Carbohydr. Res. 2012, 357, 111– 117. (20) Stewart, M. B.; Gray, S. R.; Vasiljevic, T.; Orbell, J. D. Exploring the Molecular Basis for the Metal-Mediated Assembly of Alginate Gels. Carbohydr. Polym. 2014, 102, 246–253. (21) Wolnik, A.; Albertin, L.; Charlier, L.; Mazeau, K. Probing the Helical Forms of Ca2+Guluronan Junction Zones in Alginate Gels by Molecular Dynamics 1: Duplexes. Biopolymers 2013, 99, 562–571. (22) Braccini, I.; Pérez, S. Molecular Basis of Ca2+-Induced Gelation in Alginates and Pectins: The Egg-Box Model Revisited. Biomacromolecules 2001, 2, 1089–1096. (23) The Amber Molecular Dynamics Package. http://ambermd.org/ (accessed Feb 6, 2018). ACS Paragon Plus Environment

Page 30 of 37

Page 31 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(24) Kirschner, K. N.; Yongye, A. B.; Tschampel, S. M.; González-Outeiriño, J.; Daniels, C. R.; Foley, B. L.; Woods, R. J. GLYCAM06: A Generalizable Biomolecular Force Field. Carbohydrates. J. Comput. Chem. 2008, 29, 622–655. (25) Dentini, M.; Rinaldi, G.; Barbetta, A.; Risica, D.; Anselmi, C.; Skjåk-Bræk, G. Ionic Gel Formation of a (Pseudo)alginate Characterised by an Alternating MG Sequence Produced by Epimerising Mannuronan with AlgE4. Carbohydr. Polym. 2007, 67, 465–473. (26) Arrate, M.; Aurrecoechea, J. M. Synthesis of Bicyclic Alcohols by Palladium-Catalyzed Et2Zn-Mediated Intramolecular Carbonylpropargylation. Arkivoc 2017, 257-267. (27) Guvench, O.; Mallajosyula, S. S.; Raman, E. P.; Hatcher, E.; Vanommeslaeghe, K.; Foster, T. J.; Jamison II, F. W.; MacKerell Jr., A. D. CHARMM Additive All-Atom Force Field for Carbohydrate Derivatives and Its Utility in Polysaccharide and Carbohydrate–Protein Modeling. J. Chem. Theory Comput. 2011, 7, 3162-3180. (28) Dauber-Osguthorpe, P.; Roberts, V. A.; Osguthorpe, D. J.; Wolff, J.; Genest, M.; Hagler, A. T. Structure and Energetics of Ligand Binding to Proteins: Escherichia Coli Dihydrofolate Reductase-Trimethoprim, a Drug-Receptor System. Proteins 1988, 4, 31–47. (29) Kuttel, M.; Brady, J. W.; Naidoo, K. J. Carbohydrate Solution Simulations: Producing a Force Field with Experimentally Consistent Primary Alcohol Rotational Frequencies and Populations. J. Comput. Chem. 2002, 23, 1236–1243. (30) Liu, Y.; Chipot, C.; Shao, X.; Cai, W. Solubilizing Carbon Nanotubes through Noncovalent Functionalization. Insight from the Reversible Wrapping of Alginic Acid around a SingleWalled Carbon Nanotube. J. Phys. Chem. B 2010, 114, 5783–5789. (31) Babin, V.; Sagui, C. Conformational Free Energies of Methyl-α-L-Iduronic and Methyl-β-DGlucuronic Acids in Water. J. Chem. Phys. 2010, 132, 104108. (32) Assifaoui, A.; Lerbret, A.; Uyen, H. T. D.; Neiers, F.; Chambin, O.; Loupiac, C.; Cousin, F. Structural Behaviour Differences in Low Methoxy Pectin Solutions in the Presence of Divalent Cations (Ca(2+) and Zn(2+)): A Process Driven by the Binding Mechanism of the Cation with the Galacturonate Unit. Soft Matter 2015, 11, 551–560. (33) Labus, K.; Radosiński, Ł.; Kuchta, B.; Trusek-Hołownia, A. Temperature Effect on Properties of β-Galactosidase Entrapped in Alginate Matrix: An Experimental Research Supported by Molecular Modeling. Inż. Apar. Chem. 2015, 3, 98-100. (34) Perić, L.; Pereira, C. S.; Pérez, S.; Hünenberger, P. H. Conformation, Dynamics and IonBinding Properties of Single-Chain Polyuronates: A Molecular Dynamics Study. Mol. Simul. 2008, 34, 421–446.



(35) Perić-Hassler, L.; Hünenberger, P. H. Interaction of Alginate Single-Chain Polyguluronate Segments with Mono- and Divalent Metal Cations: A Comparative Molecular Dynamics Study. Mol. Simul. 2010, 36, 778–795. (36) Anderson, R. L.; Greenwell, H. C.; Suter, J. L.; Coveney, P. V.; Thyveetil, M.-A. Determining Materials Properties of Natural Composites Using Molecular Simulation. J. Mater. Chem. 2009, 19, 7251–7262. (37) Hansen, H. S.; Hünenberger, P. H. A Reoptimized GROMOS Force Field for HexopyranoseBased Carbohydrates Accounting for the Relative Free Energies of Ring Conformers, Anomers, Epimers, Hydroxymethyl Rotamers, and Glycosidic Linkage Conformers. J. Comput. Chem. 2011, 32, 998–1032. (38) Plazinski, W.; Lonardi, A.; Hünenberger, P. H. Revision of the GROMOS 56a6CARBO Force Field: Improving the Description of Ring-Conformational Equilibria in Hexopyranose-Based Carbohydrates Chains. J. Comput. Chem. 2016, 37, 354–365. (39) Oostenbrink, C.; Villa, A.; Mark, A. E.; van Gunsteren, W. F. A Biomolecular Force Field Based on the Free Enthalpy of Hydration and Solvation: The GROMOS Force-Field Parameter Sets 53A5 and 53A6. J. Comput. Chem. 2004, 25, 1656–1676. (40) Oostenbrink, C.; Soares, T. A.; van der Vegt, N. F. A.; van Gunsteren, W. F. Validation of the 53A6 GROMOS Force Field. Eur. Biophys. J. 2005, 34, 273–284. (41) Lins, R. D.; Hünenberger, P. H. A New GROMOS Force Field for Hexopyranose-Based Carbohydrates. J. Comput. Chem. 2005, 26, 1400–1412. (42) Reif, M. M.; Winger, M.; Oostenbrink, C. Testing of the GROMOS Force-Field Parameter Set 54A8: Structural Properties of Electrolyte Solutions, Lipid Bilayers, and Proteins. J. Chem. Theory Comput. 2013, 9, 1247–1264. (43) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785–789. (44) Vosko, S. H.; Wilk, L.; Nusair, M. Accurate Spin-Dependent Electron Liquid Correlation Energies for Local Spin Density Calculations: A Critical Analysis. Can. J. Phys. 1980, 58, 1200–1211. (45) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98, 11623–11627. (46) Becke, A. D. Density-functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648–5652. (47) Hehre, W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. AB INITIO Molecular Orbital Theory; John Willey & Sons, Inc.: New York, 1986. ACS Paragon Plus Environment

Page 32 of 37

Page 33 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(48) Bayly, C. I.; Cieplak, P.; Cornell, W.; Kollman, P. A. A Well-Behaved Electrostatic Potential Based Method Using Charge Restraints for Deriving Atomic Charges: The RESP Model. J. Phys. Chem. 1993, 97, 10269–10280. (49) Horta, B. A. C.; Fuchs, P. F. J.; van Gunsteren, W. F.; Hünenberger, P. H. New Interaction Parameters for Oxygen Compounds in the GROMOS Force Field: Improved Pure-Liquid and Solvation Properties for Alcohols, Ethers, Aldehydes, Ketones, Carboxylic Acids, and Esters. J. Chem. Theory Comput. 2011, 7, 1016–1031. (50) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Petersson, G. A.; Nakatsuji, H., Li, X., et al. Gaussian 09, Revision A.02; Gaussian, Inc.: Wallingford, CT, 2016. (51) Hanwell, M. D.; Curtis, D. E.; Lonie, D. C.; Vandermeersch, T.; Zurek, E.; Hutchison, G. R. Avogadro: An Advanced Semantic Chemical Editor, Visualization, and Analysis Platform. J. Cheminformatics 2012, 4, 17. (52) Rappe, A. K.; Casewit, C. J.; Colwell, K. S.; Goddard, W. A.; Skiff, W. M. UFF, a Full Periodic Table Force Field for Molecular Mechanics and Molecular Dynamics Simulations. J. Am. Chem. Soc. 1992, 114, 10024–10035. (53) Gaweda, K.; Plazinski, W. Pyranose Ring Conformations in Mono- and Oligosaccharides: A Combined MD and DFT Approach. Phys. Chem. Chem. Phys. 2017, 19, 20760–20772. (54) Pronk, S.; Páll, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; Shirts, M. R.; Smith, J. C.; Kasson, P. M.; van der Spoel, D., et al. GROMACS 4.5: A High-Throughput and Highly Parallel Open Source Molecular Simulation Toolkit. Bioinformatics 2013, 29, 845– 854. (55) Bonomi, M.; Branduardi, D.; Bussi, G.; Camilloni, C.; Provasi, D.; Raiteri, P.; Donadio, D.; Marinelli, F.; Pietrucci, F.; Broglia, R. A., et al. PLUMED: A Portable Plugin for FreeEnergy Calculations with Molecular Dynamics. Comput. Phys. Commun. 2009, 180, 1961– 1972. (56) Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.; Hess, B.; Lindahl, E. GROMACS: High Performance Molecular Simulations through Multi-Level Parallelism from Laptops to Supercomputers. SoftwareX 2015, 1–2, 19–25. (57) Tribello, G. A.; Bonomi, M.; Branduardi, D.; Camilloni, C.; Bussi, G. PLUMED 2: New Feathers for an Old Bird. Comput. Phys. Commun. 2014, 185, 604–613. (58) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. Interaction Models for Water in Relation to Protein Hydration. In Intermolecular Forces; Pullman, E., Ed.; Reidel: Dordrecht, The Netherlands, 1981, pp 331–342.



(59) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling through Velocity Rescaling. J. Chem. Phys. 2007, 126, 014101. (60) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182–7190. (61) Hockney, R. W. Potential Calculation and Some Applications. Methods Comput Phys 1970, 9, 136-211. (62) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463–1472. (63) Hess, B. P-LINCS: A Parallel Linear Constraint Solver for Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 116–122. (64) Berendsen, H. J. C.; Van Gunsteren, W. F.; Zwinderman, H. R. J.; Geurtsen, R. G. Simulations of Proteins in Watera. Ann. N. Y. Acad. Sci. 1986, 482, 269–286. (65) Barker, J. A.; Watts, R. O. Monte Carlo Studies of the Dielectric Properties of Water-like Models. Mol. Phys. 1973, 26, 789–792. (66) Tironi, I. G.; Sperb, R.; Smith, P. E.; van Gunsteren, W. F. A Generalized Reaction Field Method for Molecular Dynamics Simulations. J. Chem. Phys. 1995, 102, 5451–5459. (67) Heinz, T. N.; van Gunsteren, W. F.; Hünenberger, P. H. Comparison of Four Methods to Compute the Dielectric Permittivity of Liquids from Molecular Dynamics Simulations. J. Chem. Phys. 2001, 115, 1125–1136. (68) Cremer, D.; Pople, J. A. General Definition of Ring Puckering Coordinates. J. Am. Chem. Soc. 1975, 97, 1354–1358. (69) Sugita, Y.; Okamoto, Y. Replica-Exchange Molecular Dynamics Method for Protein Folding. Chem. Phys. Lett. 1999, 314, 141–151. (70) Bussi, G.; Gervasio, F. L.; Laio, A.; Parrinello, M. Free-Energy Landscape for β Hairpin Folding from Combined Parallel Tempering and Metadynamics. J. Am. Chem. Soc. 2006, 128, 13435–13441. (71) Barducci, A.; Bussi, G.; Parrinello, M. Well-Tempered Metadynamics: A Smoothly Converging and Tunable Free-Energy Method. Phys. Rev. Lett. 2008, 100, 020603. (72) Haasnoot, C. A. G.; de Leeuw, F. A. A. M.; Altona, C. The Relationship between ProtonProton NMR Coupling Constants and Substituent electronegativities—I. Tetrahedron 1980, 36, 2783–2792. (73) Altona, C.; Francke, R.; de Haan, R.; Ippel, J. H.; Daalmans, G. J.; Hoekzema, A. J. A. W.; van Wijk, J. Empirical Group Electronegativities for Vicinal NMR Proton–proton Couplings along a C-C Bond: Solvent Effects and Reparameterization of the Haasnoot Equation. Magn. Reson. Chem. 1994, 32, 670–678. ACS Paragon Plus Environment

Page 34 of 37

Page 35 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(74) Huggins, M. L. Bond Energies and Polarities1. J. Am. Chem. Soc. 1953, 75, 4123–4126. (75) Sattelle, B. M.; Hansen, S. U.; Gardiner, J.; Almond, A. Free Energy Landscapes of Iduronic Acid and Related Monosaccharides. J. Am. Chem. Soc. 2010, 132, 13132–13134. (76) Sattelle, B. M.; Bose-Basu, B.; Tessier, M.; Woods, R. J.; Serianni, A. S.; Almond, A. Dependence

of

Pyranose

Ring

Puckering

on

Anomeric

Configuration:

Methyl

Idopyranosides. J. Phys. Chem. B 2012, 116, 6380–6386. (77) Mayes, H. B.; Broadbelt, L. J.; Beckham, G. T. How Sugars Pucker: Electronic Structure Calculations Map the Kinetic Landscape of Five Biologically Paramount Monosaccharides and Their Implications for Enzymatic Catalysis. J. Am. Chem. Soc. 2014, 136, 1008–1022. (78) Plazinski, W.; Drach, M.; Plazinska, A. Ring Inversion Properties of 1→2, 1→3 and 1→6Linked Hexopyranoses and Their Correlation with the Conformation of Glycosidic Linkages. Carbohydr. Res. 2016, 423, 43–48. (79) Plazinski, W.; Plazinska, A. Molecular Dynamics Simulations of Hexopyranose Ring Distortion in Different Force Fields. Pure Appl. Chem. 2017, 89, 1283–1294. (80) Rao, V. S. R. Conformation of Carbohydrates; Harwood Academic Publishers: Amsterdam, The Netherlands, 1998. (81) Angyal, S. J. Conformational Analysis in Carbohydrate Chemistry. Aust J Chem 1968, 2737– 2746. (82) Angyal, S. J. The Composition and Conformation of Sugars in Solution. Angew. Chem. Int. Ed. Engl. 1969, 8, 157–166. (83) Panczyk, K.; Plazinski, W. Pyranose Ring Puckering in Aldopentoses, Ketohexoses and Deoxyaldohexoses. A Molecular Dynamics Study. Carbohydr. Res. 2018, 455, 62–70. (84) Lemieux, R. U.; Koto, S.; Voisin, D. The Exo-Anomeric Effect. In Anomeric Effect; Szarek, W. A., Horton, D., Eds.; ACS Symposium Series: Washington, 1979; Vol. 87, pp 17–29. (85) Lonardi, A.; Oborský, P.; Hünenberger, P. H. Solvent‐Modulated Influence of Intramolecular Hydrogen‐Bonding on the Conformational Properties of the Hydroxymethyl Group in Glucose and Galactose: A Molecular Dynamics Simulation Study. Helv. Chim. Acta 2017, 100, e1600158. (86) Wang, D.; Ámundadóttir, M. L.; Gunsteren, W. F. van; Hünenberger, P. H. Intramolecular Hydrogen-Bonding in Aqueous Carbohydrates as a Cause or Consequence of Conformational Preferences: A Molecular Dynamics Study of Cellobiose Stereoisomers. Eur. Biophys. J. 2013, 42, 521–537. (87) Plazinski, W.; Drach, M. The Influence of the Hexopyranose Ring Geometry on the Conformation of Glycosidic Linkages Investigated Using Molecular Dynamics Simulations. Carbohydr. Res. 2015, 415, 17–27. ACS Paragon Plus Environment


(88) Tvaroška, I.; Bleha, T. Anomeric and Exo-Anomeric Effects in Carbohydrate Chemistry. Adv. Carbohydr. Chem. Biochem. 1989, 47, 45–123. (89) Miljkovic, M. Electrostatic and Stereoelectronic Effects in Carbohydrate Chemistry; Springer US: Boston, MA, 2014. (90) Project, E.; Nachliel, E.; Gutman, M. Parameterization of Ca+2-Protein Interactions for Molecular Dynamics Simulations. J. Comput. Chem. 2008, 29, 1163–1169. (91) Morris, E. R.; Rees, D. A.; Thom, D.; Boyd, J. Chiroptical and Stoichiometric Evidence of a Specific, Primary Dimerisation Process in Alginate Gelation. Carbohydr. Res. 1978, 66, 145– 154. (92) Grant, G. T.; Morris, E. R.; Rees, D. A.; Smith, P. J. C.; Thom, D. Biological Interactions between Polysaccharides and Divalent Cations: The Egg-Box Model. FEBS Lett. 1973, 32, 195–198.


Page 36 of 37

Page 37 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC Graphic


CARBO_R Force Field for

Recommend Documents