Molecular dynamics provides insight into how N251A and N251Y mutations in the active site of Bacillus licheniformis RN-01 levansucrase disrupt production of long-chain levan.

Produced by levansucrase, levan and levan oligosaccharides (GFn) have potential applications in food and pharmaceutical industries such as prebiotics, anti-tumor and anti-inflammatory agents. Previous study reported that Bacillus licheniformis RN-01 levansucrase could produce levan oligosaccharides and long-chain levan. However, its N251A and N251Y mutants could effectively produce short-chain oligosaccharides upto GF3, but they could not produce long-chain levan. We hypothesized that these mutations probably reduced GF3 binding affinity in levansucrase active site that contains fructosyl-Asp93 intermediate and caused GF3 to be in an unfavorable orientation for transfructosylation; therefore, levansucrase could not effectively extend GF3 by one fructosyl residue to produce GF4 and subsequently long-chain levan. However, these mutations probably did not significantly reduce binding affinity or drastically change orientation of GF2; therefore, levansucrase could still extend GF2 to produce GF3. Using this hypothesis, we employed molecular dynamics to investigate effects of these mutations on GF2/GF3 binding in levansucrase active site. Our results reasonably support this hypothesis as N251A and N251Y mutations did not significantly reduce GF2 binding affinity, as calculated by MM-GBSA technique and hydrogen bond occupations, or drastically change orientation of GF2 in levansucrase active site, as measured by distance between atoms necessary for transfructosylation. However, these mutations drastically decreased GF3 binding affinity and caused GF3 to be in an unfavorable orientation for transfructosylation. Furthermore, the free energy decomposition and hydrogen bond occupation results suggest the importance of Arg255 in GF2/GF3 binding in levansucrase active site. This study provides important and novel insight into the effects of N251A and N251Y mutations on GF2/GF3 binding in levansucrase active site and how they may disrupt production of long-chain levan. This knowledge could be beneficial in designing levansucrase to efficiently produce levan oligosaccharides with desired length.

Introduction

Levan and levan oligosaccharides (GFn) are natural fructans that contain one terminal glucopyranosyl residue and D-fructofuranosyl repeating unit linked by β-(2, 6) linkage in a main chain with some possible branching points linked by β-(2, 1) linkages [1] (Fig 1A). Properties of levan and levan oligosaccharides depend on their lengths and branching degrees [2], and they have various beneficial properties such as high-water solubility [3] and low intrinsic viscosity [4] for food, cosmetics and pharmaceutical industries. In the food industry, levan and levan oligosaccharides can be used as a prebiotic ingredient [5], encapsulating agent, emulsifier, thickener [3] and cholesterol lowering agent [6]. They can also be used as a component in cosmetics to alleviate skin irritation and moisturize skin [7]. For pharmaceutical industry, they could potentially be used as anti-tumor, anti-inflammatory and anti-viral agents [8].

Fig 1

Levan, proposed reaction mechanism and levansucrase.

(A) An example of levan structure with no branching. (B) Proposed reaction mechanism of levansucrase; (i) sucrose binds in the active site of levansucrase, (ii) the fructosyl-Asp93 intermediate is formed, while glucose is leaving the active site, (iii) sucrose binds at the acceptor binding site. The β-(2, 6) linkage is formed between sucrose and the fructosyl residue of the fructosyl-Asp93 intermediate to extend the levan chain. The bond between the fructosyl residue and Asp93 is broken, and the product is released. (C) Homology model of Bacillus licheniformis RN-01 levansucrase with the fructosyl-Asp93 intermediate.

Levan and levan oligosaccharides are synthesized by levansucrase that is mainly found in various microorganism, including Bacillus subtilis [9], Rahnella aquatilis [10], Zymomonas mobilis [11], Leuconostoc mesenteroides [12] and Leuconostoc citreum [1]. The mechanism of levansucrase was proposed to occur in two steps (Fig 1B) [13]. In the first step, sucrose is bound in the active site, and the fructosyl residue is stabilized by the transition state stabilizer (Asp256). The acid-base catalyst, Glu351, acts as a general acid, protonating the glycosidic oxygen of sucrose. Glucose is released, and oxocabenium ion of the fructosyl residue is formed. Then, a nucleophile (Asp93) attacks C2 of the oxocabenium ion, and the covalent fructosyl-enzyme intermediate is formed. In the second step, sucrose binds in the acceptor binding site. Glu351 acts as a general base that removes a proton from O6 of the non-reducing end of the acceptor. Then, this O6 attacks the fructosyl C2 of the covalent fructosyl-enzyme intermediate, creating the β-(2, 6) linkage to extend the levan chain. The bond between the fructosyl residue and Asp93 is broken, and the product is released [13].

