Replica exchange molecular dynamics simulations reveal the structural and molecular properties of levan-type fructo-oligosaccharides of various chain lengths

Background Levan and levan-type fructo-oligosaccharides (LFOs) have various potential applications in pharmaceutical and food industries due to their beneficial properties such as their low intrinsic viscosity and high water solubility. Previous studies showed that they exhibited prebiotic effects, anti-inflammatory and anti-tumor activities against Sarcoma-180 tumor cells of human. Despite their various potential applications, the structural and molecular properties of LFOs of various chain lengths are not well understood. Results We employed the replica-exchange molecular dynamics simulations method (REMD) in AMBER14 to elucidate structural and molecular properties of LFOs with chain lengths of 5 (LFO5), 10 (LFO10) and 15 (LFO15) residues in two models of generalized Born implicit solvent (GBHCT and GBOBC1). For LFO10 and LFO15, four distinct conformations (helix-like, partial helix, zig-zag and random structures) were characterized by their upper-middle and lower-middle torsions. For LFO5, two distinct conformations (partial helix and random structures) were characterized by their middle torsion and molecular angle of residues 1, 3 and 5. To determine hydrogen bonds important for the formation of helix-like structures of LFO10 and LFO15, occurrence frequencies of hydrogen bonds were analyzed, and the O6(i)--H3O(i+1) hydrogen bond was found with the highest frequency, suggesting its importance in helix formation. Among three dihedral angles between two fructosyl units [ϕ (O5’-C2’-O6-C6), ψ (C2’-O6-C6-C5) and ω (O6-C6-C5-C4)], dihedral angle distributions showed that ω was the most flexible dihedral angle and probably responsible for conformational differences of LFOs. Conclusions Our study provides important insights into the structural and molecular properties of LFOs, which tend to form helical structures as the chain length increases from 5 to 15 residues. This information could be beneficial for the selection of LFOs with appropriate lengths and properties for pharmaceutical and biological applications. Electronic supplementary material The online version of this article (doi:10.1186/s12859-016-1182-7) contains supplementary material, which is available to authorized users.


