Enhanced Conformational Sampling of N-Glycans ... - ACS Publications

Apr 11, 2017 - ABSTRACT: Molecular dynamics (MD) simulation of a N- glycan in solution is challenging because of high-energy barriers of the glycosidi...
0 downloads 15 Views 2MB Size
Subscriber access provided by MT ROYAL COLLEGE

Article

Enhanced Conformational Sampling of N-glycans in Solution with Replica State Exchange Metadynamics Raimondas Galvelis, Suyong Re, and Yuji Sugita J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.7b00079 • Publication Date (Web): 11 Apr 2017 Downloaded from http://pubs.acs.org on April 19, 2017

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 free 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 accessible to all readers and 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.

Journal of Chemical Theory and Computation 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 36

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

Journal of Chemical Theory and Computation

Enhanced Conformational Sampling of N-glycans in Solution

with

Replica

State

Exchange

Metadynamics Raimondas Galvelis†, Suyong Re†,§, and Yuji Sugita*,†,#,¶,§ †

RIKEN Theoretical Molecular Science Laboratory, 2-1 Hirosawa, Wako-shi, Saitama 351-

0198, Japan; #RIKEN iTHES, 2-1 Hirosawa, Wako-shi, Saitama 351-0198, Japan; ¶RIKEN Advanced Institute for Computational Science, Integrated Innovation Building 7F, 6-7-1 Minatojima-minamimachi, Chuo-ku, Kobe, Hyogo 650-0047, Japan; §RIKEN Quantitative Biology Center, Integrated Innovation Building 7F, 6-7-1 Minatojima-minamimachi, Chuo-ku, Kobe, Hyogo 650-0047, Japan.

Corresponding Author: Yuji Sugita, RIKEN Theoretical Molecular Science Laboratory, 2-1 Hirosawa, Wako-shi, Saitama 351-0198, Japan, TEL: +81-48-462-1407, FAX: +81-48-4674532, E-mail: [email protected]

KEYWORDS. N-glycan, molecular dynamics, enhanced conformational sampling, collective variables, conformation clustering

ACS Paragon Plus Environment

1

Journal of Chemical Theory and Computation

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

Molecular dynamics (MD) simulation of a N-glycan in solution is challenging due to highenergy barriers of the glycosidic linkages, functional group rotational barriers, and numerous intra- and inter-molecular hydrogen bonds. In this study, we apply different enhanced conformational sampling approaches, namely, metadynamics (MTD), the replica-exchange MD (REMD), and the recently proposed replica state exchange MTD (RSE-MTD), to a N-glycan in solution and compare their conformational sampling efficiencies. MTD helps to cross the highenergy barrier along the ω angle by utilizing a bias potential, but it cannot enhance sampling of the other degrees of freedom. REMD ensures moderate-energy barrier crossings by exchanging temperatures between replicas, while it hardly crosses the barriers along ω. In contrast, RSEMTD succeeds to cross the high-energy barrier along ω as well as to enhance sampling of the other degrees of freedom. We tested two RSE-MTD schemes: in one scheme, 64 replicas were simulated with the bias potential along ω at different temperatures, while simulations of 4 replicas were performed with the bias potentials for different CVs at 300 K. In both schemes, one unbiased replica at 300 K was included to compute conformational properties of the glycan. The conformational sampling of the former is better than the other enhanced sampling methods, while the latter shows reasonable performance without spending large computational resources. The latter scheme is likely to be useful when a N-glycan-attached protein is simulated.

ACS Paragon Plus Environment

2

Page 2 of 36

Page 3 of 36

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

Journal of Chemical Theory and Computation

1. INTRODUCTION Glycans consist of rigid saccharide units connected together through flexible glycosidic linkages, enabling a variety of conformers in solution.1 In contrast to the classical 'lock-and-key' mechanisms for ordinary ligand binding to proteins, glycans provide multiple 'keys' and one of them is used for specific glycan-receptor recognition.2-5 The most abundant oligosaccharide found in cellular proteins is a N-glycan, which covalently bonds to a protein via the nitrogen atom of asparagine (Asn) side chain (Figure 1a and Figure S1). The 1-6 glycosidic linkage is characterized by three dihedral angles (φ, ψ, and ω; Figure 1b), where the ω angle gives the three distinct rotamer states: gauche-gauche (gg), gauche-trans (gt) and trans-gauche (tg) (Figure 1c). The available NMR data suggest that typical N-glycans has mainly gg and gt conformations in solution.6 It has been estimated that these rotamer states are separated by significant energy barriers (~7 kcal/mol)7. In addition, the multiple hydroxyl groups form various intra- and intermolecular hydrogen-bond networks stabilizing distinct conformations. These chemical features give slow interconversion of N-glycan conformations with different rotamer states. Molecular dynamics (MD) simulations are commonly used to examine conformational fluctuation or dynamics of biomolecules.8 However, slow biological processes with the timescale of milliseconds or longer cannot be easily simulated with brute-force all-atom MD methods.9 One way to overcome the difficulty is to use enhanced conformational sampling approaches.10-11 The replica-exchange MD (REMD) was already applied to sample multiple conformers of a N-glycan in solution or gas phase.12-14 Re et al. simulated multiple conformers of a biantennary complex-type N-glycan (Bi9) in solution, such as two extended, half- and tightbackfolded conformations.13 Nishima et al. examined the conformations of other N-glycans (BiF10, BiB10, and BiBF11) by REMD to study the effect of bisecting GlcNAc and core fucosylation on the biantennary complex-type N-glycans in solution.12 Compared to the

ACS Paragon Plus Environment

3

Journal of Chemical Theory and Computation

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