Previous experimental study by Nakapong reported that Bacillus licheniformis RN-01 levansucrase could produce levan oligosaccharides and long-chain levan. However, its N251A and N251Y mutants could effectively produce short-chain oligosaccharides upto GF3, but they could not produce long-chain levan at 323 K and pH 6 [14]. In other words, the wild type could extend GF2 and GF3 by one fructosyl residue to produce GF3 and GF4, respectively, while the mutants could extend GF2 by one fructosyl residue to produce GF3, but they could not effectively extend GF3 to produce GF4 and subsequently long-chain levan. However, the molecular-level understanding on how these mutations cause production of short-chain products is lacking.

Molecular dynamics (MD) simulations is widely used to elucidate conformational changes of enzymes over a period of time and to gain insights into the interactions between enzymes and substrates that may not be accessible by experiments [15-20]. This method can also be used to calculate the binding free energy of ligand binding in macromolecules [17, 21]. However, to our knowledge, MD technique has not been employed to investigate levan oligosaccharides binding in the active site of levansucrase and the effects of mutations on the binding of these substrates.

In this study, MD simulations were performed at experimental temperature and pH on six complexes: GF2-wild-type levansucrase (GF2-LSwt), GF2-N251A mutant levansucrase (GF2-LSN251A), GF2-N251Y mutant levansucrase (GF2-LSN251Y), GF3-wild-type levansucrase (GF3-LSwt), GF3-N251A mutant levansucrase (GF3-LSN251A) and GF3-N251Y mutant levansucrase (GF3-LSN251Y) to elucidate the effects of N251A and N251Y mutations on the binding of GF2/GF3 in the active site of Bacillus licheniformis RN-01 levansucrase that contains the fructosly-Asp93 intermediate (fru-Asp93). This molecular-level understanding on GF2/GF3 binding in Bacillus licheniformis RN-01 levansucrase might be beneficial for designing mutants that can produce levan oligosaccharides with desired lengths.

Materials and methods

Structure preparation

The structures of GF2 and GF3 were constructed using the LEaP module in AMBER14 [22] and the GLYCAM06j-1force field parameters [23]. To remove unfavorable interactions, these structures were minimized by 2,500 steps of steepest descent and 2,500 steps of conjugate gradient. The target sequence of Bacillus licheniformis RN-01 levansucrase (GenBank ID: ACI15886.1) was obtained from the National Center for Biotechnology Information (NCBI). SWISS-MODEL server [24-27] was used to construct the homology model of levansucrase from Bacillus licheniformis RN-01 based on the crystal structure of Bacillus subtilis levansucrase (PDB ID: 1OYG [28]), which has the highest sequence identity to the target sequence. The quality of the homology model was evaluated by Ramachandran plot produced by the RAMPAGE server [29]. S1 Fig shows that a majority of its residues are in favored region (96.0%) and allowed region (3.3%), indicating reasonable quality of this homology model. Moreover, the catalytic residues (Asp93, Asp256 and Glu351) of this homology model were found in positions, where they should be able to catalyze the transfructosyation (Fig 1C). All ionizable amino acids were protonated at pH = 6, using the H++ server [30]. To construct the structure of fru-Asp93, the initial structure of Asp86 and fructosyl residue were taken from the crystal structure of Bacilus subtilis levansucrase in complex with sucrose (PDB ID: 1PT2 [28]). GaussView05 program [31] was used to create a bond between OD2 of Asp and C2 of the fructosyl residue. The atomic charges and the electrostatic potential (ESP) charges of fru-Asp93 were calculated using the HF/6-31G* basis set in the Gaussian09 program [32]. Using Antechamber module in AMBER14, the ESP charges of fru-Asp93 intermediate was converted into restrained ESP (RESP) charges, and other force filed parameters of fru-Asp93 intermediate were generated from general AMBER force field (GAFF). The LEaP module was then used to construct the structure of levansucrase with fru-Asp93 in its active site, using ff14SB force field (Fig 1C).

Identification of catalytically competent binding conformations and molecular dynamics