Background
Levan-type fructo-oligosaccharides (LFOs) are short chain fructans that contain D-fructofuranosyl residues and are predominantly linked by β- (2,6) linkages in a main chain with some β-(2, 1) linked branching points (Fig. 1a). Produced by levansucrase, levan and LFOs are found in various microorganisms such as Bacillus subtilis [1], Zymomonas mobilis [2], and Leuconostoc mesenteroides [3], and play important roles as sources for energy utilization and biofilm formation [4,5]. The properties of levan and LFOs depend on their chain lengths and branching degrees [6], and they have various desirable properties such as their unusually low intrinsic viscosity [7] and high water solubility [8]. These properties are very beneficial for various industrial applications, especially in food and pharmaceutical industries. For example, levan and LFOs showed prebiotic effects, stimulating the growth of beneficial intestinal bacteria, and also could potentially act as cholesterol lowering agents (MW 2000 kDa) [9,10]. Furthermore, they could also be served as carbon sources for probiotics such as four strains of Bifidobacterium sp. that produce short chain fatty acids, lactate and acetate (MW < 3600 Da) [11]. Moreover, they exhibited antiinflammatory and anti-tumor activities against Sarcoma-180 tumor cells of human (MW 380-710 kDa) [12,13]. Despite their various potential applications, the knowledge on the structural and molecular properties of levan and LFOs of various chain lengths is still limited.
Replica exchange molecular dynamics (REMD) method provides an extensive conformational sampling at various temperatures by allowing each replica to exchange their configurations through parallel tempering [14]. Raising the temperature can increase the probability of the system to overcome the energy barrier, consequently enhancing the probability of attaining the global minimum and allowing the sampling of large volumes of phase space. Therefore, the incorporation of higher temperature systems allow the lower temperature systems to access a representative set of the low free energy minima that are accessible by the higher temperature systems [15]. This technique has been used to investigate the properties of oligosaccharides in solution. For example, Re et al. employed REMD to elucidate the structural diversity and the changes in conformational equilibria of biantennary The parameters used in the characterization of the conformations of LFOs: the upper-middle (χ 6-7 ) and lower-middle (χ 9-10 ) torsions for LFO 15 , (c) the upper-middle (χ 4-5 ) and lower-middle (χ 7-8 ) torsions for LFO 10 , (d) the molecular angles (θ a ) and middle torsions (χ 3-4 ) for LFO 5 . (e) Examples of three dihedral angles between four fructosyl residues of LFO 5 , ω (C4-C5-C6-O6), ψ (C5-C6-O6-C2') and ɸ (C6-O6-C2'-O5'). The carbon, oxygen and hydrogen atoms were colored in green, red and white, respectively. The six-membered ring of the glucosyl residue and the five-membered rings of fructosyl residue are shown in green and yellow, respectively complex-type N-glycans [16]. Moreover, Nishima et al. used this technique to investigate the effects of bisecting GlcNAc and core fucosylation on conformational properties of biantennary complex-type N-glycans [17]. Recently, Jo et al. employed this technique to examine the conformational freedom of the N-glycan core pentasaccharide moiety in solution and found that the conformational variability of the pentasaccharide in solution was more restricted than the N-glycan on the protein surface [18]. This method was also employed to investigate conformational flexibility of cellulose oligomers as well as their chain length and temperature dependence [19]. However, to our knowledge, REMD method has not been used in the elucidation of the properties of LFOs of various chain lengths.
In this study, we performed REMD on the models of LFOs with the chain lengths of 5, 10 and 15 residues in two models of generalized Born implicit solvent (GB HCT and GB OBC1 ) to elucidate their structural and molecular properties as well as the relationship between these properties and the chain length. Such information would be beneficial for the selection of LFOs with appropriate lengths and properties for pharmaceutical and biological applications.

Structure preparation and minimization
The structures of LFO 5 , LFO 10 and LFO 15 were constructed using the LEaP module in AMBER14 [20], and their atom types and force field parameters were assigned based on GLYCAM06j-1 [21]. Two implicit solvent models (GB HCT and GB OBC1 ) were used in the minimization and simulations of each system. All systems were minimized with 2500 steepest-descent minimization cycles and 2500 conjugate-gradient minimization cycles [22,23].

Replica exchange molecular dynamics simulations
Sixteen replicas of each system were initially equilibrated for 500 ps to reach the desired temperature range from 262 to 802 K. REMD of all systems were performed using the SANDER module in AMBER14. Langevin dynamics with a collision frequency of 1 ps −1 were used to control the temperatures in all systems. Initial velocity of each system was reseeded by the random number generator [24]. A cut off of 999 Å was used to truncate nonbonded pairs, and the maximum distance of 999 Å between atom pairs was employed to compute the pairwise summation involved in the effective Born radii calculation. All bond-stretching freedoms associated with hydrogen were eliminated by SHAKE algorithm, allowing a time step of 0.002 ps [25]. Each replica was simulated for 100 ns and exchanged every 2 ps. The replicas at 298 K were employed for the analyses of the structural and molecular properties of LFOs with different chain lengths.
To measure the sizes of all systems, their average radii of gyration (ROG) were determined. To determine possible representative structures of LFOs of each chain length, K-means clustering algorithm, as implemented in MMTSB tool sets [26], was employed to cluster the structures from their 100 ns trajectories based on their structural similarities, calculated from their heavy-atom root-mean-square-deviation. To determine a reasonable representative of each cluster, a structure that is most similar to the average structure of all members of each cluster was chosen to be a "centroid;" i.e., a "centroid" is a structure with the lowest heavy-atom root-meansquare-deviation to the average structure. Based on their shapes, these "centroids" were further classified into helix-like, partial helix, zig-zag or random structures as major representative conformers. Helix-like structures were defined as conformations that had more than 1 helical turn, while partial helix structures were defined as conformations that had 1 helical turn. Zig-zag structures were defined as conformations that had zigzag shapes. Random structures were defined as structures that were not classified as helix-like, partial helix or zig-zag structures.

Reliability of REMD simulations
To determine whether the temperatures were optimally distributed and the number of replicas was sufficient, the acceptance ratios of replica exchange were calculated. The acceptance ratios of the simulations of LFO 15 in the GB HCT model were almost constant around 28 %, implying a free random walk in the replica (temperature) space (Additional file 1: Figure S1a). Moreover, a free random walk both in the replica space (Additional file 1: Figure  S1b) and the temperature space (Additional file 1: Figure  S1c) were also confirmed. Furthermore, the canonical probability distribution of the total potential energy at each temperature had sufficient overlap with those of neighbors (Additional file 1: Figure S1d). The results of the REMD simulations of LFO 10 and LFO 5 in the GB HCT model were also similar, and their average acceptance ratios were almost constant around 37 and 50 % for LFO 10 and LFO 5 , respectively. For the systems simulated in the GB OBC1 model, the results of REMD simulations were also similar to those simulated in the GB HCT model, and their average acceptance ratios were almost constant around 28, 36 and 50 % for LFO 15 , LFO 10 and LFO 5 , respectively. These results indicate good reliability of the REMD simulations of all systems.

Sizes of LFOs
The sizes of LFOs were determined by measuring their radii of gyration. Figure 2 shows that the trends of the radii of gyration of LFOs simulated in the GB HCT model and those simulated in the GB OBC1 model are similar. The radii of gyration of LFOs tended to increase as their chain lengths increased from 5 to 15 residues. These results suggest the extension of the structures of LFOs as their chain lengths increase. 15 , LFO 10 and LFO 5 Figure 3 shows the free-energy maps of LFO 15 , LFO 10 and LFO 5 as simulated in GB HTC and GB OBC1 models as well as their major representative conformers and their population sizes from clustering analysis and centroid classification. For LFO 15 , four major conformations such as helix-like (a), partial helix (b), zig-zag (c) and random (d) structures were observed after clustering analysis and centroid classification (Fig. 3a, b and Additional file 1: Figure S2 and Figure S5), and they were characterized by their upper-middle and lower-middle torsions (χ 6-7 and χ 9-10 ). Helix-like structures were found with the highest population of 54.1 and 63.2 % for those simulated in GB HCT and GB OBC1 models, respectively. Helix-like structures took up conformations of left-handed 3-fold helices and tended to have their upper-middle and lower-middle torsions in the similar range of around 240-315°. The conformations with the second highest population were partial helix structures, and their population sizes were 33.9 and 22.3 % for systems simulated in GB HCT and GB OBC1 models, respectively. The other two conformations were zig-zag and random structures. Zig-zag structures were found with the population sizes of 2.8 and 6.7 % for systems simulated in GB HCT and GB OBC1 models, respectively. The population sizes of random structures simulated in GB HCT and GB OBC1 models were 9.2 and 7.8 %, respectively.

Conformations of LFO
Similar to the conformations of LFO 15 , four major conformations such as helix-like (a), partial helix (b), zig-zag (c) and random (d) structures were found for LFO 10 after clustering analysis and centroid classification (Fig. 3c, d and Additional file 1: Figure S3 and Figure S6). These conformations were characterized by their upper-middle and lower-middle torsions (χ 4-5 and χ 7-8 ). The conformation with the highest population sizes of 50.5 and 57.5 % was helix-like structures for those simulated in GB HCT and GB OBC1 models, respectively. Partial helix structures occurred with the second highest population sizes of 34.5 and 25.4 % for those simulated in GB HCT and GB OBC1 models, respectively. The population sizes of zig-zag structures were 6.8 and 8.7 % and those of random structures were 8.2 and 8.4 % for systems simulated in GB HCT and GB OBC1 models, respectively.
For LFO 5 , two major conformations such as partial helix (b) and random (d) structures were observed after clustering analysis and centroid classification, probably due to its shorter chain length as compared to those of Fig. 2 The average radii of gyration calculated from the heavy atoms of LFOs simulated in GB HCT (square) and GB OBC1 (dot) models LFO 10 and LFO 15 (Fig. 3e, f and Additional file 1: Figure  S4 and Figure S7). These conformations were characterized by their molecular angles (θ a ) and middle torsion (χ 3-4 ). Partial helix structures were observed with the population sizes of 92.8 and 92.5 % for those simulated in GB HCT and GB OBC1 models, respectively. Random structures were also found with the population sizes of 7.2 and 7.5 % for those simulated in GB HCT and GB OBC1 models, respectively. Table 1 shows the populations of major representative conformers of LFO 15 , LFO 10 and LFO 5 simulated in GB HCT and GB OBC1 models as determined from clustering analysis and centroid classification. As the chain length increased, the population of the helix-like structures tended to increase. These results may suggest that LFOs have tendencies to form helices as their chain lengths are extended.

Hydrogen bonds important for the formation of helix-like structures
To elucidate the hydrogen bonds important for the formation of helix-like structures, the occurrence  15 and LFO 10 with the occurrence frequencies of at least 1 % were analyzed. For the systems simulated in the GB HCT model, the O6 (i) -H3O (i+1) hydrogen bonds (between residue i and i + 1) were found with the highest frequency, and their glycosidic oxygens acted as important hydrogen bond acceptors that interacted with the hydroxyl groups of C3 atoms of the furanose rings and probably helped stabilize the helix-like structures ( Table 2 and Fig. 4). The hydrogen bonds with the second and third highest occurrence frequencies for both LFO 15 and LFO 10 were the O1 (i) -H3O (i) and O5 (i) -H1O (i) hydrogen bonds, which were the hydrogen bonds within the same residue (Table 2 and Fig. 4). The trends of the occurrence frequencies of the hydrogen bonds of LFO 15 and LFO 10 in the GB OBC1 model were also similar to those in the GB HCT model ( Table 2). These three hydrogen bonds (O6 (i) -H3O (i+1), O1 (i) -H3O (i) and O5 (i) -H1O (i) hydrogen bonds), especially the O6 (i) -H3O (i+1) hydrogen bond that was found with the highest frequency, are probably important for the formation of helix-like structures of LFO 15 and LFO 10 as their occurrence frequencies are higher than other hydrogen bonds.

Conformational flexibilities
To investigate the conformational flexibilities of LFO 15, LFO 10 and LFO 5 , the occurrence frequencies of ω, ψ and ɸ of all glycosidic bonds were measured. For the systems simulated in the GB HCT model, ψ and ɸ of all glycosidc linkages of all LFOs exhibited single major peaks around 173°and -63°, respectively (Fig. 5). However, ω was more flexible than ψ and ɸ as it exhibited one major peak and two minor peaks (Fig. 5). The results from the systems simulated in the GB OBC1 model were similar (Additional file 1: Figure S8); ω exhibited more peaks and was more flexible than ψ and ɸ. These results suggest that the flexibility of ω may be responsible for the conformational diversity of LFOs since this dihedral angle has more possibilities in rotating and changing the conformations of LFOs.

Conclusions
To elucidate the structural and molecular properties of LFOs as well as the relationship between these properties and their chain lengths, REMD were performed on systems of LFO 5 , LFO 10 and LFO 15 in GB HCT and GB OBC1 solvent models. We found that as the chain length increased, the radii of gyration tended to increase, suggesting the extension of the conformations as the chain length increases. After clustering analysis and centroid classifications, four major representative conformations (helix-like, partial helix, zig-zag and random structures) were found for LFO 15 and LFO 10, while two conformations (partial helix and random structures) were identified for LFO 5 . The free energy maps show that the four conformations of LFO 15 and LFO 10 were characterized by their upper-middle and lower-middle torsions, whereas the two conformations of LFO 5 were characterized by their molecular angles and middle torsions. As the chain length increased from 5 to 15 residues, the conformation populations of the helix-like structures tended to increase, suggesting the possible tendency of LFOs to form helices as their chain lengths are extended. Moreover, the O6 (i) -H3O (i+1) hydrogen bond was found with the highest frequency, suggesting its importance in helix formation of LFO 15 and LFO 10 . Furthermore, ω was found to be more flexible than ψ and ɸ and probably responsible for the conformational diversity of LFOs. This study gives important insights into the structural and molecular properties of LFOs; they tend to form helical structures as the chain length increases from 5 to 15 residues. Our findings may be useful in the selection of LFOs with appropriate chain lengths and structural properties for pharmaceutical and biological applications.