conventional MD, these studies extended the conformational space of N-glycans, but the freeenergy landscapes of N-glycans in solution were not perfectly converged within the simulation times. An alternative generalized-ensemble approach is the hamiltonian replica-exchange molecular dynamics (H-REMD) (or the replica-exchange umbrella sampling (REUS)).15-16 Patel et al. used a set of two-dimensional grid-based bias potentials in H-REMD to improve the rotamer sampling of the glycosidic linkages.17-18 However, the separate umbrella sampling (US) simulations of model disaccharides are required for each glycosidic linkage to parameterize the bias potentials before the production simulations. In this study, we examine conformational sampling efficiency of several enhanced sampling methods, namely, metadynamics (MTD)19-20, temperature-REMD21, and the recently proposed replica state exchange metadynamics (RSE-MTD)22 for N-glycan simulations in solution. Since each method enhances sampling of a conformational space in a different way, it is meaningful to compare their features in the case of N-glycans. In the next section, we introduce a target Nglycan system and then summarize the theories of MTD, REMD, and RSE-MTD respectively. The simulation results are analyzed in terms of the time-series of collective variables (CVs), the conformational entropy, the free-energy landscapes, and structure clustering of internal coordinates, and so on. We propose two RSE-MTD schemes: one uses 64 replicas with the bias potentials for the same CV at different temperatures and another reduces the number of replicas to 4 with three bias potentials for different CVs at the room temperature. The former shows the best sampling efficiency, while the latter gives a reasonable performance and could be used in simulations of larger N-glycan-attached proteins.

2. METHODS

ACS Paragon Plus Environment

4

Page 4 of 36

Page 5 of 36

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

Journal of Chemical Theory and Computation