To determine whether Autodock vina [33] and its parameters were appropriate for the studied systems, the crystal sucrose was redocked into the active site of the crystal structure of Bacillus subtilis levansucrase (1PT2). The best docked and crystal binding conformations were compared and found to be reasonably similar with the RMSD value of 0.64 Å (S2 Fig), indicating that Autodock Vina and its parameters were appropriate for this system. To determine catalytically competent binding conformations, Autodock Vina was employed to dock GF2/GF3 in the active site of the homology model of wild-type levansucrase with fru-Asp93 to create GF2-LSwt and GF3-LSwt complexes. A grid box of 40 Å x 40 Å x 40 Å with a grid spacing of 1 Å was employed. 20 independent docking runs were performed for each ligand, where each run gave nine possible binding conformations. In order for the wild type to be able to extend the levan chain, GF2/GF3 should bind in catalytically competent orientations, where O6 of the non-reducing end of GF2/GF3 turns toward C2 of the fructosyl residue of fru-Asp93. Employing this assumption, only binding conformations that have O6 of the non-reducing end of GF2/GF3 turns toward C2 of the fructosyl residue of fru-Asp93 were selected. The binding conformations that passed this criterion were later clustered by MMTSB tool set [34] based on their structural similarities as measured by the RMSD values of heavy atoms. To identify a reasonable representative binding conformation of each cluster, a binding conformation that is most similar to the average structure of all members of each cluster was chosen to be a centroid. The centroid of each cluster was immersed in an isomeric truncated octahedral box of TIP3P water molecules with the buffer distance of 13 Å using the LEaP module. Chloride ions (Cl-) were added to neutralize all systems. To reduce unfavorable interactions, the complexes were minimized with the five step procedure. All steps include 5,000 steps of steepest descent and 5,000 steps of conjugate gradient with different restraints on the proteins. Initially, to relax each system, the hydrogen atoms and water molecules were minimized, while heavy atoms of protein were restrained with a force constant of 5 kcal/ (mol Å2). The backbone of the protein was subsequently restrained with force constants of 10, 5 and 1 kcal/ (mol Å2), respectively. Finally, the entire system was minimized without any restraining force. The GPU (CUDA) version of PMEMD module of the AMBER14 was employed to simulate all systems under the periodic boundary condition [35-37]. The SHAKE algorithm [38] was used to constrain all bonds involving hydrogen atoms, allowing a simulation time step of 0.002 ps. A cutoff distance of 12 Å was used for non-bonded interactions, and the particle mesh Ewald method was applied to calculate the long-range electrostatic interaction [38]. The Langevin dynamics technique [39] was employed to control the temperature with a collision frequency of 1.0 ps-1. All systems were heated from 0 K to the experimental temperature of 323 K (50°C) for 200 ps in the NVT ensemble, while the backbone of proteins were restrained with a force constant of 10 kcal/ (mol Å2). Subsequently, all systems were equilibrated for 300 ps with no restraint in the NVT ensemble. These systems were further simulated for 80 ns in the NPT ensemble at 323 K and 1 atm. With the assumption that catalytically competent binding conformations should have the position of O6 of the non-reducing end of GF2/GF3 that is not too far from that of C2 of the fructosyl residue of fru-Asp93 after simulations, the distances between O6 of the non-reducing end of GF2/GF3 and C2 of the fructosyl residue of fru-Asp93 (O6-C2 distance) of all centroids were measured. The centroids with the O6-C2 distances greater than 5 Å were eliminated. One centroid of GF2 binding conformations and one centroid of GF3 binding conformations passed this criterion. Since it may still be possible that the selected centroid may not necessary be the most stable binding conformation of the cluster, similar setup, minimization and MD procedure were also performed on two additional binding conformations of GF2/GF3 that are in the same cluster as the selected centroid. These binding conformations were second and third most similar to the average structure of all members of each cluster. With the assumption that catalytically competent binding conformations of GF2/GF3 in the active site of the wild type should be the ones, where GF2/GF3 stably binds in its active site, the binding conformations, whose heavy-atom RMSD of GF2/GF3 during 60–80 ns (the last 20 ns of the simulation) have the lowest values and fluctuation out of the three binding conformations, were chosen to be the catalytically competent binding conformations of the wild type complex (GF2-LSwt and GF3-LSwt complexes). To construct the structures of the mutant complexes, Asn251 of the selected binding conformations of GF2-LSwt and GF3-LSwt complexes were mutated to Ala251 to build GF2-LSN251A and GF3-LSN251A complexes, and it was mutated to Tyr251 to build GF2-LSN251Y and GF3-LSN251Y complexes. Similar setup, minimization and MD were performed on the mutant systems.

In terms of analyses, the RMSD values with respect to the minimized structure were calculated to monitor the stability of all systems. Since the RMSD values of all systems were stable around 60–80 ns, these trajectories were used for further analyses. To measure the proximity between atoms necessary for transfructosylation, the O6-C2 distances of all systems were measured. To measure binding affinity between GF2/GF3 and levansucrase, total binding free energies and decomposition of free energies per residue were calculated by Molecular Mechanics/Generalized Born Surface Area (MM/GBSA) method. MM/GBSA technique [40, 41] is widely employed to approximate the binding affinities, as calculated by binding free energies, of small ligands to macromolecules [42]. This method is stable, reproducible [42] and giving promising results in correctly ranking the molecules with known affinity to their target proteins [43-50]. This technique was also employed for rigorous free energy decomposition into contributions from different groups of atoms or types of interaction in various studies to determine important binding residues [51-54].

