All-Atom Internal Coordinate Mechanics (ICM) Force Field for Hexopyranoses and Glycoproteins

Abstract

We present an extension of the all-atom internal-coordinate force field, ICMFF, that allows for simulation of heterogeneous systems including hexopyranose saccharides and glycan chains in addition to proteins. A library of standard glycan geometries containing α- and β-anomers of the most common hexapyranoses, i.e., d-galactose, d-glucose, d-mannose, d-xylose, l-fucose, N-acetylglucosamine, N-acetylgalactosamine, sialic, and glucuronic acids, is created based on the analysis of the saccharide structures reported in the Cambridge Structural Database. The new force field parameters include molecular electrostatic potential-derived partial atomic charges and the torsional parameters derived from quantum mechanical data for a collection of minimal molecular fragments and related molecules. The ϕ/ψ torsional parameters for different types of glycosidic linkages are developed using model compounds containing the key atoms in the full carbohydrates, i.e., glycosidic-linked tetrahydropyran-cyclohexane dimers. Target data for parameter optimization include two-dimensional energy surfaces corresponding to the ϕ/ψ glycosidic dihedral angles in the disaccharide analogues, as determined by quantum mechanical MP2/6-31G** single-point energies on HF/6-31G** optimized structures. To achieve better agreement with the observed geometries of glycosidic linkages, the bond angles at the O-linkage atoms are added to the internal variable set and the corresponding bond bending energy term is parametrized using quantum mechanical data. The resulting force field is validated on glycan chains of 1-12 residues from a set of high-resolution X-ray glycoprotein structures based on heavy atom root-mean-square deviations of the lowest-energy glycan conformations generated by the biased probability Monte Carlo (BPMC) molecular mechanics simulations from the native structures. The appropriate BPMC distributions for monosaccharide-monosaccharide and protein-glycan linkages are derived from the extensive analysis of conformational properties of glycoprotein structures reported in the Protein Data Bank. Use of the BPMC search leads to significant improvements in sampling efficiency for glycan simulations. Moreover, good agreement with the X-ray glycoprotein structures is achieved for all glycan chain lengths. Thus, average/median RMSDs are 0.81/0.68 Å for one-residue glycans and 1.32/1.47 Å for three-residue glycans. RMSD from the native structure for the lowest-energy conformation of the 12-residue glycan chain (PDB ID 3og2) is 1.53 Å. Additionally, results obtained for free short oligosaccharides using the new force field are in line with the available experimental data, i.e., the most populated conformations in solution are predicted to be the lowest energy ones. The newly developed parameters allow for the accurate modeling of linear and branched hexopyranose glycosides in heterogeneous systems.

Figure 1

Distribution of numbers of observed monosaccharide moieties at the N-glycosylation sites in PDB X-ray structures. Data was generated by querying PDB using the ICM program.

Figure 2

(a) Distribution of the amide angle of N-acetyl glucosamine in PDB structures. Peaks at ±180° correspond to the expected flat configuration. (b) Example of poor geometry: N-acetylglucosamine at N58 in PDB 1OW0; internal amide bond has a torsion angle at 85°, and the amino acid/sugar amide bond is twisted out of plane by 61°; electron density contouring (brown mesh) reveals that experimental data does not seem to warrant assignment of such a strained conformer.

Figure 3

Atom notation and torsional definitions for the hexopyranose fragment.

Figure 4

Model molecules used for deriving ϕ/ψ torsional parameters for O-linkage: (1) equatorial–equatorial (eq–eq), (2) axial–equatorial (ax–eq), (3) axial–axial (ax–ax), and (4) equatorial–axial (eq–ax) linkages.

Figure 5

Distribution of C₁–O–C valence angles (in degrees) in CSD saccharide structures with R < 10%.

Figure 6

Conformational preferences of different C–O–C disaccharide linkages: (a–d) eq–eq, ax–eq, ax–ax, and eq–ax C–O–C disaccharide linkages, respectively. (Left) Distribution of ϕ/ψ angles in high-resolution PDB structures of glycoproteins with each type of linkage (6319 linkages from 1810 structures for eq–eq, 2303 linkages from 920 structures for ax–eq, 434 linkages from 200 structures for ax–ax, and 239 linkages from 101 structures for eq–ax). (Middle, right) QM and total ICMFF energy surfaces, respectively, for model molecules 1–4 (Figure 4). The color code from purple to red of the energy maps corresponds to the 0–8 kcal/mol range. Contours are drawn with 1 kcal/mol step.

Figure 7

Conformations of the eq–eq model molecule corresponding to the three minima on the ϕ/ψ energy map (Figure 6a).

Figure 8

Distribution of (a) ϕ, (b) ψ, and (c) ω torsional angles in PDB structures of glycoproteins with 1–6 disaccharide linkages. ψ histogram is offset into the 0 to 360° range rather than −180 to 180° to better show the major peak at 180°.

Figure 9

Distribution of ϕ/ψ torsional angles in PDB structures of proteins with N-linkages: (a) α-d-GlcNAc-ASN and (b–d) β-d-GlcNAc-ASN.

Figure 10

Distribution of ϕ/ψ torsional angles in PDB structures of glycoproteins with O-linkages: (a) α-*-Ser and (b) α-*-Thr structures.

Figure 11

Model molecules used for parametrization of ω torsional potential in (a) glucopyranoside and (b) galactopyranoside fragments.

Figure 12

Overlay of the experimental (magenta), the lowest-energy BPMC (green), and the lowest-energy no BPMC (yellow) conformations of the nine-residue glycan chain of PDB 1gai.

Figure 13

Progression of the lowest energy achieved with the time of simulation for PDB 1gai using (a) BPMC steps and ϕ, ψ, ω, and side chain torsional angles as search variables, (b) evenly distributed MC steps and ϕ, ψ, ω, and side chain torsional angles as search variables, and (c) BPMC steps with only ϕ, ψ, and ω angles as search variables.

Figure 14

CH-pi stacking in PDB 3c1u. The experimental BNag conformation is shown in cyan.

Figure 15

Overlay of the experimental (magenta) and the lowest-energy (green) conformations of the 12-residue glycan chain in PDB 3og2. Heavy-atom RMSD between the two structures is 1.76 Å.