2.1. A target molecule and simulation methods A biantennary complex-type N-glycan, named as Bi9, (Figure 1a and Figure S1) contains a branched penta-saccharide core sequence, Manα1-6(Manα1-3)Manβ1-4GlcNAcβ1-4GlcNAcβ1Asn-X-Ser/Thr. This form is shared by all N-glycans. The conformational space of N-glycans is, in general, determined by the flexibility of glycosidic linkages. In Bi9, the α1-6 linkage, which is characterized by three dihedral angles φ, ψ, and ω, largely determines the overall conformational space (Figure 1b). We constructed two initial configurations corresponding to the gg and gt conformations using GLYCAM Web Carbohydrate Builder23 and solvated each with 3437 water molecules using VMD24. Bi9 consists of eight linkages connecting nine saccharides. To characterize its global motions, the eight pseudo-dihedral angles are defined: χ(12): 2O5-2C11C4-1C3, χ(23): 3O5-3C1-2C4-2C3, χ(34): 4O5-4C1-3C4-3C3, χ(45): 5O5-5C1-4C4-4C3, χ(56): 6O5-6C1-5C4-5C3, χ(34'): 4'O5-4'C1-3C6-3C5, χ(4'5'): 5'O5-5'C1-4'C2-4'C1, and χ(5'6'): 6'O5-6'C15'C4-5'C3 (Figure S1). In addition to the ω angle, these pseudo-dihedral angles are used as collective variables and/or for the trajectory analysis. All the simulations were carried out using the GENESIS software package25. The CHARMM36 additive force field for carbohydrates26-28 and TIP3P model for water29 were used. The particle mesh Ewald (PME) method30 was used to evaluate long-range electrostatic interactions. Non-bonded interactions were truncated at 10.0 Å with an atom-based switching function that was effective at 8.0 Å. The pairlist for non-bonded interactions was truncated at 11.5 Å and updated every 10 steps. All bonds involving hydrogen atoms in Bi9 were constrained using SHAKE31 and water molecules were treated as rigid with SETTLE32. Langevin thermostat33 with a damping coefficient of 5 ps-1 was used to simulate the NVT ensemble. Leapfrog integrator was used with an integration time-step of 2 fs. All rings in Bi9 were

ACS Paragon Plus Environment

5

Journal of Chemical Theory and Computation

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

Page 6 of 36

restrained to prevent inversions by adding the flat-bottom harmonic potentials to the 6 dihedral angles (O5-C1-C2-C3, C1-C2-C3-C4, C2-C3-C4-C5, C3-C4-C5-O5, C4-C5-O5-C1, and C5-O5-C1-C2):    Vrest (φ ) =    

1 2 k (φ − φlb ) 2 0 1 2 k (φ − φub ) 2

φ < φlb φlb ≤ φ ≤ φub , φub < φ

(1)

where φlb = 30°, φub = 90°, and k = 20.0 (kcal/mol)/deg2. The two initial Bi9 conformations, gg and gt, with explicit water molecules were minimized and equilibrated with the two MD simulations: 1 ns NPT dynamics at 300 K and 1 atm, and another 1 ns NVT dynamics at 300 K before the production runs. We use the notations of s(gg) and s(gt) to specify different simulations starting from gg or gt conformations. Additionally, two independent conventional MD simulations (1 µs), s(gg) and s(gt), were performed for comparison. Four Graphic Processing Units (GPUs) were used for each simulation.34 In order to further speed up the simulation, the multiple time step integrator (rRESPA, reference system propagator algorithm)35 was used with 2.5 fs time step for the shortrange interactions and 5 fs time step for the long-range interactions.

2.2. MTD, REMD, and RSE-MTD for enhancing conformational samplings. 2.2.1. Metadynamics (MTD) In MTD, history-dependent bias potential, Vbias, is added to the system potential V0. ,

(2)

The bias potential is a function of a few collective variables (CV),

, which describe motions

related with high-energy barriers on the potential energy surface. During the simulation, every time interval τ, Vbias is updated with a small Gaussian-shaped kernel:

ACS Paragon Plus Environment

6

Page 7 of 36

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

Journal of Chemical Theory and Computation

, where

(3)

is an instantaneous CV value at time t. w and σ represent the height and width of the

Gaussian kernels, respectively. The kernel height is rescaled with the well-tempered MTD framework.20, 36 We selected the ω angle in Bi9 as a CV, because it was identified as crucial for the conformational flexibility of Bi9.13 For constructing Vbias in MTD, we used the parameters in eq. (2) and (3): τ = 1.0 ps (500 steps), w = 1.0 kcal/mol, σ = 0.1 rad, ∆T = 2700 K. The system temperature, T, was kept at 300 K. We performed two-replicas MTD simulations: one with the bias and another unbiased for trajectory analysis. These two replicas were exchanged with the Metropolis-Hastings algorithm37-38.

2.2.2. Replica-exchange molecular dynamics (REMD) In REMD21, copies of the original target systems, replicas, are simulated with canonical MD simulations at different temperatures. The temperatures are exchanged with the MetropolisHastings algorithm37-38 between the neighboring replicas or among all the replicas by the SuwaTodo algorithm39 (this method is often called as the replica-permutation method40). Theoretically, the former is designed to impose the detailed balance conditions for the replica exchanges and the latter satisfies only the global balance condition. MD simulations of replicas at high temperatures allow us to sample wider conformational spaces of the target molecules by avoiding from trapping at local energy minimum states, while those at low temperatures search stable conformations near the minimum energy states. By exchanging temperatures, we can realize a random walk in the temperature space, which, in turn, induces a random walk in the potential energy space to overcome energy barriers in the systems with rugged energy surfaces.

ACS Paragon Plus Environment

7

Journal of Chemical Theory and Computation

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

In the current study, we used 64 replicas, which cover a temperature range of 300-497 K41. Replica exchanges were attempted every 2.0 ps (1000 steps) based on the Metropolis-Hasting algorithm with alternating even/odd pairs of replicas. The simulations show 41~55% of the acceptance ratios of replica exchanges (Figure S2). We found random walks in both replica (Figure S3) and temperature (Figure S4) spaces. The highest temperature might not be sufficient to overcome high-energy barriers along with the ω angles in Bi9, but the number of replicas in this study was determined due to the resource limitations.

2.2.3. Replica state exchange metadynamics (RSE-MTD) RSE-MTD22 is one of the generalized-ensemble algorithms, which was recently developed to combine the merits of MTD and REMD algorithms. Two similar methods42-43 were already proposed by developers of MTD. In one method, multiple MTDs are performed at different temperatures, and the temperatures and bias potentials are exchanged with the MetropolisHasting algorithm like temperature-REMD. In another method, MTD with the bias potentials for different CVs are carried out and the bias potentials are exchanged with the Metropolis-Hasting algorithm like Hamiltonian-REMD (or REUS). The former is called as PT-MTD42, while the latter is named as BE-MTD43. RSE-MTD generalizes these methods in the following two points. First, we allow the replicas with different bias potentials and at multiple temperatures. Second, we use the Suwa-Todo algorithm for exchanging all the replicas at once, by imposing the global balance condition. To emphasize the exchanges of replica states, we named it as the replica state exchange metadynamics (RSE-MTD). We tested the conformational sampling efficiency of this method with alanine polypeptides in gas and water.22 Note that similar approaches combining the different types of enhanced sampling methods have also been reported recently.44

ACS Paragon Plus Environment

8

Page 8 of 36

Page 9 of 36

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

Journal of Chemical Theory and Computation

The construction of bias potentials in RSE-MTD was almost the same as that described in the MTD section. We carried out simulations with two RSE-MTD schemes. One simulation (RSEMTD(PT)) includes 64 replicas that cover a temperature range of 300-493 K. One replica is simulated without any bias, while another 63 replicas have the Vbias for the ω angle in Bi9. The parameters and construction in Vbias are the same as those in MTD. In RSE-MTD(PT) scheme, the exchanges were performed with the partial Suwa-Todo algorithm, where 64 replicas are divided into 15 overlapping groups: A groups (1-8, 9-16, 17-24, 25-32, 33-40, 41-48, 49-56, and 57-64) and B groups (5-12, 13-20, 21-28, 29-36, 37-44, 45-52, and 53-60). The exchanges within the A and B groups were attempted in alternating order (i.e. A, B, A, B, …). In another scheme of RSE-MTD (RSE-MTD(BE)), one unbiased and another three replicas with different Vbias for each collective variables (ω, χ(23), and χ(3'4')) were simulated at 300 K. The exchanges of replicas were done using the Suwa-Todo algorithm. The replica conditions for all the methods used in this study are summarized in Table 1.

2.3. Conformational Analysis of the Simulation Trajectories The simulation length for each replica in MTD, REMD, RSE-MTD(PT), and RSE-MTD(BE) was 50 ns. We analyzed 20000 snapshots taken from the last 40 ns of the unbiased replica in each simulation. The trajectories are analyzed in terms of the time-series of collective variables (CVs), the conformational entropy45, the free-energy landscapes, and structure clustering. Regarding CVs, we focus on the ω angle which is enhanced with a bias potential in MTD as well as RSE-MTDs, while the χ(23), and χ(3'4') angles are selected as examples of unbiased CVs in MTD, REMD, and RSE-MTD(PT). The conformational entropy, S, is defined using the equations (4) and (5):

ACS Paragon Plus Environment

9

Journal of Chemical Theory and Computation

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

Page 10 of 36



Sk = − ∫ 0 Pk ( q) ln ( Pk ( q) ) dq

,

(4)

K

S = ∑ ( Sk − S0 ) k=1

,

(5)

where Pk(q), Sk, and K are the probability distribution function of the k-th degree of freedom, the conformational entropy originating from the k-th degree of freedom, and the number of freedom, respectively. We use all the flexible dihedral angles in Bi9 to computing the entropy. This quantity was previously used in protein folding simulations as an estimator of the simulation convergence.45 The free-energy landscapes of Bi9 are drawn using the spherical coordinates (η and θ angles), which were proposed by Nishima and coworkers.12 η represents the swing motion of the 1-6 branch, whereas θ shows the up-down motion with respect to the xy-plane defined with the core hexagon (see Figure 1 in Ref. [12]). They can separate major conformers of Bi9. In our structure clustering, we use the ω angle and the eight pseudo-dihedral angles defined in Section 2.1: χ(12), χ(23), χ(34), χ(45), χ(56), χ(34'), χ(4'5'), and χ(5'6') for investigating the linkage flexibility. The density-based spatial clustering of applications with noise (DBSCAN) algorithm46 (as implemented in Scikit-learn47) is used, where the minimum size of cluster and the spatial resolution (ε) are set to 1% of samples and 15º, respectively. We consider the periodicity of the dihedral and pseudo-dihedral angles in the structure clustering.

3. RESULTS AND DISCUSSION First, we examine the conformational sampling features of MTD, REMD and RSE-MTD(PT). After investigating the results of RSE-MTD(PT) simulations, three essential CVs are found and used in RSE-MTD(BE). Finally, we discuss the sampling efficiency of all the four simulation methods and compare the results with experimental data.

ACS Paragon Plus Environment

10

Page 11 of 36

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

Journal of Chemical Theory and Computation

3.1. Validation of MTD, REMD, and RSE-MTD First, we verify each of the enhanced conformational sampling methods used in this study. For MTD and RSE-MTD, the convergence of the bias potential is necessary for estimating freeenergy changes. In Figure 2, we plot the root mean square deviations (RMSDs) of the freeenergy changes along the ω angle with respect to the final values at 50 ns. Large deviations were observed in the first 10 ns of MTD s(gg) and s(gt), while RSE-MTD(PT) s(gg) and s(gt) show quick convergences toward the final value with small errors of less than 1 kcal/mol. After 10 ns, the RMSDs for RSE-MTD(PT) s(gg) and s(gt) are much smaller than those of MTD. Figure S2 shows the free-energy changes along the ω angle. The free energy is converged within ~1 kcal/mol in all simulation methods.

3.2. Time series of ω and other degrees of freedom Next, we study the conformational fluctuations of Bi9 in different simulations. Figure 3 shows the time series of the ω angle in the selected two replicas of the simulations (replicas 1 and 2 for MTD, REMD, and RSE-MTD(PT)). In all the cases, the transitions between the gg and gt conformations were observed frequently, implying good random walks in the space. These transitions rarely happen in the normal MD simulations at the room temperatures. Due to the bias potentials, MTD and RSE-MTD(PT) show transitions more often. In contrast, the transitions are less frequent in REMD. For 10-20 ns, the simulation stayed in one of the states (either ω ≈ 60º or -60º) and suddenly jumped into the other stable state. This slow transition happened due to the random walks in temperature space, ensuring all the moderate-energy barrier crossings.

ACS Paragon Plus Environment

11

Journal of Chemical Theory and Computation

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

Page 12 of 36

Since the frequent structural transitions along the CV (ω) in MTD and RSE-MTD(PT) were expected in theory, we examine the other coordinates, which were not biased in all the three simulations. We selected the pseudo-dihedral angles, χ(34') and χ(23), for this purpose and plot their time-series in Figure 4 (see also Figure S4). We used the same replicas for REMD and RSE-MTD(PT) as those in Figure 3. In MTD s(gg) and s(gt), those two angles did not show frequent transitions: χ(34') stayed at -120º in most of the time and rarely jumped into 120º. In MTD s(gg), no transition was observed from χ(23) = 45º, while only one transition event occurred in MTD s(gt) for the two replicas. In contrast, REMD and RSE-MTD(PT) show frequent transitions in both χ(34') and χ(23). The frequency of the transitions between REMD and RSE-MTD(PT) are similar with each other. In Figure 5, the conformational entropy45 shows that the MTD s(gg) and s(gt) converge slowly to different values, suggesting that the sampling is non-ergodic. In contrast, REMD and RSEMTD(PT) rapidly converge to the equilibrium values within the error of ~1 kBT just after the first 10 ns. The result suggests that REMD and RSE-MTD(PT) sample approximately the same number of micro-states. Those differences result from the methodology of each enhanced sampling method. In MTD and RSE-MTD(PT), simulations were enhanced to lower the critical energy barriers of local coordinates, ω (Figure 3), with the bias potentials. However, as shown in Figure 4, other degrees of freedoms were not sufficiently sampled just by the bias potentials in MTD. REMD and RSEMTD(PT) utilized the exchanges of temperatures among the replicas. The high temperature simulations in those methods make the fluctuations of global motions of Bi9 larger than those at low temperatures. Moderate-energy barriers along the pseudo-dihedral angles in Figure 4 could be overcome by the increase of fluctuations. The conformational entropy in Figure 5 shows that

ACS Paragon Plus Environment

12

Page 13 of 36

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

Journal of Chemical Theory and Computation

these two methods realize more ergodic samplings. In terms of both local and global conformational sampling features, RSE-MTD(PT) seems to be the best candidate for efficient simulations of N-glycans in solution.

3.3. Free-energy landscapes and structural clustering at room temperatures For an ergodic simulation, the free energy estimate has to be reproducible and independent from the initial conformation. The free energy plot along the spherical coordinates (η and θ angles) has been shown that it can show distributions of the important conformers of Bi9, namely, Extended A (Ex-a), Extended B (Ex-b), Tight backfold (TBF), and Backfold (BF).12 Their positions in the free-energy plot and the representative conformations are shown in Figure 6 (in the plot of RSE-MTD(PT) s(gg)) and Figure 7, respectively. In MTD, noticeable differences between s(gg) and s(gt) were observed and TBF structures were rarely sampled. In REMD, the sampled spaces are more consistent between s(gg) and s(gt) and all the four conformations were observed as independent basins. However, the populations of local minima are different, i.e. the s(gg) favoured the gg conformations (θ ≈ 40°), but the s(gt) sampled the gt conformations (θ ≈ -20°) more often. In other words, REMD s(gg) oversampled Ex-a, because of the initial conformation dependency. The RSE-MTD(PT) results look the most converged. The populations of Ex-a and Ex-b seem to be the same regardless of the starting conformations. We further examine the equilibrium solution structures of Bi9 in terms of 9 internal coordinates: the ω angle and 8 pseudo-dihedral angles. If we neglect the internal conformational changes in each saccharide unit, this should be a complete description of Bi9 structure through all the linkage connectivity. The trajectories of unbiased replicas from the RSE -MTD(PT) s(gg) and s(gt) were merged for this clustering analysis based on the DBSCAN algorithm. DBSCAN

ACS Paragon Plus Environment

13

Journal of Chemical Theory and Computation

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

Page 14 of 36

detected 5 structures in Figure 7. The first four conformers were already assigned in our previous study12, such as Ex-a, Ex-b, TBF, and BF, while the last one, Extended C (Ex-c), was found only in this study. In Figure 8, we project each conformer to the space spanning with the η and θ angles (same as the free energy map in Figure 6). We found that the distributions of Ex-b and Ex-c are overlapped in the space. Projection of the five representative clusters to the 9 internal coordinates reveal that just 3 interresidue angles, ω, χ(34') and χ(23), are necessary to distinguish all the clusters (Figure 8 and Figure S3). The ω angle makes a distinction between the gg (Ex-a and BF) and gt (Ex-b, TBF, and Ex-c) conformers. χ(34') distinguishes the extended (Ex-a, Ex-b, and Ex-c) with compact (TBF and BF) conformations, while χ(23) distinguishes Ex-a with Ex-c. For quantitative analysis, conformations from unbiased trajectory in each simulation are assigned to the clusters using the nearest-neighbour classifier. The estimated cluster populations in Table 2 show the same trends with the free energy plots (Figure 6), i.e. the RSE-MTD(PT) simulations are better converged and more consistent than MTD and REMD.

3.4. RSE-MTD(BE) simulations The key information from the clustering is that three inter-residue angles, ω, χ(34') and χ(23), are able to account for the N-glycan flexibility and all the representative conformers of Bi9 in solution. Assuming that the linkage flexibility is the only source of high-energy barriers, we propose RSE-MTD(BE), where the three angles are used as CVs. In addition to one unbiased replica, the four replicas are simulated at 300 K. Since the bias potentials could enhance sampling of the CVs, good random walks in the coordinates were observed in Figures 3 and 4. Furthermore, the conformational entropy obtained with RSE-MTD(BE) converges to the same

ACS Paragon Plus Environment

14

Page 15 of 36

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

Journal of Chemical Theory and Computation

values as REMD and RSE-MTD(PT), suggesting that these three methods sampled almost the same micro states. However, the convergences of free-energy map (Figure 6) and populations of the representative clusters with RSE-MTD(BE) are slightly worse compared to RSE-MTD(PT) simulation: the differences between s(gg) and s(gt) are larger than those of RSE-MTD(PT). Considering that the total number of replicas in RSE-MTD(BE) is significantly smaller than RSE-MTD(PT), this method can be a reasonable choice for simulating larger biological systems, such as singly or multiply glycosilated proteins in solution or in membranes.

3.5. Comparison between simulations and experiments The comparison with experimental data cannot be done without touching force field accuracy and experimental data interpretation issues. Nevertheless, it serves as an indicator how realistic are these simulations. First, we compare the cluster populations with those estimated from the NMR pseudo-contact shift (PCS) data.48 In general, there is a reasonably good agreement in terms of the Ex-a, Ex-b, and BF populations (Table 2), e.g. in case of RSE- MTD(PT), the estimate error are below 5%. However, we find the discrepancy in the HBF and TBF populations between the simulated and the experimentally derived results. RSE-MTD(PT) gives the gg:gt ratio of 52:48 on average, while the experiment suggests the ration of 65:35. Second, the simulation results are compared with the NMR scalar 3J-coupling constants (J56 and J56’)49, which are related with the ω angle by the modified Karplus equations.12, 50 The J56 values (Table 2) are not sensitive to the sampling protocols, i.e. there is just ~0.2 Hz variations and ~1 Hz systematic overestimation. Meanwhile, there is a larger variation of the J56' values in the different simulations. The RSE-PT-MTD simulations show the best agreement (~0.2 Hz error) with the experiment.

ACS Paragon Plus Environment

15

Journal of Chemical Theory and Computation

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

Page 16 of 36

3.6. Comparison with long time conventional MD simulation Finally, we compare the enhanced sampling approaches with conventional MD simulations (1 µs). The transitions between the gg and gt conformations occur every few hundreds ns (three times in s(gg) and six times in s(gt) simulations, Figure S8) due to the high-energy barrier separating the two states. This underscores the importance of the enhanced sampling approach to obtain the converged conformational sampling of N-glycans in solution.

4. CONCLUSION AND FUTURE PERSPECTIVES In the current study, we examine the conformational sampling efficiency of enhanced sampling schemes, MTD, REMD, RSE-MTD(PT), and RSE-MTD(BE). The simple MTD does not show good performance except for the sampling of the biased CV. Among the four methods, RSE-MTD(PT) shows the best sampling performance in terms of the frequency of structural transitions, the accuracy of free-energy landscapes, and the cluster populations. However, considering the computational cost, RSE-MTD(BE) could be a reasonable choice for larger Nglycan-attached proteins. The penta-saccharide part, which is crucial for overall conformation of Bi9, is shared with the most of N-glycans. Although we simulated only Bi9, RSE-MTD(BE) methods, applying the bias potential to the three inter-residue angles: ω, χ(34') and χ(23), would be generally useful for all N-glycans. In other systems, including proteins, nucleic acids, or other biomolecules, if we could find essential CVs a proiri, or after short REMD simulations, it is possible to reduce the computational cost in the same way as we propose with RSE-MTD(BE).

ACS Paragon Plus Environment

16

Page 17 of 36

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

Journal of Chemical Theory and Computation

Figure 1. (a) The topological structure of a N-glycan, Bi9, with the branch names. (b) The central part of Bi9 with the φ, ψ, and ω angles, and the three pseudo-dihedral angles, χ (23), χ(34), and χ(34').

ACS Paragon Plus Environment

17

Journal of Chemical Theory and Computation

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

Page 18 of 36

Figure 2. Root mean square deviations of the free energy along the ω angle with respect to the final values at 50 ns in the MTD s(gg), MTD s(gt), RSE-MTD(PT) s(gg), and RSE-MTD(PT) s(gt) simulations, respectively.

ACS Paragon Plus Environment

18

Page 19 of 36

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

Journal of Chemical Theory and Computation

Figure 3. Time series of the ω angle in MTD, REMD, RSE-MTD(PT), and RSE-MTD(BE) simulations. s(gg) and s(gt) are different simulations starting from gg and gt conformations of Bi9, respectively.

ACS Paragon Plus Environment

19

Journal of Chemical Theory and Computation

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

Page 20 of 36

Figure 4. Time series of the two pseudo-dihedral angles, χ(34') and χ(23) in MTD, REMD, RSEMTD(PT), and RSE-MTD(BE) simulations. s(gg) and s(gt) are different simulations starting from gg and gt conformations of Bi9, respectively.

ACS Paragon Plus Environment

20

Page 21 of 36

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

Journal of Chemical Theory and Computation

Figure 5. Conformational entropy from the MTD, REMD, RSE-MTD(PT), and RSE-MTD(BE) simulations.

ACS Paragon Plus Environment

21

Journal of Chemical Theory and Computation

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

Page 22 of 36

Figure 6. The free-energy plots in the spherical coordinates by Nishima et al. (see text for details). Each plot is obtained from the trajectory of MTD, REMD, RSE-MTD(PT), and RSEMTD(BE). For each method, two independent simulations with different initial conformations were carried out as s(gg) and s(gt). The positions of major conformers found by the clustering analysis were indicated in the plot of RSE-MTD(PT) s(gg).

ACS Paragon Plus Environment

22

Page 23 of 36

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

Journal of Chemical Theory and Computation

Figure 7. The five representative clusters are shown with top ten representative structures in each cluster. The core, 1-3 branch, and 1-6 branch are colored in blue, green, and red, respectively. Both top and front views are shown.

ACS Paragon Plus Environment

23

Journal of Chemical Theory and Computation

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

Page 24 of 36

Figure 8. The projections of the five clusters to (a) (ω, χ(23)) (b) (ω, χ(34')) and (c) (η, θ) spaces. Red, green, blue, yellow, and magenta dots represent Ex-a, Ex-b, TBF, BF, and Ex-c conformers, respectively.

ACS Paragon Plus Environment

24

Page 25 of 36

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

Journal of Chemical Theory and Computation

Figure 9. The ratio of gg:gt conformations estimated from MTD, REMD, RSE-MTD(PT), and RSE-MTD(BE). The corresponding ratio derived from NMR PCS data (Ref. [48]) is also shown for comparison.

ACS Paragon Plus Environment

25

Journal of Chemical Theory and Computation

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

Page 26 of 36

Table 1. Replica conditions (temperature and bias) for the MTD, REMD, RSE-MTD(PT), and RSE -MTD(BE) methods.

Replica condition

MTD

REMD

RSE-MTD(PT)

RSE-MTD(BE)

T [K]

CV*

T [K]

CV*

T [K]

CV*

T [K]

CV*

1

300.00

none

300.00

none

300.00

none

300.00

none

2

300.00

ω

302.49

none

300.00

ω

300.00

ω

3

305.00

none

302.49

ω

300.00

χ(34')

4

307.52

none

305.00

ω

300.00

χ(23)











63

493.04

none

489.20

ω

64

496.90

none

493.04

ω

* CV “none” refers to unbiased replica.

ACS Paragon Plus Environment

26

Page 27 of 36

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

Journal of Chemical Theory and Computation

Table 2. Experimental and simulated cluster populations† and scalar 3J-coupling constants. Each row corresponds to an independent simulation with the specified method (MTD, REMD, RSEMTD(PT), or RSE-MTD(BE)) and initial condition (s(gg) or s(gt)). 3

gg conformers [%]

Method/ Simulation

gt conformers [%]

J-coupling [Hz]

Ex-a

BF

HBF

Ex-b

TBF

Ex-c

J56

J56'

45*

10*

10*

35*

0*



2.1#

5.8#

s(gg)

55.9

4.4

0.0§

32.9

5.9

0.0

2.9

4.8

s(gt)

61.9

3.9

0.0§

30.0

1.7

2.1

2.8

4.4

s(gg)

56.0

8.5

0.0§

23.0

10.8

0.6

2.9

4.6

s(gt)

37.9

4.6

0.0§

36.9

14.8

4.9

3.0

6.5

s(gg)

42.7

6.9

0.0§

30.3

16.9

2.2

3.0

5.9

s(gt)

46.9

6.7

0.0§

33.9

10.6

1.3

2.9

5.6

s(gg)

43.1

10.5

0.0§

32.1

10.5

3.7

2.9

5.6

s(gt)

37.8

7.8

0.0§

44.2

10.1

0.1

2.9

6.3

Experimental data

MTD

REMD

RSE-MTD(PT)

RSE-MTD(BE) †

DBSCAN Ref. [46-47] classifies some samples as noise, so the cluster populations do not add up to 100%. *

The experimental cluster populations obtained from the PCS data fitting Ref. [48].

#

The experimental NMR scalar 3J-coupling constants taken from Ref. [49].



Ex-c was not taken into account in the PCS data fitting Ref. [48].

§

In our simulations, HBF was not found as a separate cluster, but rather as a part of BF.

ACS Paragon Plus Environment

27

Journal of Chemical Theory and Computation

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

ASSOCIATED CONTENT Supporting Information. The supporting information is available free of charge on the ACS publications website. Four supplementary figures (S1–S4; PDF) AUTHOR INFORMATION Corresponding Author *E-mail: [email protected] (Y. S.). Author Contributions The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. Notes The authors declare no competing financial interest. ACKNOWLEDGMENT This research was supported in part by the RIKEN basic science interdisciplinary research projects (dynamic structural biology, molecular systems research, iTHES, and integrated lipidology research), Innovative Drug Discovery Infrastructure through Functional Control of Biomolecular Systems, Priority Issue 1 in Post-K Supercomputer Development (Project ID:hp150270), MEXT Grant-in-Aid for Scientific Research on Innovative Areas Grant Number 26119006 (to Y.S.), and MEXT/JSPS KAKENHI Grant Number 26220807 (to Y.S.), 16K00415 (to S. R.), Center of innovation Program from Japan Science and Technology Agency, JST (to Y.

ACS Paragon Plus Environment

28

Page 28 of 36

Page 29 of 36

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

Journal of Chemical Theory and Computation

S.). This research used computational resources of the HPCI system provided by the University of Tokyo through the HPCI System Research Project (Project ID: hp140157), the RIKEN Integrated Cluster of Clusters (RICC), and MEXT SPIRE Supercomputational Life Science (SCLS).

REFERENCES (1)

Werz, D. B.; Ranzinger, R.; Herget, S.; Adibekian, A.; von der Lieth, C. W.; Seeberger, P. H., Exploring the structural diversity of mammalian carbohydrates ("glycospace") by statistical databank analysis. ACS Chem. Biol. 2007, 2 (10), 685-91.

(2)

Gabius, H. J., The magic of the sugar code. Trends Biochem. Sci. 2015, 40 (7), 341.

(3)

Gabius, H. J.; Andre, S.; Jimenez-Barbero, J.; Romero, A.; Solis, D., From lectin structure to functional glycomics: principles of the sugar code. Trends Biochem. Sci. 2011, 36 (6), 298-313.

(4)

Hardy, B. J., The glycosidic linkage flexibility and time-scale similarity hypotheses. J. Mol. Struct. -Theochem 1997, 395, 187-200.

(5)

Weiss, A. A.; Iyer, S. S., Glycomics aims to interpret the third molecular language of cells. Microbe 2007, 2, 489-497.

(6)

Wormald, M. R.; Petrescu, A. J.; Pao, Y. L.; Glithero, A.; Elliott, T.; Dwek, R. A., Conformational studies of oligosaccharides and glycopeptides: complementarity of NMR, X-ray crystallography, and molecular modelling. Chem. Rev. 2002, 102 (2), 371-86.

ACS Paragon Plus Environment

29

Journal of Chemical Theory and Computation

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

(7)

Page 30 of 36

Peric-Hassler, L.; Hansen, H. S.; Baron, R.; Hunenberger, P. H., Conformational properties of glucose-based disaccharides investigated using molecular dynamics simulations with local elevation umbrella sampling. Carbohydr. Res. 2010, 345 (12), 1781-801.

(8)

Karplus, M.; McCammon, J. A., Molecular dynamics simulations of biomolecules (vol 9, pg 646, 2002). Nat. Struct. Biol. 2002, 9 (10), 788-788.

(9)

Dror, R. O.; Dirks, R. M.; Grossman, J. P.; Xu, H.; Shaw, D. E., Biomolecular simulation: a computational microscope for molecular biology. Annu. Rev. Biophys. 2012, 41, 429-52.

(10) Bernardi, R. C.; Melo, M. C.; Schulten, K., Enhanced sampling techniques in molecular dynamics simulations of biological systems. Biochim. Biophys. Acta 2015, 1850 (5), 872-7. (11) Mitsutake, A.; Sugita, Y.; Okamoto, Y., Generalized-ensemble algorithms for molecular simulations of biopolymers. Biopolymers 2001, 60 (2), 96-123. (12) Nishima, W.; Miyashita, N.; Yamaguchi, Y.; Sugita, Y.; Re, S., Effect of bisecting GlcNAc and core fucosylation on conformational properties of biantennary complex-type N-glycans in solution. J. Phys. Chem. B 2012, 116 (29), 8504-12. (13) Re, S.; Miyashita, N.; Yamaguchi, Y.; Sugita, Y., Structural diversity and changes in conformational equilibria of biantennary complex-type N-glycans in water revealed by replica-exchange molecular dynamics simulation. Biophys. J. 2011, 101 (10), L44-6. (14) Yamaguchi, Y.; Nishima, W.; Re, S.; Sugita, Y., Confident identification of isomeric Nglycan structures by combined ion mobility mass spectrometry and hydrophilic interaction liquid chromatography. Rapid Commun. Mass Spectrom. 2012, 26 (24), 2877-84.

ACS Paragon Plus Environment

30

Page 31 of 36

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

Journal of Chemical Theory and Computation

(15) Fukunishi, H.; Watanabe, O.; Takada, S., On the Hamiltonian replica exchange method for efficient sampling of biomolecular systems: Application to protein structure prediction. J. Chem. Phys. 2002, 116 (20), 9058-9067. (16) Sugita, Y.; Kitao, A.; Okamoto, Y., Multidimensional replica-exchange method for freeenergy calculations. J. Chem. Phys. 2000, 113 (15), 6042-6051. (17) Patel, D. S.; Pendrill, R.; Mallajosyula, S. S.; Widmalm, G.; MacKerell, A. D., Jr, Conformational properties of α- or β-(1→6)-linked oligosaccharides: Hamiltonian replica exchange MD simulations and NMR experiments. J. Phys. Chem. B 2014, 118 (11), 285171. (18) Yang, M.; MacKerell, A. D., Jr, Conformational sampling of oligosaccharides using Hamiltonian replica exchange with two-dimensional dihedral biasing potentials and the weighted histogram analysis method (WHAM). J. Chem. Theory Comput. 2015, 11 (2), 788-99. (19) Laio, A.; Parrinello, M., Escaping free-energy minima. Proc. Natl. Acad. Sci. U. S. A. 2002, 99 (20), 12562-6. (20) Barducci, A.; Bussi, G.; Parrinello, M., Well-tempered metadynamics: a smoothly converging and tunable free-energy method. Phys. Rev. Lett. 2008, 100 (2), 020603. (21) Sugita, Y.; Okamoto, Y., Replica-exchange molecular dynamics method for protein folding. Chem. Phys. Lett. 1999, 314 (1-2), 141-151. (22) Galvelis, R.; Sugita, Y., Replica state exchange metadynamics for improving the convergence of free energy estimates. J. Comput. Chem. 2015, 36 (19), 1446-55.

ACS Paragon Plus Environment

31

Journal of Chemical Theory and Computation

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

Page 32 of 36

(23) Woods, R. J., GLYCAM Web. Complex Carbohydrate Research Center, University of Georgia: Athens, GA: 2005-2011. (24) Humphrey, W.; Dalke, A.; Schulten, K., VMD: visual molecular dynamics. J. Mol. Graph. 1996, 14 (1), 33-8, 27-8. (25) Jung, J.; Mori, T.; Kobayashi, C.; Matsunaga, Y.; Yoda, T.; Feig, M.; Sugita, Y., GENESIS: a hybrid-parallel and multi-scale molecular dynamics simulator with enhanced sampling algorithms for biomolecular and cellular simulations. Wiley Interdiscip Rev. Comput. Mol. Sci. 2015, 5 (4), 310-323. (26) Guvench, O.; Greene, S. N.; Kamath, G.; Brady, J. W.; Venable, R. M.; Pastor, R. W.; MacKerell, A. D., Jr, Additive empirical force field for hexopyranose monosaccharides. J. Comput. Chem. 2008, 29 (15), 2543-64. (27) Guvench, O.; Hatcher, E. R.; Venable, R. M.; Pastor, R. W.; MacKerell, A. D., Jr, CHARMM Additive All-Atom Force Field for Glycosidic Linkages between Hexopyranoses. J. Chem. Theory Comput. 2009, 5 (9), 2353-2370. (28) Guvench, O.; Mallajosyula, S. S.; Raman, E. P.; Hatcher, E.; Vanommeslaeghe, K.; Foster, T. J.; Jamison, F. W., II; MacKerell, A. D., Jr, CHARMM additive all-atom force field for carbohydrate derivatives and its utility in polysaccharide and carbohydrate-protein modeling. J. Chem. Theory Comput. 2011, 7 (10), 3162-3180. (29) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L., Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79 (2), 926-935.

ACS Paragon Plus Environment

32

Page 33 of 36

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

Journal of Chemical Theory and Computation

(30) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G., A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103 (19), 8577-8593. (31) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C., Numerical-Integration of Cartesian Equations of Motion of a System with Constraints - Molecular-Dynamics of N-Alkanes. J. Comput. Phys. 1977, 23 (3), 327-341. (32) Miyamoto, S.; Kollman, P. A., Settle - an Analytical Version of the Shake and Rattle Algorithm for Rigid Water Models. J. Comput. Chem. 1992, 13 (8), 952-962. (33) Adelman, S. A.; Doll, J. D., Generalized Langevin Equation Approach for Atom-SolidSurface Scattering - General Formulation for Classical Scattering Off Harmonic Solids. J. Chem. Phys. 1976, 64 (6), 2375-2388. (34) Jung, J.; Naurse, A.; Kobayashi, C.; Sugita, Y., Graphics Processing Unit Acceleration and Parallelization of GENESIS for Large-Scale Molecular Dynamics Simulations. J. Chem. Theory Comput. 2016, 12 (10), 4947-4958. (35) Tuckerman, M.; Berne, B. J.; Martyna, G. J., Reversible Multiple Time Scale MolecularDynamics. J. Chem. Phys. 1992, 97 (3), 1990-2001. (36) Dama, J. F.; Parrinello, M.; Voth, G. A., Well-tempered metadynamics converges asymptotically. Phys. Rev. Lett. 2014, 112 (24), 240602. (37) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E., Equation of State Calculations by Fast Computing Machines. J. Chem. Phys. 1953, 21 (6), 10871092.

ACS Paragon Plus Environment

33

Journal of Chemical Theory and Computation

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

Page 34 of 36

(38) Hastings, W. K., Monte-Carlo Sampling Methods Using Markov Chains and Their Applications. Biometrika 1970, 57 (1), 97-109. (39) Suwa, H.; Todo, S., Markov chain Monte Carlo method without detailed balance. Phys. Rev. Lett. 2010, 105 (12), 120603. (40) Itoh, S. G.; Okumura, H., Replica-Permutation Method with the Suwa-Todo Algorithm beyond the Replica-Exchange Method. J. Chem. Theory Comput. 2013, 9 (1), 570-81. (41) Patriksson, A.; van der Spoel, D., A temperature predictor for parallel tempering simulations. Phys. Chem. Chem. Phys. 2008, 10 (15), 2073-7. (42) Bussi, G.; Gervasio, F. L.; Laio, A.; Parrinello, M., Free-energy landscape for beta hairpin folding from combined parallel tempering and metadynamics. J. Am. Chem. Soc. 2006, 128 (41), 13435-41. (43) Piana, S.; Laio, A., A bias-exchange approach to protein folding. J. Phys. Chem. B 2007, 111 (17), 4553-9. (44) Xie, L.; Shen, L.; Chen, Z. N.; Yang, M., Efficient free energy calculations by combining two complementary tempering sampling methods. J. Chem. Phys. 2017, 146 (2), 024103. (45) Li, D. W.; Bruschweiler, R., In silico relationship between configurational entropy and soft degrees of freedom in proteins and peptides. Phys. Rev. Lett. 2009, 102 (11), 118108. (46) Ester, M.; Kriegel, H.-p.; S, J.; Xu, X. In A density-based algorithm for discovering clusters in large spatial databases with noise, 1996; AAAI Press: 1996; pp 226-231.

ACS Paragon Plus Environment

34

Page 35 of 36

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

Journal of Chemical Theory and Computation

(47) Pedregosa, F.; Varoquaux, G.; Gramfort, A.; Michel, V.; Thirion, B.; Grisel, O.; Blondel, M.; Prettenhofer, P.; Weiss, R.; Dubourg, V.; Vanderplas, J.; Passos, A.; Cournapeau, D.; Brucher, M.; Perrot, M.; Duchesnay, E., Scikit-learn: Machine Learning in Python. J. Mach. Learn. Res. 2011, 12, 2825-2830. (48) Canales, A.; Mallagaray, A.; Perez-Castells, J.; Boos, I.; Unverzagt, C.; Andre, S.; Gabius, H. J.; Canada, F. J.; Jimenez-Barbero, J., Breaking Pseudo-Symmetry in Multiantennary Complex N-Glycans Using Lanthanide-Binding Tags and NMR Pseudo-Contact Shifts. Angew. Chem. Int. Edit. 2013, 52 (51), 13789-13793. (49) Homans, S. W.; Dwek, R. A.; Boyd, J.; Mahmoudian, M.; Richards, W. G.; Rademacher, T. W., Conformational transitions in N-linked oligosaccharides. Biochemistry 1986, 25 (20), 6342-50. (50) 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 (19), 2783-2792.

ACS Paragon Plus Environment

35

Journal of Chemical Theory and Computation

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

Page 36 of 36

For Table of Contents Only

Enhanced Conformational Sampling of N-glycans in Solution with Replica State Exchange Metadynamics, Raimondas Galvelis, Suyong Re, and Yuji Sugita

ACS Paragon Plus Environment

36