Hydrogen bond interactions between GF2/GF3 and levansucrase were determined by calculating hydrogen bond occupations between amino acid residues and GF2/GF3. In this study, a hydrogen bond occurred if the following criteria were met: (i) a proton donor-acceptor distance ≤ 3.5 Å and (ii) a donor-H-acceptor bond angle ≥ 120°. Strong and medium hydrogen bonds were defined as hydrogen bonds with occupation > 75% and 50–75%, respectively. Weak hydrogen bonds were defined as hydrogen bonds with occupation < 50% but ≥ 25%.

Results and discussion

System stability

Using the minimized structures as references, the RMSD values of all atoms, backbone atoms and ligand atoms of all systems were calculated to determine the stabilities of these systems and identify appropriate trajectories for further analyses (Fig 2). As shown by these plots, the simulations of all systems were likely to reach equilibrium around 80 ns. As a result, the 60–80 ns trajectories of all systems were employed for further analyses.

Fig 2

RMSD plots of A) GF2-LSwt, B) GF2-LSN251A, C) GF2-LSN251Y, D) GF3-LSwt, E) GF3-LSN251A and F) GF3-LSN251Y complexes. The RMSD values of all atoms, backbone atoms and ligand atoms are shown in red, black and green, respectively.

The proximity between atoms necessary for transfructosylation

With the assumption that the system that allow transfructosylation to occur should be the one that has O6 of the non-reducing end of GF2/GF3 turning toward C2 of the fructosyl residue of fru-Asp93, and these two atoms should not be too far from each other, the O6-C2 distances of all systems were measured as shown in Fig 3 and S3 Fig. The O6-C2 distances of GF2-LSwt, GF2-LSN251A and GF2-LSN251Y are reasonable and quite stable during the 60–80 ns simulations. These findings suggest that transfructosylation should be able to occur in these systems, i.e., the wild type, the N251A and N251Y mutants should be able to extend GF2 by one fructosyl residue to create GF3. These results support the previous experimental findings that the wild type, the N251A and N251Y mutants could produce GF3.

Fig 3

The distance between O6 of the non-reducing end of GF2/GF3 and C2 of the fructosyl residue of fru-Asp93: A) GF2-LSwt, B) GF2-LSN251A, C) GF2-LSN251Y, D) GF3-LSwt, E) GF3-LSN251A and F) GF3-LSN251Y complexes.

Superimpositions between the crystal structure of Erwinia amylovora levansucrase in complex with fructose and glucose (PDB ID: 4D47 [55]) and the homology model of Bacillus licheniformis RN-01 levansucrase with docked sucrose, between 4D47 and the homology model of Bacillus licheniformis RN-01 levansucrase containing fru-Asp93 intermediate as well as between 4D47 and the homology model of Bacillus licheniformis RN-01 levansucrase containing fru-Asp93 intermediate with catalytically competent binding conformation of GF2 are shown Fig 4. These results show that the fructosyl reside of sucrose and fructosyl residue of fru-Asp93 are in similar position to that of fructose in the crystal structure of Erwinia amylovora levansucrase (Fig 4A and 4B). Their orientations are slightly different probably because fructose in the crystal structure of Erwinia amylovora levansucrase is the hydrolysis product of sucrose and it does not connect to other residue; therefore, it has more flexibility in terms of orientation than the fructosyl residue of sucrose/fru-Asp93. Moreover, Fig 4C shows that the position of the fructosyl residue of the non-reducing end of GF2 is close to that of glucose in the crystal structure of Erwinia amylovora levansucrase.

Fig 4

Superimpositions between the crystal structure of Erwinia amylovora in complex with fructose and glucose (grey) and (A) the homology model of Bacillus licheniformis RN-01 levansucrase (orange) with docked sucrose (black), (B) the homology model of Bacillus licheniformis RN-01 levansucrase containing fru-Asp93 intermediate (orange) and (C) the homology model of Bacillus licheniformis RN-01 levansucrase containing fru-Asp93 intermediate (orange) with catalytically competent binding conformation of GF2 (black).

For the systems involving GF3, the O6-C2 distance of GF3-LSwt is reasonable and stable, suggesting that transfructosylation should be able to occur. However, N251A and N251Y mutations drastically increase the O6-C2 distances to around 10.8 Å and 8.6 Å for GF3-LSN251A and GF3-LSN251Y complexes, respectively. For these two mutant complexes, O6 of the non-reducing end of GF3 is too far from C2 of the fructosyl residue of fru-Asp93 for transfructosylation to occur. Moreover, O6 of the non-reducing end of GF3 also points away from C2 of the fructosyl residue of fru-Asp93. Therefore, the orientations of GF3 in these two mutant systems are not favorable for transfructosylation to occur, and these mutants should not be able to effectively extend GF3 by one fructosyl residue to produce GF4. These findings support the previous experimental results that the wild type could produce GF4, while the N251A and N251Y mutants could not effectively produce GF4 and long-chain levan.

Binding free energies

To determine if binding affinity is an important factor associated with the experimental results that the wild type, N251A and N251Y mutants could extend GF2 to produce GF3, while only the wild type could effectively extend GF3 to produce GF4, MM-GBSA method was employed to calculate the binding free energies of GF2-LSWT, GF2-LSN251A, GF2-LSN251Y, GF3-LSwt, GF3-LSN251A and GF3-LSN251Y complexes during the 60–80 ns trajectories (Table 1). The binding free energies of GF2-LSwt, GF2-LSN251A and GF2-LSN251Y complexes are -4.7±0.9, -10.1±0.8 and -4.5±0.8 kcal/mol, respectively. These results suggest that N251A and N251Y mutations did not significantly reduce the binding affinities of GF2 in the active site of levansucrase. Since the distances between atoms necessary for transfructosylation of these systems are also reasonable, the wild type, N251A and N251Y mutants should all be able to bind GF2, extend it by one fructosyl residue and produce GF3, supporting the experimental results. In terms of GF3 binding, the binding free energies of GF3-LSWT, GF3-LSN251A and GF3-LSN251Y complexes are -20.5±0.7, 1.1±0.9 and -8.7±0.8 kcal/mol, respectively. These results show that the binding affinities of GF3-LSN251A and GF3-LSN251Y complexes are worse than that of the GF3-LSWT complex, suggesting that these mutations reduce the binding affinities of GF3 in the active site of the mutants as compared to that of the wild type. Since the distance between atoms necessary for transfructosylation is reasonable only for the wild-type complex, these results suggest that only the wild type could potentially bind GF3, extend it by one fructosyl residue and produce GF4, while the N251A and N251Y mutants could not tightly bind GF3 to effectively produce GF4, supporting the experimental results.

Table 1

The Binding free energies (kcal/mol) and their components of GF2-LSwt, GF2-LSN251A, GF2-LSN251Y, GF3-LSwt, GF3-LSN251A, and GF3-LSN251Y complexes.

System	ΔEvdW	ΔEele	ΔGpol	ΔGnp	a)ΔGsolv	-TΔStot	b)ΔGbind	Standard error of the mean of ΔGbind
GF2-LSwt	-35.3	-55.9	66.4	-5.6	60.8	25.7	-4.7	0.9
GF2-LSN251A	-38.2	-66.0	71.8	-5.7	66.1	28.0	-10.1	0.8
GF2-LSN251Y	-38.0	-47.8	60.3	-5.6	54.7	26.5	-4.5	0.8
GF3-LSwt	-47.4	-103.3	101.1	-7.9	93.2	37.0	-20.5	0.7
GF3-LSN251A	-45.7	-36.3	62.6	-6.2	56.4	26.7	1.1	0.9
GF3-LSN251Y	-54.7	-50.3	74.7	-7.4	67.3	29.0	-8.7	0.8

a) ΔGsolv = ΔGpol + ΔGnp

b) ΔGbind = ΔEvdW + ΔEele + ΔGsolv—TΔStot

Levansucrase from Gram-positive bacteria generally produce long-chain levan polymer, while that from Gram-negative bacteria produce short-chain levan oligosaccharides [55]. Bacillus licheniformis RN-01 levansucrase is from Gram-positive bacteria; therefore, it usually produces long-chain levan polymer. However, its N251A and N251Y mutants could effectively produce short-chain oligosaccharides, and they could not produce long-chain levan. The multiple sequence alignment of B. licheniformis RN-01 levanscurase and levansucrase from Gram-positive bacteria such as B. megaterium [56], B. amyloliquefaciens [57], B. atrophaeus [58], B. subtilis [9], and B. stearothermophilus [59], and from Gram-negative bacteria such as G. diazotrophicus [60], Z. mobilis [11], P. chlororaphis [61], R. aquatilis [62] and E. amylovora [63] shows that Asn251 of levansucrase from Gram-positive bacteria is generally conserved, while that from Gram-negative bacteria is mutated to other residues such as Phe or Tyr (Fig 5 and S4 Fig). These mutations could potentially reduce the binding affinity to substrates such as GF3 in the active site of levansucrase from Gram-negative bacteria and could potentially disrupt the production of long-chain levan polymer, similar to the N251A and N251Y mutations of Bacillus licheniformis RN-01 levansucrase.

Fig 5

The multiple sequence alignment of the sequence near Asn251 of B. licheniformis RN-01 levanscurase and levansucrase from Gram-positive bacteria such as B. megaterium [56], B. amyloliquefaciens [57], B. atrophaeus [58], B. subtilis [9], and B. stearothermophilus [59], and from Gram-negative bacteria such as G. diazotrophicus [60], Z. mobilis [11], P. chlororaphis [61], R. aquatilis [62] and E. amylovora [63].

Moreover, PDBeFold [64] was employed to perform structural superposition of our model with available crystal structures. S1 Table shows crystal structures of enzymes that are most similar to our model. Fig 6 shows the superimposition of our model with most similar levansucrase structures from B. sutilis (PDB ID: 1OYG [28]) and B. megaterium (PDB ID: 3OM2 [65]). These results show that the overall structure of our model is similar to the crystal structures of levansucrase from B. sutilis and B. megaterium. Additionally, the positions and orientations of N251 of our model and these most similar structures are very similar.

Fig 6

The superimposition of the homology model of Bacillus licheniformis RN-01 levansucrase (grey) with most similar levansucrase structures from B. sutilis (PDB ID: 1OYG [28], green) and B. megaterium (PDB ID: 3OM2 [65], blue) from PDBeFold.

In terms of binding free energy components of GF2 binding, the main components contributing to the substrate binding affinities of GF2-LSwt, GF2-LSN251A, GF2-LSN251Y complexes are the electrostatic interaction terms (ΔEele) as they have the most favorable values that are in the range of -66.0 –-47.8 kcal/mol. Other terms favor GF2 binding are the van der Waals energy terms (ΔEvdw), which are in the range of -38.2 –-35.3 kcal/mol, and the non-polar solvation terms (ΔGnp), which are in the range of -5.7 –-5.6 kcal/mol. The polar solvation terms (ΔGpol) have unfavorable contribution to GF2 binding, and they are in the range of 60.3–71.8 kcal/mol.

In terms of GF3 binding, the main component contributing to the substrate binding affinity of GF3-LSwt complex is ΔEele with the value of -103.3 kcal/mol. ΔEvdw and ΔGnp are also favorable with the values of -47.4 and -7.9 kcal/mol, respectively. ΔGpol is unfavorable with the value of 101.1 kcal/mol. However, the main component contributing to GF3 binding affinities of GF3-LSN251A and GF3-LSN251Y complexes are ΔEvdw with the values of -45.7 and -54.7 kcal/mol, respectively. The values of ΔEele of GF3-LSN251A (-36.3 kcal/mol) and GF3-LSN251Y (-50.3 kcal/mol) complexes are significantly worse than that of GF3-LSwt complex (-103.3 kcal/mol). These results were probably caused by the fact that GF3 could not bind in a favorable orientation in the active sites of the N251A and N251Y mutants. In this case, GF3 could not form as many favorable interactions with residues in the active sites of the mutants as with those of the wild type. Additionally, ΔGnp of GF3-LSN251A and GF3-LSN251Y are favorable with the values of -6.2 and -7.4 kcal/mol. Their ΔGpol values are unfavorable with the values of 62.6 and 74.7 kcal/mol.

Per residue substrate-enzyme interactions

To identify important binding residues that make major contributions to the calculated binding free energies as well as the effects of N251A and N251Y mutations on the binding residues, the values of free energy decomposition on a per residue basis () were calculated as shown in Fig 7. In this study, an importance binding residue was defined to be a residue with the total energy contribution better than -1.0 kcal/mol. For GF2-LS complexes, residues with energy contribution better than -1 kcal/mol for all three complexes are Trp92, fru-Asp93, Val123, Arg369 and Arg442, indicating their importance in GF2 binding in the active sites of wild-type and mutant levansucrase. However, there are also residues with total energy contribution better than -1 kcal/mol in the wild-type complex, but not in the mutant complexes such as Trp170, Arg255 and Glu351, suggesting their importance in GF2 binding only in the active site of wild-type levansucrase. For GF3-Ls complexes, Trp92, fru-Asp93, Trp170, Asn441 and Arg442 have energy contribution better than -1 kcal/mol for all three systems, suggesting their importance in GF3 binding in the active sites of wild-type and mutant levansucrase. Thr126, Gln168, Arg255, Arg369 and Tyr438 have energy contribution better than -1 kcal/mol in the wild-type complex but not in the mutant complexes, suggesting their importance in GF3 binding only in the active site of wild-type levansucrase.

Fig 7

Per-residue decomposition of binding free energy contributions of A) GF2-LSwt, B) GF2-LSN251A, C) GF2-LSN251Y, D) GF3-LSwt, E) GF3-LSN251A and F) GF3-LSN251Y complexes.

In terms of the mutated residue 251, we found that the N251A and N251Y mutations did not cause significant changes to the total energy contribution of residue 251. However, these mutations caused significant changes to the total energy contributions of other residues, especially Arg255. For GF2-LS complexes, the value of the total energy contribution of Arg255 was changed from -4.7 kcal/mol in the wild-type complex to -0.1 kcal/mol in the N251A mutant complex and to -0.3 kcal/mol in the N251Y mutant complex. For GF3-LS complexes, the value of the total energy contribution of Arg255 was changed from -3.5 kcal/mol in the wild-type complex to -0.4 kcal/mol in the N251A mutant complex and to -1.2 kcal/mol in the N251Y mutant complex. These results suggest the importance of Arg255 in GF2/ GF3 binding.

Hydrogen bond interactions

To identify hydrogen bonds important for GF2/GF3 binding, hydrogen bond occupations of all systems were calculated as shown in Table 2, S2 and S3 Tables. In terms of GF2 binding, the N251A and N251Y mutations did not reduce the number of strong and medium hydrogen bonds between GF2 and binding residues. On the contrary, the total number of strong and medium hydrogen bonds of the mutant complexes is slightly more than that of the wild-type complex. These results suggest that these mutations may not significantly reduce the binding affinity of GF2, supporting the binding free energy results of GF2. In terms of GF3 binding, the N251A and N251Y mutations drastically reduce the number of strong and medium hydrogen bonds between GF3 and binding residues. These results suggest that these mutations probably cause the reduction of GF3 binding affinity, supporting the binding free energy results of GF3.

Table 2

Number of strong and medium hydrogen bonds formed between GF2/GF3 and binding residues in the GF2-LSwt, GF2-LSN251A, GF2-LSN251Y, GF3-LSwt, GF3-LSN251A, and GF3-LSN251Y complexes.

Complex	Number of strong and medium hydrogen bond	Binding residues that form hydrogen bonds with GF2/GF3
GF2-LSwt	5 (4S*, 1M**)	Arg255, Glu349, Glu351
GF2-LSN251A	6 (4S, 2M)	fru-Asp93, Val123, Arg369, Arg442
GF2-LSN251Y	6 (4S, 2M)	fru-Asp93, Val123, Gln168, Tyr413, Arg442
GF3-LSwt	14 (9S, 5M)	Trp92, fru-Asp93, Thr126, Arg255, Glu349, Glu351, Arg369, Tyr438, Arg442
GF3-LSN251A	5 (3S, 2M)	Trp92, fru-Asp93, Glu351, Tyr413, Arg442
GF3-LSN251Y	4 (2S, 2M)	Trp92, fru-Asp93, Arg442

*S; Strong hydrogen bond

**M; medium hydrogen bond

To determine the importance of Arg255 in GF2/GF3 binding in the active site of wild-type levansucrase, hydrogen bond networks involving Arg255 were identified as shown in Fig 8. At the beginning of the 80 ns MD simulations, GF2/GF3 formed hydrogen bond networks with Asn251, Glu349 and Arg255 in the wild-type complexes. The N251A and N251Y mutations disrupt these hydrogen bond networks in the mutant complexes; therefore, Arg255 could not effectively form hydrogen bonds with GF2/GF3 during the 80 ns MD simulations. However, there were other residues, instead of Asn251, Glu349 and Arg255, that later formed hydrogen bonds with GF2, still keeping it in a reasonable binding affinity and orientation for transfructosylation (Fig 3). These binding residues are fru-Asp93, Val123, Arg369, Arg442 for the N251A mutant, and fru-Asp93, Val123, Gln168, Tyr413 and Arg442 for the N251Y mutant (Table 2). These residues are different from the binding residues that formed hydrogen bonds with GF2 in the wild-type complex, and the total number of strong and medium hydrogen bonds formed between GF2 and the binding residues in the N251A and N251Y mutants is slightly more than that of the wild type. As a result, the binding conformations of GF2 in the active site of the mutants are slightly different from that of the wild type, but they are still in reasonable orientations and distances for transfructosylation. In terms of GF3 binding, the binding conformations of GF3 in the active sites of the N251A and N251Y mutants are drastically different from that of the wild type (Fig 3). Arg255 could not effectively form hydrogen bonds with GF3, and there are significantly less number of strong and medium hydrogen bonds formed between GF3 and the binding residues of the N251 A and N251Y mutants than that of the wild type (Table 2). Therefore, GF3 was not able to bind and stay in a favorable orientation for transfructosylation (Fig 3).

Fig 8

Hydrogen bond networks involving Arg255 and GF2/GF3 in A) GF2-LSwt, B) GF2-LSN251A, C) GF2-LSN251Y, D) GF3-LSwt, E) GF3-LSN251A and F) GF3-LSN251Y complexes at the beginning of the 80 ns MD simulations.

Conclusions

In this work, MD was performed on the GF2-LSwt, GF2-LSN251A, GF2-LSN251Y, GF3-LSwt, GF3-LSN251A and GF3-LSN251Y complexes to gain insight into the effects of N251A and N251Y mutations on the binding of GF2/GF3 in the active site of Bacillus licheniformis RN-01 levansucrase. Our results of binding free energies and hydrogen bond occupations as well as the distances between atoms necessary for transfructosylation of GF3-LSwt, GF3-LSN251A and GF3-LSN251Y complexes support the hypothesis that these mutations reduced GF3 binding affinity in active site of levansucrase with fructosyl-Asp93 intermediate and caused GF3 to be in an unfavorable orientation for transfructosylation; therefore, transfructosylation could not occur in GF3-LSN251A and GF3-LSN251Y complexes. As a result, only the wild type should be able to extend GF3 by one fructosyl residue to produce GF4, supporting the experimental results that the wild type can produce GF4, but the N251A and N251Y mutants cannot effectively produce GF4. However, these mutations did not drastically change binding affinity or orientation of GF2 as shown by the binding free energy and hydrogen bond occupation results as well as the distances between atoms necessary for transfructosylation of GF2-LSwt, GF2-LSN251A and GF2-LSN251Y complexes. Therefore, the wild type, the N251A and N251Y mutants should be able to extend GF2 by one fructosyl residue to produce GF3, supporting the experimental results that the wild type, the N251A and N251Y mutants can produce GF3. Moreover, the free energy decomposition results also suggest the importance of Arg255 in GF2/GF3 binding in the active site of the wild type. Our results also show that Arg255 formed hydrogen bond networks with GF2/GF3, Asn251 and Glu349 in the wild-type complexes at the beginning of the 80 ns MD simulations, and the N251A and N251Y mutations disrupted these hydrogen bond networks. Although these hydrogen bond networks were disrupted in the GF2-LsN251A and GF2-LsN251Y complexes, GF2 could still bind in a favorable orientation for transfructosylation in the active sites of these mutants probably because there were other residues binding and forming hydrogen bonds with GF2, and these interactions helped prevent misorientation of GF2. However, GF3 could not bind in a favorable orientation for transfructosylation in the active sites of these mutants because there was significantly less number of residues binding and forming hydrogen bonds with GF3 in the mutant complexes than that in the wild-type complex. Our study provides important and novel insight into the binding of GF2/ GF3 in the active site of Bacillus licheniformis RN-01 levansucrase and into how N251A and N251Y mutations may disrupt production of long-chain levan.

Ramachandran plot of the homology model of Bacillus licheniformis RN-01 levansucrase.

(TIF)

Click here for additional data file.

Superimposition of the crystal binding conformation (black) and best docked conformation (pink).

(TIF)

Click here for additional data file.

The distance between O6 of the non-reducing end of GF2/GF3 and C2 of the fructosyl residue of fru-Asp93 during 0–80 ns: A) GF2-LSwt, B) GF2-LSN251A, C) GF2-LSN251Y, D) GF3-LSwt, E) GF3-LSN251A and F) GF3-LSN251Y complexes.

(TIF)

Click here for additional data file.

The multiple sequence alignment of B. licheniformis RN-01 levanscurase and levansucrase from Gram-positive bacteria such as B. megaterium [56], B. amyloliquefaciens [57], B. atrophaeus [58], B. subtilis [9], and B. stearothermophilus [59], and from Gram-negative bacteria such as G. diazotrophicus [60], Z. mobilis [11], P. chlororaphis [61], R. aquatilis [62] and E. amylovora [63].

(TIF)

Click here for additional data file.

Structural alignment results of Bacillus licheniformis RN-01 levansucrase using PDBeFold.

(DOCX)

Click here for additional data file.

Hydrogen bond occupations of GF2-LSwt, GF2-LSN251A and GF2-LSN251Y complexes.

(DOCX)

Click here for additional data file.

Hydrogen bond occupations of GF3-LSwt, GF3-LSN251A and GF3-LSN251Y complexes.

(DOCX)

Click here for additional data file.