Computing Nonequilibrium Conformational Dynamics of Structured

Nov 30, 2015 - ... of Biological Engineering, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, United States. ¶ Department of P...
0 downloads 0 Views 4MB Size
Subscriber access provided by KUNGL TEKNISKA HOGSKOLAN

Article

Computing Nonequilibrium Conformational Dynamics of Structured Nucleic Acid Assemblies Reza Sharifi Sedeh, Keyao Pan, Matthew Ralph Adendorff, Oskar Hallatschek, Klaus-Jürgen Bathe, and Mark Bathe J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.5b00965 • Publication Date (Web): 30 Nov 2015 Downloaded from http://pubs.acs.org on December 6, 2015

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 49

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

Computing Nonequilibrium Conformational Dynamics of Structured Nucleic Acid Assemblies Reza Sharifi Sedeh,†,# Keyao Pan,‡,# Matthew Ralph Adendorff,‡ Oskar Hallatschek,¶ Klaus-Jürgen Bathe,*,† and Mark Bathe*,‡



Department of Mechanical Engineering, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, United States



Department of Biological Engineering, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, United States ¶

Department of Physics, University of California, Berkeley, Berkeley, California 94720, United States

Abstract. Synthetic nucleic acids can be programmed to form precise three-dimensional structures on the nanometer-scale. These thermodynamically stable complexes can serve as structural scaffolds to spatially organize functional molecules including multiple enzymes, chromophores, and force-sensing elements with internal dynamics that include substrate reaction-diffusion, excitonic energy transfer, and force-displacement response that often depend critically on both the local and global conformational dynamics of the nucleic acid assembly. However, high molecular weight assemblies exhibit long time-scale and large length-scale

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

Page 2 of 49

motions that cannot easily be sampled using all-atom computational procedures such as molecular dynamics. As an alternative, here we present a computational framework to compute the over-damped conformational dynamics of structured nucleic acid assemblies and apply it to a DNA-based tweezer, a nine-layer ring origami object, and a pointer-shaped DNA origami object, which consist of 204, 3,600, and over seven thousand basepairs, respectively. The framework employs a mechanical finite element model for the DNA nanostructure combined with an implicit solvent model to either simulate the Brownian dynamics of the assembly or alternatively compute its Brownian modes. Computational results are compared with an all-atom molecular dynamics simulation of the DNA-based tweezer. Several hundred microseconds of Brownian dynamics are simulated for the nine-layer ring origami object to reveal its long time-scale conformational dynamics, and the first ten Brownian modes of the pointer-shaped structure are predicted.

Keywords. DNA Origami; Structural Nucleic Acid Assemblies; Brownian Dynamics; CoarseGrained Modeling

1

Introduction

Programmed self-assembly of synthetic nucleic acids enables the design of high molecular weight assemblies with complex three-dimensional architectures for diverse applications in biomolecular science and nanotechnology.1-7 For example, a nanometer-scale tweezer with two distinct conformational sub-states that are alternatively stabilized by addition of an extrinsic fuel strand was recently used to probe enzyme reaction kinetics at the single-molecule level.4 In the closed-state, the two arms of the tweezer are brought into close proximity of several nanometers, with approximately six-fold enhancement of the reaction product over the open-state.

ACS Paragon Plus Environment

2

Page 3 of 49

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

Alternatively, the principle of scaffolded DNA origami has been used in diverse contexts to program stable architectures of 1D, 2D, and 3D DNA objects for light-harvesting and nanoscale excitonic energy transfer.8, 9 RNA has similarly been assembled into large-scale structures using, for example, the principle of architectonics10 to form RNAi delivery vehicles or origami to cotranscriptionally fold high molecular weight assemblies in vitro.11 In the preceding applications of DNA origami, long time-scale conformational dynamics of the DNA assembly are coupled with short time-scale reaction-diffusion or excitonic energy transfer dynamics that can be simulated using a variety of stochastic and deterministic particle-based or continuum simulation procedures. Because the long time-scale conformational dynamics of DNA nanostructures are highly intensive to simulate computationally, there is a need for simulation procedures that enable the generation of conformational trajectories that can subsequently be coupled to complementary physics-based models of, for example, reaction-diffusion12, 13 or excitonic energy transfer dynamics.9, 14 Structured DNA assemblies that serve as three-dimensional scaffolds are typically programmed to adopt a single, dominant ground-state conformation that is fully determined by their Watson-Crick basepairing. This programmed basepairing can then be used to predict the target, three-dimensional structure of the assembly using coarse-grained finite element (FE) modeling.5, 15 Thermal fluctuations of these high molecular weight nanostructures, which span from several hundreds to many thousands of basepairs (bps), in a viscous solvent environment typically result in large-scale conformational dynamics that may range from tens to hundreds of nanoseconds (ns). While all-atom molecular dynamics (MD) is powerful in its versatile ability to investigate the mechanical, thermodynamic, and dynamical properties of DNA motifs,16-18 its first

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

Page 4 of 49

application to DNA origami structures was in 2013 by Yoo et al.19 and required a high performance supercomputer to facilitate its all-atom simulation consisting of 3,000–6,000 nucleotides on the 30–100 ns time-scale. While highly valuable to reveal both local and global conformational dynamics, such all-atom simulations of high molecular weight nanostructures is accessible to only a limited computational community, and alternative coarse-grained procedures are therefore needed to estimate large length-scale and long time-scale conformational dynamics more broadly and efficiently. While vacuum normal mode analysis (NMA) is useful to simulate the equilibrium conformational properties of high molecular weight macromolecular assemblies20-22 including DNA5, 15 and RNA-protein assemblies,23 it does not account for the effects of solvent on the structural mode shapes and their physical time-scales, rendering it inapplicable to the coarse-grained simulation of time-dependent conformational properties. Similar to models employed for NMA, a harmonic internal elastic energy is also assumed in recently proposed rigid base and bp models that, respectively, represent DNA bases and bps as rigid bodies with relative rotations and displacements of the bases and bps defining the DNA conformation.24, 25 While these models have proven successful in describing DNA structure and flexibility on the intermediate scale of 10–100 bps, they ignore effects of solvent damping on DNA and therefore do not compute physical time-scales of its motions. For this reason, a number of coarse-grained DNA models that include solvent effects have been developed.26-31 In these models, a reduced representation of DNA is typically employed instead of a conventional all-atom representation and force-field such as CHARMM3234

or AMBER35, 36, thereby dramatically reducing computational cost by grouping atoms into

pseudo-atoms37 or beads.38 Each pseudo-atom or bead may represent any number of atoms based on the level of coarse-graining desired, typically ranging from several pseudo-atoms per

ACS Paragon Plus Environment

4

Page 5 of 49

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

nucleotide37 to one bead per tens of bps.38 DNA hybridization and strand displacement kinetics can also be simulated using coarse-grained techniques.39, 40 Recently, He et al.41 also developed the physics-based coarse-grained DNA model NARES-2P and applied it to find that mean-field interactions between nucleic-acid-base dipoles are the main driving force for the formation of double-stranded DNA. In general, the degree of coarse-graining is typically kept homogeneous throughout the entire model so that fast internal degrees of freedom (DOF) are uniformly eliminated.37, 42 For example, in a one-bead-per-nucleotide model42 only inter-nucleotide interactions are considered and all atomic interactions within each nucleotide are ignored. These models with homogeneous degree of coarse-graining therefore retain a similar level of detail throughout the model. In contrast, Pantano and coauthors43 recently developed a hybrid model that employs varying levels of coarse-graining, from six pseudo-atoms per nucleotide37 to all-atom, for the multi-scale simulation of DNA. This model retains atomic detail in important regions of double-stranded DNA while treating the remainder of the system using six pseudo-atoms per nucleotide. For the efficient simulation of the coarse-grained dynamics of high molecular weight DNA assemblies such as DNA origami, oxDNA244, a recently extended version of oxDNA45, was developed. In this model, DNA is treated as a string of rigid nucleotides with one interaction site for the backbone and three per base, where finite-extensible nonlinear elastic springs are used to connect backbone sites. A “top-down” approach is used to derive the energy function and parameters such that the model captures accurately generic properties of DNA that are essential for modeling self-assembly.26,

44-48

This model has proven successful in reproducing duplex

hybridization,40 single-stranded stacking, and hairpin formation, as well as in capturing structural

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 49

and mechanical aspects of double-stranded DNA,49 and has been applied to DNA nanostructures including nanotweezers45 and DNA walkers48. Brownian dynamics (BD) is an alternative, well established stochastic model that uses implicit solvent modeling with white thermal noise to simulate dynamics due to fluctuationdissipation that are captured explicitly by classical MD. A landmark paper by Ermak and McCammon50 developed a generalized procedure to simulate the BD of N particles accounting

for hydrodynamic interactions, which were described by a 3N × 3N diffusion tensor. To compute

random displacements at each step in the Ermak-McCammon algorithm, the diffusion tensor is Cholesky-factorized, which results in a computational time that scales with O(N 3 ).51 Thus far, several approaches have been developed to reduce this computational cost52-55 in order to make application of the algorithm to large biomacromolecules for long time-scales feasible. For example, Fixman52 proposed the use of Chebyshev polynomials for the approximation of the square-root of the diffusion tensor, resulting in a computation-time scaling of O(N 2.25 ).56 Geyer and Winter presented yet another approach, termed the Truncated Expansion Ansatz,53 which approximates the Brownian displacement vector by truncating the expansion of the hydrodynamic multiparticle correlations as two-body contributions at the second-order and scales as O(N 2 ). Finally, Ando et al.54 recently proposed a new algorithm based on Krylov subspaces to approximate the random Brownian displacements, where computation-time also scales as O(N 2 ). Instead of approximating the square-root of the diffusion tensor, one may also keep the tensor unchanged for several sequential time-steps in order to speed up BD simulations57-59. In the past three decades, numerous BD models have been developed based on the original work of Ermak and McCammon50 to analyze the dynamics of linear and circular DNA.38, 57, 60-64 One of the first successful general-purpose BD packages for the simulation of

ACS Paragon Plus Environment

6

Page 7 of 49

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

linear and circular DNA was developed by Klenin and coauthors57. In their program, DNA was modeled as a chain of rigid segments connected by bending, torsion, and stretching potentials, and the Rotne-Prager tensor65 is used to model hydrodynamic interactions.57 This model was further adapted by others to study the dynamical behavior of DNA.58, 66-69 Recently, the scalable BD package BD_BOX was developed by Trylska and coauthors70 for the multi-scale simulation of systems comprised of various interacting macromolecules using flexible, coarse-grained bead models, where the package has been parallelized using MPI (message passing interface) and OpenMP (open multiprocessing) libraries to speed up the BD simulation of multicomponent systems.70-72 Also, Marenduzzo and coauthors73 recently developed a coarse-grained BD model for twistable elastic polymers using continuum rod theory based on a previous work of Chirico et al.38, and the model was implemented in the highly-scalable multi-purpose MD software package LAMMPS.74 Most of these BD models were not developed for the simulation of the overdamped conformational dynamics of highly structured DNA assemblies such as scaffolded DNA origami. Thus, the development of a coarse-grained, numerically stable, and efficient BD model designed specifically for simulations of structured DNA assemblies is needed. Here, we present such a specialized BD procedure for this purpose. The approach employs an underlying FE representation for the nucleic acid 3D structure,5, 15 in which DNA duplexes and crossover motifs are treated as elastic beams or wormlike chains with empirically prescribed ground-state geometric and mechanical properties. Electrostatic and van der Waals interactions or excluded volume between DNA duplexes are ignored together with any solventmediated effects that impact these interactions, consistent with previous assumptions.5,

15

Because our model captures the highly nonlinear and complex coupling between stretching, bending, and twisting DOF in megaDalton-scale DNA assemblies, it is capable of predicting

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

Page 8 of 49

complex 3D structure and flexibility. Over-damped dynamics of structured DNA assemblies are simulated using a bead-based model to capture the effects of viscous solvent drag on thermal fluctuations, which have previously been assumed to be at equilibrium and in vacuum. This is similar in spirit to previous BD simulation procedures for DNA that have been applied to long linear strands in solution, yet adapting these approaches to our FE based structural modeling approach.57,

75

We apply our procedure to compute the conformational dynamics of a DNA

nanotweezer, a nine-layer ring origami object, and a pointer-shaped origami object and compare the results with all-atom MD simulation of the tweezer and vacuum NMA. While the procedure is applied here only to three DNA-based assemblies, it is generally applicable to the simulation of long time-scale conformational dynamics of highly structured, thermodynamically stable nucleic acid assemblies that persist in a single dominant ground-state conformation. The computational procedure is available to the broader community as MATLAB code available as Supporting Information and online at http://cando-dna-origami.org.

2

Theory

FE modeling has previously been shown to be a useful coarse-grained approach to predict the 3D equilibrium structure and flexibility of programmed DNA assemblies.5, 15 Briefly, B-form DNA is treated as an elastic beam with stretching, bending, and twisting stiffnesses that are known empirically. Crossover strands that constrain neighboring duplexes at multi-way junctions have ground-state geometries and flexibilities that are assumed to adopt unique, ground-state values. Sequence effects on duplex and crossover mechanical properties, as well as non-bonded electrostatic and van der Waals interactions are ignored. The latter modeling assumptions are clearly only valid when the overall 3D structure is determined by the numerous crossovers that

ACS Paragon Plus Environment

8

Page 9 of 49

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

are used to constrain DNA duplexes to reside in well-defined 3D organizations, as is the case for a large set of programmed DNA (and RNA) assemblies.2,

5, 15

In cases where overall

nanostructure conformation is determined by complex non-bonded and excluded volume interactions, alternative coarse-grained modeling approaches should be used.44, 45 The FE model of the DNA nanostructure then contains N nodes whose ground-state structure can be computed using nonlinear FE analysis, where each “node” in the FE model corresponds to a single bp and is coincident with the bp reference point given by the software 3DNA. In addition, each FE node has three right-handed orthonormal reference axes that are collinear with the three reference directions of the bp given by the software 3DNA.76 Each FE node then models the three translational and three rotational DOF of a single bp, which allows conversion between the FE model of a DNA assembly and the corresponding all-atom model. The configuration of the FE model is then fully represented by a 6N-dimensional column vector x that contains these DOF. Because our BD framework is linear, infinitesimal displacements about the mechanical equilibrium, ground-state structure are assumed. Because the rotational DOF about the global coordinate axes are assumed to be infinitesimal and are therefore fully decoupled, the order in which these rotations are applied to each FE node is immaterial.51 A bead-based model for solvent drag is incorporated into our FE model to simulate effects of solvent-damping on thermally induced fluctuations that impact the long time-scale conformational dynamics of DNA nanostructures. In order to model the conformational dynamics of DNA nanostructures using the coarse-grained FE approach, discretized mass, M, stiffness, K, and drag, Z, matrices are needed together with governing equations of motion, which we describe next.

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

2.1

Page 10 of 49

FE Stiffness, Mass and Friction Matrices

The 6N × 6N diagonal mass matrix M and the 6N × 6N stiffness matrix K of the FE model of the DNA nanostructure are obtained directly from the ground-state structure that is computed using our previously introduced modeling approach (Figure 1).5, 15 In this model, DNA duplexes are assumed to adopt a canonical B-form conformation with axial rise of 0.34 nm per bp and righthanded helicity of 10.5 bp per turn in their ground-state conformation, although alternative DNA conformations such as A-form DNA can also be treated using the same FE framework. Each Bform DNA duplex is modeled as a beam with stretching modulus of 1,100 pN, bending modulus of 230 pN nm2, and twisting modulus of 460 pN nm2. Each bp in the duplex is modeled as a FE node with mass of 652 Da, which represents the average molecular weight of a single bp, and two consecutive bps in the same duplex are connected by a two-node Hermitian beam element with the aforementioned geometric and mechanical properties. Four-way junctions are assumed to exist only in “stacked-X” conformations that occur with divalent ions present, in which coaxial stacking occurs between two pairs of arms,77 consistent with previous modeling assumptions in which each pair of stacked arms is modeled as a continuous DNA double-helix.5 Each “stacked-X” four-way junction (or equivalently, double-stranded crossover) is therefore modeled using a single torsional spring, assuming a ground-state conformation with interhelical distance of 1.85 nm and between 0 and 60 degree right-handed scissor-like interhelical angle Jtwist. The torsional stiffness associated with Jtwist is assumed here to be 135 pN nm rad–1.5 Single crossovers are modeled similarly except with a zero torsional stiffness. The stiffness and mass matrices of the DNA nanostructures are computed for this FE model in the ground-state conformation using the commercial FE analysis software ADINA Version 9.0 (ADINA R&D, Inc., Watertown, MA, USA) (Figure 1).5, 15 Single-stranded DNA elements such as hairpins and

ACS Paragon Plus Environment

10

Page 11 of 49

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

loops are ignored. Nicks are treated as regions with the same elastic stiffness as B-form DNA duplexes.

The 6N × 6N friction matrix Z is obtained from a bead-based model assuming 1-nm

hydrodynamic radius beads78 that are positioned at the centroids of the bp positions of the DNA nanostructure, which are coincident with the FE model nodes (Figure 1). Bead models have previously been used to simulate the diffusivity of rigid dsDNA fragments and proteins79-82 by accounting for their molecular structure on solvent-induced drag. To validate our approach, we first used the bead model to calculate the impact of solvent drag on the diffusivity of rigid DNA nanostructures, and subsequently integrated this hydrodynamic model with a structurally compliant FE representation of the DNA nanostructure to develop a coarse-grained BD framework to simulate long time-scale dynamics of structured DNA assemblies.

The 3N × 3N Rotne-Prager tensor, T,65 is used to compute the hydrodynamic interactions

between the beads, 6πµσ 1 I i=j  T T rij rij rij rij σ2 1

1  8πµrij I + + 2  I−  i ≠ j and rij ≥ 2σ 2 2 rij rij 3 rij 2 Tij =  9 rij 3 rij rij T

1  6πµσ 1 − I +  i ≠ j and rij < 2σ 32 σ 32 σrij 

(1)

where σ is the bead hydrodynamic radius, µ = 1 mPa s is the viscosity of water at 20 °C, rij is the vector from the reference point of bp i to the reference point of bp j, rij is its magnitude, and I is

 , which only contains hydrodynamic the 3 × 3 identity matrix. The 3N × 3N friction matrix Z damping coefficients corresponding to the 3N translational DOF, is obtained by inverting T.83

 and T only include translational DOF, we need to compute the 6N × 6N friction Because Z matrix Z , which contains both rotational and translational DOF. In order to do so, the

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 49

components of Z that correspond to the translational DOF are defined to be identical to those of

 matrix, and all other components of Z are set to zero. This modeling approximation the Z assumes that only point hydrodynamic forces occur on beads, whereas hydrodynamic torques are absent. Notwithstanding, overall hydrodynamic torques are present on nanostructures due to the coupling between structural DOF.

2.2

FE Vacuum Normal Mode Analysis

The FE-based model predicts the 3D equilibrium structure of programmed DNA assemblies with a single dominant mechanical ground-state.5,

15

Vacuum NMA yields the root-mean-square

fluctuations (RMSF) of each bp, which is based on the following governing equation,15 Mx + Kx = 0

(2)

where x is the 6N-dimensional acceleration vector of the DNA nanostructure. Eq. 2 can be transformed to the frequency domain to yield the generalized eigenvalue problem, Kϕ ϕi = ωi 2 Mϕ ϕi ϕi T Mϕ ϕj = δij

(3) (4)

where ϕi is the i th normal mode (NM), ωi is the square root of the i th eigenvalue, δij is the Kronecker delta, and the eigenvalues are sorted in ascending order. The RMSF associated with the kth DOF, 〈xk 2 〉 , is computed using the Equipartition Theorem of statistical thermodynamics,20, 84

ACS Paragon Plus Environment

12

Page 13 of 49

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

〈xk 2 〉 =  6N

i7

kB T 2 ϕ ωi 2 ik

(5)

in which xk and ϕik are the kth element of x and ϕi , respectively, kB is the Boltzmann constant, and T is temperature. While vacuum NMA yields the correct magnitude of the local and global equilibrium thermal fluctuations of macromolecular assemblies,20, 85 the frequencies and modeshapes associated with vacuum NMA are unphysical because solvent damping is ignored.

2.3

FE Langevin and Brownian Dynamics

In contrast to vacuum NMA, Langevin Dynamics and BD model the effects of solvent frictional forces or damping, embodied in Z, and thermally induced fluctuations due to the solvent heat bath using the stochastic force vector f(t). The time-dependent conformational dynamics of DNA nanostructures can be modeled using the Langevin equation,86 Mx + Zx + Kx = f(t)

(6)

where the velocity vector of the DNA assembly is denoted x. The external stochastic force vector satisfies the following conditions, 〈f(t)〉 = 0

〈f(t) f(t')T 〉 = 2kB T Z δt − t'

(7) (8)

where δ is the Dirac delta function. The governing equations for BD are obtained by neglecting the inertial term in Eq. 6, which is valid in the over-damped limit, Zx + Kx = f(t),

(9)

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 49

when inertial effects can be ignored.75, 87 This modeling approximation is expected to be valid for long time-scale and large length-scale thermally induced motions. Similar to vacuum NMs, the Brownian modes (BMs) are computed from the generalized eigenvalue problem, Zφi = τi Kφi

φi T Kφj = δij

(10) (11)

in which φi is the ith BM and τi is the ith eigenvalue and relaxation time associated with φi . The governing equations for the BD model in Eq. 9 can also be cast as, Zx + Kx = Cdw  C = !C 0" 0 0

(12) (13)

 is a 3N × 3N matrix resulting from the Cholesky decomposition where C is a 6N × 6N matrix, C  ,50 of 2kB T Z

C  T = 2kB T Z , C

(14)

and dw is a 6N-dimensional random vector whose components are obtained from a real-valued

Gaussian distribution with zero mean and variance 1⁄∆t,88 where ∆t is the time-step size. The BD trajectory of the nucleic acid nanostructure is simulated by solving Eq. 12 in time using the implicit time-integration trapezoidal rule.51

$, To implement the time-integration scheme we first compute K $ = K + 2⁄∆t Z K

(15)

ACS Paragon Plus Environment

14

Page 15 of 49

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

which is triangularized using,51 $ = LDLT K

(16)

and initialize the displacement and velocity vectors at time 0. We then compute the displacement

vector xt%∆t and velocity vector xt%∆t at time t + ∆t using the BD governing equations at time t + ∆t,

Zxt%∆t + Kxt%∆t = Cdwt%∆t

(17)

where dwt%∆t is the value of the stochastic forcing vector at time t + ∆t. Calculation of xt%∆t and xt%∆t is comprised of the following three steps:

$ t%∆t at time t + ∆t using, i) Calculate the effective load vector R $ t%∆t = Cdwt%∆t + Z&2⁄∆t xt + xt '. R

(18)

ii) Solve for the displacement vector xt%∆t at time t + ∆t using,51 $ LDLT xt%∆t = R

t%∆t

.

(19)

iii) Calculate the velocity vector xt%∆t at time t + ∆t using, xt%∆t = 2⁄∆t xt%∆t − xt − xt .

(20)

This time integration scheme is an unconditionally stable implicit time-integration scheme for the linear BD dynamical problem solved here, so that the time-step size does not impact solution stability and therefore only needs to be chosen to obtain sufficient accuracy, which is a subjective choice of the modeler.51 Ensuring that the time-step size is at least five times smaller

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 49

than the relaxation time of a given BM of interest will in general account for its dynamical contribution in the BD simulation.

3

Methods

3.1

Brownian Dynamics

The ground-state three-dimensional structure of a previously published DNA-based “tweezer” was computed using the off-lattice CanDo model5, 15 assuming default geometric and mechanical properties of B-form DNA4, 89 (Figure 2a and Table S1). Double-stranded crossovers J1, J2, J3, J5, and J6 in the tweezer adopt angles that deviate from their ground-state values that were prescribed from all-atom simulations in which they were computed to be –13.2°, –23.3°, –21.3°, –23.3°, and 1.0°, respectively, where a positive sign indicates right-handed twist (Figure 2a). The coarse-grained FE model contains N = 212 FE nodes and 6N = 1,272 DOF.

The first 20 non-rigid-body BMs and corresponding relaxation times were calculated from the K and Z matrices by solving the generalized eigenvalue problem in Eqs. 10 and 11. The time-step size was chosen to be 4.0 picoseconds (ps), which is 10–4 times the longest relaxation time of the model corresponding to the first BM. The total simulation time was chosen to be 150 times the longest relaxation time, which corresponds to 6.0 µs, to ensure the dynamics of the first BM of interest is captured in the simulation. Ten replicate simulations were performed and sampled every ten time-steps to compute statistical quantities, where each replicate was run on a single CPU core (Intel Xeon CPU E5-2670 v2) in approximately 32 hours. The over-damped conformational dynamics of a previously published nine-layer ring origami object that contains 3,600 bps was also simulated using BD.5,

90

The structure was

modeled using 3,740 FE nodes total, with 22,440 total DOF in the FE model. The integration

ACS Paragon Plus Environment

16

Page 17 of 49

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

time-step was chosen to be 180 ps, which is 10–4 times the longest relaxation time of the model, and the total BD simulation time was set to 270 µs in order to capture the longest relaxation time. Ten replicates were simulated, with each replicate computed on a single CPU core (Intel Xeon CPU E5-2670 v2) in 190 hours. Conformations were sampled every ten time-steps to produce trajectories of 150,000 frames each. The last system simulated is a pointer-like origami object that consists of 82 doublehelices total and consists of the full-length M13mp18 phage genomic DNA as scaffold strand.91 The ground-state structure of the pointer was computed using the same model employed previously.2,

15

Each single-stranded or double-stranded crossover in the pointer origami was

modeled as a rigid beam with length of 2.25 nm, as previously published.15, 91 The FE model

consists of N = 7,249 FE nodes total, corresponding to 6N = 43,494 translational and rotational DOF. The first ten BMs were computed and presented in results, with the BD trajectory left to future work.

3.2

Molecular Dynamics

The all-atom PDB file for the tweezer model was generated from CanDo assuming default geometric and mechanical properties of B-form DNA, as in the BD simulations, and subsequently ionized and solvated in preparation for MD. The all-atom system was immersed in TIP3P92 water (109,204 water molecules) and explicit Mg2+ and Cl– ions were added to neutralize DNA charges and to bring the Mg2+ ion concentration to 12.5 mM in the simulation cell while maintaining electro-neutrality of the simulation cell. Energy minimization was performed for 10,000 steps with the conjugate gradient and line search algorithm implemented in NAMD2 prior to MD equilibration and production.

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 49

All simulations were performed using NAMD293 employing the CHARMM27 force field,32, 33 and Allnér Mg2+ parameters94 with an integration time-step size of 2 femtoseconds (fs). A periodic orthogonal cell was used with dimensions 105 Å × 232 Å × 126 Å. Van der Waals energies were computed using a 12 Å cut-off with a switching function (10–12 Å). The Particle Mesh Ewald method95 was used to compute full electrostatics with a minimum grid-point density of 1 Å–1. Full electrostatic forces were computed every two time-steps (equal to every 4 fs) and non-bonded forces were calculated each time-step (2 fs). Simulations were performed in the NpT ensemble using the Nose-Hoover Langevin piston method96, 97 for pressure control (1 atm), with an oscillation period of 200 fs and a damping time of 100 fs. Langevin forces98 were applied to heavy atoms for temperature control (300 K) with coupling coefficients of 5 ps–1. All hydrogens were constrained to their equilibrium lengths and system snapshots were recorded every 1 ps. Values for the six in-plane scissor motions, Jtwist,99 of the system’s four-way junctions J1– J6 were calculated at each frame of the trajectory, from which a probability density was obtained for each of the six junctions. For each junction, the duplex arms of the four-way junction were defined to be the ten bps flanking the strand crossover site on each side of the junction. In a “stacked-X” four-way junction, the four duplex arms form two helices, denoted by H1 and H2. The program Curves+100 calculates an axial curve segment for each helix. A line segment, denoted by L1 or L2, respectively, was fitted to the axial curve segment for H1 or H2. Three vectors were then computed, with the first defined as the vector from the endpoint of L1 at the 5′-end of the non-crossing strand to the midpoint of L1; the second defined as the vector from this midpoint of L1 to the midpoint of L2; and the third defined as the vector from the midpoint of L2 to the endpoint of L2 at the 3′-end of the non-crossing strand. The dihedral angles were then used to define the Jtwist angles, where positive values denote right-handed configurations

ACS Paragon Plus Environment

18

Page 19 of 49

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

and negative values denote left-handed configurations. Single Gaussians were fitted to each distribution of angles and the mean values used as the ground-state angles for the BD simulations. The atomic model of the DNA tweezer was simulated using explicit water and ions, and consisted of 338,556 frames with a time difference of 1.0 ps between each consecutive frame. Two bps at each end were neglected in the downstream analysis because of bp fraying in the MD simulation (Supporting Text). Thus each MD frame has 1,176 DOF after coarse-graining. The frames corresponding to the first 50 ns of the MD trajectory were removed in the downstream analysis for equilibration purposes. The MD trajectory was computed using 100 CPU cores (Intel Xeon CPU E5-2670 v2) plus ten GPU units (NVIDIA K20 with 2,496 CUDA cores each) in 36 days.

3.3

Conversion between All-Atom and FE Models In order to compare the MD and BD trajectories, each frame in the MD trajectory was

coarse-grained by calculating the displacements corresponding to the three nodal translational DOF and three nodal rotational DOF. First, the program Curves+100 was used with 3DNA convention76 to coarse-grain the j th bp in frame i in the MD trajectory by calculating the reference point as a 3D column vector bij and three right-handed orthonormal reference axes as

three columns of a 3 × 3 matrix Bij . Second, a reference FE structure of the same DNA tweezer

in the mechanical ground state was generated using the lattice-free FE model,5 which gave the

reference point of the jth bp as a 3D column vector rj and reference axes as a 3 × 3 orthogonal

matrix Rj . Finally, the Kabsch algorithm was applied to align the reference points rj , j =

1, 2, 3, … to those bij , j = 1, 2, 3, … in frame i in the MD trajectory, yielding the rigid-body

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 49

translation and rotation from rj , j = 1, 2, 3, … to bij , j = 1, 2, 3, … as a 3D column vector ai and an orthogonal matrix Ai , respectively. The displacements corresponding to the three translational DOF of the jth bp in frame i were given by, uij = Ai T bij − ai − rj .

(21)

The displacements corresponding to the three nodal rotational DOF of the same bp were given by first representing the rotation matrix, Uij = Ai T Bij Rj T

(22)

as a normalized rotation axis e and a rotation angle θ and then multiplying θ with the three elements of e, u( ij =θe1

θe2

θe3 T .

(23)

This procedure enables coarse-graining the atomic positions in frame i to a 6N-dimensional

displacement vector ui in which N is the number of FE nodes. The 6N × 6N covariance matrix

was computed using all frames u1 , u2 , …, uM, of which the eigenvalues σ1 2 , σ2 2 , …, σ6N 2 were sorted in descending order and the corresponding orthonormal eigenvectors ψ1 , ψ2 , …, ψ6N were obtained. The eigenvector ψi is defined as the ith essential mode (EM), and the eigenvalue σi 2 is the variance of the trajectory projected onto ψi . In addition, each frame in the BD trajectory was used to generate the corresponding allatom model using the same approach presented in previous works5, 9. Briefly, four standard reference atomic structures of the Watson-Crick bps A-T, T-A, G-C, and C-G, including backbone atoms, were built using the software Accelrys DS Visualizer Version 3.5 (Accelrys,

ACS Paragon Plus Environment

20

Page 21 of 49

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

Inc., Burlington, MA). For each of these atomic structures, the software 3DNA76 gave the bp reference point and three orthogonal reference directions that point from the minor groove to the major groove, from one base to its complementary base, and perpendicular to the bp-plane, respectively. Without loss of generality, a right-handed orthonormal reference frame e0 , e1 , e2 , e3 was generated for each standard reference atomic structure using the above 3DNA

convention. For example, the reference frame of the standard reference atomic structure of the AT bp has the center e0 located at the bp reference point, axis e1 pointing to the major groove, axis e2 pointing to the adenine, axis e3 perpendicular to the bp plane in the 5′ to 3′ direction of the DNA strand containing the adenine. Each bp in the all-atom model was generated by applying rigid-body translation and rotation to the specific standard reference atomic structure, including atoms in the backbones, such that its bp reference point and reference axes coincide with those of the FE node.

4.

Results and Discussion

The DNA tweezer consists of two “arms” that are predicted by the BMs to exhibit out-of-plane and in-plane scissoring motions, with relaxation times corresponding to 42 and 32 ns (Figures 2b, 2c and S1). The autocorrelation function (ACF) of the BD motion projected onto each BM is expected theoretically to decay exponentially with a time-constant that is equal to the relaxation time corresponding to the BM. The time constants of these exponential functions were computed from BD to be similar to their corresponding BM relaxation times, providing support for the validity of the BD trajectory (Figures 2b, 2c, 3a and S6). Interestingly, the corresponding time constants of the ACFs of the projections of the BD trajectory onto the corresponding NMs (Figure S7) and MD-EMs (Figure S8) are distinct from those of the BMs (Figure S6), which is

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 49

likely attributable to the differences between their shapes (Figures 3b, S1–S5). In contrast, the time-constants obtained by projecting the MD trajectory onto the MD EMs (28 ns and 30 ns, respectively, for the first two modes) largely agree with those of the corresponding BD trajectory and BMs (36 ns and 28 ns, respectively, for the first two modes), but are distinct for higher modes. Differences in the relaxation times of higher modes may be attributable to a number of factors including the fact that higher mode-shapes may involve conformational dynamics that the simplified coarse-grained FE model does not treat, such as non-bonded electrostatic and van der Waals interactions, as well as nonlinear DNA dynamics involving four-way junctions and interstrand hybridization (Figure 3a). In equilibrium, the RMSFs of the coarse-grained BD model are expected to agree quantitatively with vacuum NMA results by the equipartition theorem (Figure 3c). The translational and rotational diffusion coefficients of the tweezer are also computed to be 41 µm2 s 1 and 1.1 × 106 s 1 , respectively, using the Z matrix and the bead positions, as previously performed.101 A basic utility of the BD simulation is to obtain realistic dynamics of geometric properties of DNA nanostructures. The practical value of the BD trajectory of the DNA tweezer is demonstrated here by studying the time evolution of the distance between two atoms in the DNA tweezer, to which an enzyme and its co-factor, respectively, are attached. In particular, the function of the DNA tweezer is to scaffold the enzyme and co-factor at the ends of each of its two arms4 to bring them into proximity in an extrinsically controlled reaction. In the closed-state of the tweezer, the distance between the two atoms, chain C, residue 1, atom O5′; and chain H, residue 1, atom O5′, can be monitored to compare the coarse-grained and all-atom conformational results (Figure 4a). The Euclidean distance d between the two atoms is somewhat lower in the BD trajectory, with a mean of 3.3 nm and a standard deviation of 1.5 nm compared

ACS Paragon Plus Environment

22

Page 23 of 49

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

with the corresponding distance from the MD trajectory that has a mean of 4.7 nm and a standard deviation of 0.8 nm (Figure 4b–4d). The blocking method of Flyvbjerg and Petersen102 estimates the standard errors of the mean distances from the BD and MD trajectories to be 0.05 nm and 0.2 nm, respectively (Figure 4f). This disagreement in mean distance and fluctuation magnitude is not unexpected given the coarse-grained nature of the FE model employed. Notwithstanding, the relaxation time of the distance d is predicted from the BD trajectory to be 21.8 ns compared with 15.8 ns predicted by MD (Figure 4e). Similarly, deviation is present between the probability distribution functions for the six junction angles Jtwist calculated from BD and MD (Figure 5a and 5b), as well as their relaxation times (Figure S11). It is likely that the equilibrium angles of the flexible four-way junctions computed using the BD model deviate from the input, ground-state values obtained from MD because of mechanical coupling between the junctions that is mediated by intervening DNA duplexes. Interestingly, cross-correlations between distinct junction angles are predicted to be slightly greater in the all-atom model, which may be due to nonlinear motions captured by the latter, or another underlying difference between the two models (Figure 5c). Taken together, large-scale dynamical motions of this relatively small DNA construct are reasonably well captured both in their conformational nature and also temporal dynamics by the coarse-grained FE model. Next we applied the BD model to a considerably higher molecular weight DNA object whose dynamical simulation is out of reach of all-atom representations using standard computing resources. The translational and rotational diffusion coefficients of the nine-layer ring origami are, respectively, calculated to be 10 µm2 s 1 and 1.1×104 s 1 . The rigid-body rotational diffusivity of the origami is therefore predicted to be two orders of magnitude slower than that of the tweezer, whereas the translational diffusion coefficient of the origami is only four-fold lower

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 49

than that of the tweezer. The longest relaxation time of the nine-layer ring corresponding to its first BM is 1.8 µs, which is over 40 times longer than that of the tweezer. Similar to the DNA tweezer, the Brownian motion of the nine-layer ring origami is dominated by the first two BMs, which correspond to out-of-plane twisting and bending of this disk-like structure, respectively (Figure 6). Overall fluctuations of the structure measured by the RMSFs of the bps calculated from the BD simulation are also similar to those obtained using NMA (Figure 7c). Because the effects of solvent damping are distributed uniformly on the object and in the same manner as its mass, the first several NMs and BMs are anticipated to be similar (Figures 7b, S12–S14). This is in contrast to the tweezer (Figures 3b, S1–S2, and S4), for which the impact of solvent damping on the first modes is expected to be considerably more complex due to its smaller size and irregular geometry. As a result, the first ten relaxation times of the nine-layer ring object computed from its BD trajectory projected onto its BMs are highly similar to those projected onto the corresponding NMs (Figures 7a, S15 and S16). Full relaxation of the first 20 modes is observed on the 270 µs time-scale simulated (Figure S15). In practice, the nanopore in the center of the nine-layer ring origami may have useful applications in nano-biotechnology,103 and a major geometric feature of the nanopore is its diameter. To gain insight into its dynamics using the BD trajectory, we characterized its diameter over time. The pore diameter, D, is defined as the distance between the reference points of two bps separated by 100 bps in the innermost layer of the ring origami, which contain nucleotides number 1,846 and 1,946 of the scaffold strand (Figure 8a). Here the counting of nucleotides is from the first paired nucleotide in the 5′-direction and starts from one. Note that the pore diameter may be defined by the distance between any pair of bps separated by 100 bps in the innermost layer without significantly changing the results presented below. The temporal

ACS Paragon Plus Environment

24

Page 25 of 49

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

evolution of the pore diameter is seen to fluctuate rapidly between approximate minimum and maximum values of 15.9 and 26.9 nm, with a mean value of 21.3 nm that is similar to the mechanical ground-state value of 21.3 nm, and a standard deviation of 1.2 nm (Figure 8b–8d). The standard error of the mean is estimated to be 0.02 nm using the blocking method (Figure 8f). The ACF of the pore diameter D is fitted using a double exponential model, which reports a slow relaxation time of 1.3 µs, and a faster relaxation time of 10.7 ns (Figure 8e). The slow relaxation time is similar to those corresponding to the first two BMs, which affect the pore diameter, D (Figure 7a). Finally, in application of the BD framework to a full-scale scaffolded origami object we computed the relaxation times of the first ten BMs of the previously published pointer origami object (Figure 9a). Results indicate that the object’s BD is dominated by its first two BMs that consist of global twisting and bending (Figure 9b), with highly similar shapes indicated by the overlap matrix (Figure 9c). In contrast to the tweezer and nine-layer ring origami in which the DNA helices are arranged in 2D planes, the pointer comprises helices arranged in a 3D solid lattice arrangement (Figure 9a). The effect of solvent damping is therefore distributed uniformly on the object, with interior duplexes shielded hydrodynamically. This is in contrast to the preceding origami object in which all duplexes are subject to solvent damping, and results in vacuum NMs and BMs of the pointer that are dissimilar (Figure 9c). The translational and rotational diffusion coefficients of the pointer were predicted from the model to be 14 µm2 s 1

and 4.5 × 104 s 1, respectively. While the number of bps in the pointer is twice that of the ninelayer ring object, its rigid-body rotational diffusivity is almost four times faster due to its

compact nature. The longest relaxation times corresponding to the first two BMs are τ1 = 93 ns

and τ2 = 91 ns, which are only approximately twice the longest relaxation time of the

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 49

considerably smaller DNA tweezer and approximately 20 times lower than that of the nine-layer ring origami. The relatively fast relaxation times of the slowest BMs of the pointer are likely due to the fact that it is a dense 3D solid object with a high density of internal crossovers that results in a compact object of high rigidity. While the temporal dynamics of the pointer structure are not simulated explicitly here, the BM relaxation times that describe its internal motions could be used to alter or optimize these or other specific motions of interest. For example, crossover density might be reduced to alter twisting time-scales of the structure, or its cross-section altered to achieve the same. Application of the present BD framework to this and other practical design applications of DNA nanostructures will be of interest to examine in future work.

5.

Concluding Discussion

Structured nucleic acids serve as versatile nanoscale scaffolds to control biomolecular assemblies including multi-enzyme synthesis cascades,4 energy transport cascades,9,

14

inorganic particle

synthesis molds,6 and cellular sensors.3 The FE method is a practical coarse-grained modeling approach to predict the 3D structure of these high molecular weight assemblies in a computationally efficient manner.2, 5, 15 However, previous FE models are unable to compute the temporal dynamics of these assemblies, limiting their application to the analysis of equilibrium conformational fluctuations. Here, we overcome this limitation by incorporation of a bead-based model that has previously been applied to simulate long time-scale conformational dynamics of linear DNA chains.57, 58, 75, 87 This results in a general computational methodology that enables the long time-scale simulation of up to microseconds of structured nucleic acid assemblies that are thermodynamically stable about a single, dominant ground-state conformation. In addition to calculation of time-course trajectories using BD, the framework enables rapid calculation of the

ACS Paragon Plus Environment

26

Page 27 of 49

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

rigid-body translational and rotational diffusivity of DNA assemblies, as well as their BMs that indicate the conformational nature and relaxation times associated with their long time-scale dynamics. While large length-scale and long time-scale dynamics are reasonably well captured by the coarse-grained model, more subtle conformational dynamics involving inter-duplex distance and junction angle correlations are only approximately captured, requiring all-atom descriptions to treat properly. Ignoring non-bonded electrostatic and van der Waals interactions additionally limits application of the model to highly structured assemblies with high crosslink density as typically designed in DNA origami.2,

5, 15

Notwithstanding these limitations, the

present computational tool enables rapid feedback on conformational dynamics. For example, a time-scale of 6.0 μs is simulated for the tweezer in only 32 hours on a single CPU core, which

compares with 36 days computation time using 100 CPU cores plus ten GPU units for the equivalent MD model. As always, the level of modeling detail required to address a specific question of interest needs to be balanced with associated computational cost, and the present coarse-grained framework should offer a viable alternative for specific applications in which overall conformational dynamics of high molecular weight DNA assemblies are sought without resorting to computationally intensive MD. Future applications of this framework may include hybrid approaches that benefit from all-atom descriptions for short length-scale and fast timescale motions together with use of this coarse-grained FE model for large length-scale and slow time-scale motions. Such approach may be useful to enable the simulation of coupled, multiscale conformational dynamics that impact dynamical processes such as excitonic energy transfer9,

14

and multi-enzyme reactions4. The current FE-based BD framework can also be

applied with suitable modeling modifications to proteins and protein assemblies,20, 104 which may be the subject of future work. The source code for the BD simulation is available for use by the

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

Page 28 of 49

broader community as MATLAB run-time files available in the Supporting Information and at http://cando-dna-origami.org.

Associated Content

Supporting Information. Sequence design of the DNA tweezer, discussion on the bp fraying in the MD simulation, comparison of analytical and numerical Brownian and vacuum eigensolutions, and other supporting figures. Also included are the movies showing the BD and MD trajectories and the source code for the BD simulation. This material is available free of charge via the Internet at http://pubs.acs.org.

Author Information Corresponding Authors *Email: [email protected] and [email protected]

Author Contributions #

These authors contributed equally to this work. R.S.S., O.H., K.J.B., and M.B. conceived of the

coarse-grained BD approach for DNA assemblies. R.S.S. developed the BD simulation framework and verified its performance against analytical and other numerical solutions. R.S.S. implemented the pre-simulation and BD simulation scripts. R.S.S. and K.P. implemented the scripts for post-simulation data analysis. K.P. built the FE models and atomic models for the DNA assemblies, performed the final BD simulations, and processed the simulation results to generate the figures and movies. M.R.A. performed and analyzed the MD simulation for the

ACS Paragon Plus Environment

28

Page 29 of 49

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

DNA tweezer. All authors discussed the results and their presentation. R.S.S., K.P., M.R.A., and M.B. wrote the manuscript text. All authors commented on and edited the manuscript.

Notes The authors declare no competing financial interest.

Acknowledgements Funding from the Army Research Office Multi-University Research Initiative W911NF1210420 to RSS, MRA, and MB, from the Office of Naval Research N000141210621 to KP and MB, and from the Simons Foundation to OH are gratefully acknowledged.

Figure Legends Figure 1. BD simulation framework for DNA nanostructures. The all-atom DNA structure is modeled using a FE representation5, 15 that yields the FE stiffness and mass matrices of the DNA assembly. Beads positioned at bps are used to model hydrodynamic drag and interactions. These matrices provide input to the Langevin and BD models used to simulate long time-scale conformational dynamics of the DNA object in solvent.

Figure 2. BM calculation of the DNA tweezer. (a) (Upper) Topological routing of the nine single-stranded DNAs of the DNA tweezer. The nucleotide sequences are available in Table S1. The 3′-end of each strand is represented by the arrowhead, and the identity of each strand is marked near its 5′-end. Six four-way junctions are denoted by J1–J6, in which J4 is a singlestranded crossover and the remaining five are double-stranded crossovers. (Lower) Three

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

Page 30 of 49

orthogonal views of the atomic model of the DNA tweezer in its mechanical ground state. (b) Analysis of the first BM. (Upper) Analytical and numerical ACFs of the BD in the first BM. This panel presents the analytical ACF (BM), the numerical ACF calculated from the simulated BD trajectories (BD) in which the shaded area represents the 95% confidence interval, and an exponential model fitted to the numerical ACF (BD Fit). (Lower) The mode shape is visualized using the same rendering scheme as in (a). (c) Analysis of the second BM. (Upper) Analytical and numerical ACFs of the BD in the second BM. (Lower) The corresponding mode shape.

Figure 3. Comparison between BM, BD, NMA, and MD simulations of the DNA tweezer. (a) Relaxation times of the first ten BMs, including the analytical results (BM), the numerical results obtained by fitting exponential functions to the ACFs of the simulated BD trajectories projected onto the first ten BMs (BD Fit), NMs (NM), and MD EMs (EM). Also shown are the numerical results obtained by fitting exponential functions to the ACFs of the simulated MD trajectory projected onto the first ten MD EMs (MD Fit). (b) Overlap matrix between the first ten BMs and NMs, where the ijth element of the overlap matrix is calculated by the dot-product of the ith normalized BM with the jth normalized NM. (c) RMSFs of each FE node of the DNA tweezer calculated from NMA and the simulated BD trajectories. The Pearson correlation coefficient is R = 1.0.

Figure 4. Comparison between BD and MD trajectories of the distance d between the arms of the DNA tweezer. (a) Definition of the distance d. The variable d denotes the Euclidean distance between two atoms: chain C, residue 1, atom O5′ (orange sphere) and chain H, residue 1, atom O5′ (green sphere). (b) Time evolution of the distance d in the entire 6.0-µs BD trajectory (blue) and the 338.6-ns MD trajectory (red). (c) A plot focusing in on the region between 300 ns and

ACS Paragon Plus Environment

30

Page 31 of 49

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

330 ns in the BD and MD trajectories in (b). (d) Probability density functions of the distance d calculated from the entire simulated BD (blue) and MD (red) trajectories. The vertical lines denote the mean values. (e) Numerical ACFs of the distance d calculated from the simulated BD (BD) and MD (MD) trajectories. Also included are two exponential models fitted to the numerical ACFs of BD (BD Fit) and MD (MD Fit). (f) The standard errors, σ(d), of the mean distances are estimated using the blocking method, in which n denotes the number of block transformations, and the shaded area represents the 95% confidence interval.

Figure 5. Comparison between BD and MD trajectories of the junction angles Jtwist of the DNA tweezer. (a) Definition of the junction angle Jtwist. (Upper) Topological routing of a four-way junction in a “stacked-X” conformation. (Lower) Two orthogonal views of the atomic model of the junction in its mechanical ground state. The junction angle Jtwist is defined as the angle between the two helices, where positive and negative signs indicate right-handed and left-handed twist, respectively. (b). Probability density functions of the six junction angles denoted by Jtwist at Junctions J1–J6. (c) Pearson correlation coefficients between the two sets of the density functions of the six junction angles calculated from the BD and MD trajectories.

Figure 6. BM calculation of the nine-layer ring origami. (a) Three orthogonal views of the atomic model of the nine-layer ring origami in the mechanical ground state. (b,c) Mode shapes of the first (b) and second (c) BMs. The mode shapes are visualized using the same rendering scheme as in (a). (d) Analytical and numerical ACFs of the BD trajectories in the first (upper) and second (lower) BMs. Each plot presents the analytical ACF (BM), the numerical ACF calculated from the simulated BD trajectories (BD) in which the shaded area represents the 95% confidence interval, and an exponential function which is fitted to the numerical ACF (BD Fit).

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 49

Figure 7. Comparison between BM, BD, and NMA calculations of the nine-layer ring origami. (a) Relaxation times of the first ten BMs, including the analytical results (BM), the numerical results obtained by fitting exponential functions to the ACFs of the simulated BD trajectories projected onto the first ten BMs (BD Fit) and NMs (NM). (b) Overlap matrix between the first ten BMs and NMs. (c) RMSFs of each FE node of the nine-layer ring origami calculated from NMA and the simulated BD trajectories, with a Pearson correlation coefficient of R = 1.0.

Figure 8. Analysis of the pore diameter, D, of the nine-layer ring origami. (a) Definition of the pore diameter as the distance between the reference points of two bps containing Nucleotides 1,846 (green sphere) and 1,946 (orange sphere) of the scaffold strand. (b) Time evolution of the pore diameter, D, in the entire 270-µs BD trajectory. (c) A plot focusing in on the region between 269.5 µs and 270 µs in the BD trajectory in (b). (d) Probability density function of the pore diameter computed from the simulated BD trajectory. The vertical line denotes the mean value. (e) The ACF of the pore diameter D computed from the simulated BD (BD) trajectory. Also included is a double exponential model fitted to the numerical ACF of BD (BD Fit). (f) The standard error, σ(D), of the mean pore diameter is estimated using the blocking method, in which n denotes the number of block transformations, and the shaded area represents the 95% confidence interval.

Figure 9. BM calculation of the 82-helix pointer origami. (a) Three orthogonal views of the atomic model of the pointer origami in the mechanical ground state. The scaffold strand is colored in blue and all staple strands in gray. (b) Shape of the first BM aligned with the groundstate structure. The mode shape is visualized using the same rendering scheme as in (a), and the ground-state structure is in semi-transparent rendering. (c) Analysis of the first ten BMs and

ACS Paragon Plus Environment

32

Page 33 of 49

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

NMs. (Upper) Analytical relaxation times of the first ten BMs. (Lower) Two overlap matrices between the first ten BMs and NMs.

References 1.

Rothemund, P. W. K. Nature 2006, 440, 297-302.

2.

Castro, C. E.; Kilchherr, F.; Kim, D. N.; Shiao, E. L.; Wauer, T.; Wortmann, P.; Bathe, M.; Dietz, H. Nat. Methods 2011, 8, 221-229.

3.

Krishnan, Y.; Bathe, M. Trends Cell Biol. 2012, 22, 624-633.

4.

Liu, M.; Fu, J.; Hejesen, C.; Yang, Y.; Woodbury, N. W.; Gothelf, K.; Liu, Y.; Yan, H. Nat. Commun. 2013, 4, 2127.

5.

Pan, K.; Kim, D. N.; Zhang, F.; Adendorff, M. R.; Yan, H.; Bathe, M. Nat. Commun. 2014, 5, 5578.

6.

Sun, W.; Boulais, E.; Hakobyan, Y.; Wang, W. L.; Guan, A.; Bathe, M.; Yin, P. Science 2014, 346, 1258361.

7.

Castro, C. E.; Su, H. J.; Marras, A. E.; Zhou, L.; Johnson, J. Nanoscale 2015, 7, 59135921.

8.

Dutta, P. K.; Varghese, R.; Nangreave, J.; Lin, S.; Yan, H.; Liu, Y. J. Am. Chem. Soc. 2011, 133, 11985-11993.

9.

Pan, K.; Boulais, E.; Yang, L.; Bathe, M. Nucleic Acids Res. 2014, 42, 2159-2170.

10.

Jaeger, L.; Chworos, A. Curr. Opin. Struct. Biol. 2006, 16, 531-543.

11.

Geary, C.; Rothemund, P. W.; Andersen, E. S. Science 2014, 345, 799-804.

12.

Fu, J.; Yang, Y. R.; Johnson-Buck, A.; Liu, M.; Liu, Y.; Walter, N. G.; Woodbury, N. W.; Yan, H. Nat. Nanotechnol. 2014, 9, 531-536.

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 49

13.

Idan, O.; Hess, H. Nat. Nanotechnol. 2012, 7, 769-770.

14.

Sawaya, N. P.; Huh, J.; Fujita, T.; Saikin, S. K.; Aspuru-Guzik, A. Nano Lett. 2015, 15, 1722-1729.

15.

Kim, D. N.; Kilchherr, F.; Dietz, H.; Bathe, M. Nucleic Acids Res. 2012, 40, 2862-2868.

16.

Laughton, C. A.; Harris, S. A. WIREs Comput. Mol. Sci. 2011, 1, 590-600.

17.

Perez, A.; Luque, F. J.; Orozco, M. Acc. Chem. Res. 2012, 45, 196-205.

18.

Maffeo, C.; Yoo, J.; Comer, J.; Wells, D. B.; Luan, B.; Aksimentiev, A. J. Phys. Condens. Matter. 2014, 26, 413101.

19.

Yoo, J.; Aksimentiev, A. Proc. Natl. Acad. Sci. USA 2013, 110, 20099-20104.

20.

Bathe, M. Proteins 2008, 70, 1595-1609.

21.

Sedeh, R. S.; Bathe, M.; Bathe, K. J. J. Comput. Chem. 2010, 31, 66-74.

22.

Miyashita, O.; Gorba, C.; Tama, F. J. Struct. Biol. 2011, 173, 451-460.

23.

Tama, F.; Valle, M.; Frank, J.; Brooks, C. L., 3rd. Proc. Natl. Acad. Sci. USA 2003, 100, 9319-9323.

24.

Drsata, T.; Lankas, F. WIREs Comput. Mol. Sci. 2013, 3, 355-363.

25.

Lankas, F.; Gonzalez, O.; Heffler, L. M.; Stoll, G.; Moakher, M.; Maddocks, J. H. Phys. Chem. Chem. Phys. 2009, 11, 10565-10588.

26.

Doye, J. P. K.; Ouldridge, T. E.; Louis, A. A.; Romano, F.; Sulc, P.; Matek, C.; Snodin, B. E. K.; Rovigatti, L.; Schreck, J. S.; Harrison, R. M.; Smith, W. P. J. Phys. Chem. Chem. Phys. 2013, 15, 20395-20414.

27.

Potoyan, D. A.; Savelyev, A.; Papoian, G. A. WIREs Comput. Mol. Sci. 2013, 3, 69-83.

28.

Hinckley, D. M.; Freeman, G. S.; Whitmer, J. K.; de Pablo, J. J. J. Chem. Phys. 2013, 139, 144903.

ACS Paragon Plus Environment

34

Page 35 of 49

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

29.

Hsu, C. W.; Fyta, M.; Lakatos, G.; Melchionna, S.; Kaxiras, E. J. Chem. Phys. 2012, 137, 105102.

30.

Morriss-Andrews, A.; Rottler, J.; Plotkin, S. S. J. Chem. Phys. 2010, 132, 035105.

31.

DeMille, R. C.; Cheatham, T. E.; Molinero, V. J. Phys. Chem. B 2011, 115, 132-142.

32.

Foloppe, N.; MacKerell, A. D. J. Comput. Chem. 2000, 21, 86-104.

33.

MacKerell, A. D.; Banavali, N. K. J. Comput. Chem. 2000, 21, 105-120.

34.

Hart, K.; Foloppe, N.; Baker, C. M.; Denning, E. J.; Nilsson, L.; MacKerell, A. D. J. Chem. Theory Comput. 2012, 8, 348-362.

35.

Cheatham, T. E.; Cieplak, P.; Kollman, P. A. J. Biomol. Struct. Dyn. 1999, 16, 845-862.

36.

Perez, A.; Marchan, I.; Svozil, D.; Sponer, J.; Cheatham, T. E., 3rd; Laughton, C. A.; Orozco, M. Biophys. J. 2007, 92, 3817-3829.

37.

Dans, P. D.; Zeida, A.; Machado, M. R.; Pantano, S. J. Chem. Theory Comput. 2010, 6, 1711-1725.

38.

Chirico, G.; Langowski, J. Biopolymers 1994, 34, 415-433.

39.

Srinivas, N.; Ouldridge, T. E.; Sulc, P.; Schaeffer, J. M.; Yurke, B.; Louis, A. A.; Doye, J. P.; Winfree, E. Nucleic. Acids. Res. 2013, 41, 10641-10658.

40.

Ouldridge, T. E.; Sulc, P.; Romano, F.; Doye, J. P.; Louis, A. A. Nucleic. Acids. Res. 2013, 41, 8886-8895.

41.

He, Y.; Maciejczyk, M.; Oldziej, S.; Scheraga, H. A.; Liwo, A. Phys. Rev. Lett. 2013, 110, 098101.

42.

Trovato, F.; Tozzini, V. J. Phys. Chem. B 2008, 112, 13197-13200.

43.

Machado, M. R.; Dans, P. D.; Pantano, S. Phys. Chem. Chem. Phys. 2011, 13, 1813418144.

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

44.

Page 36 of 49

Snodin, B. E. K.; Randisi, F.; Mosayebi, M.; Sulc, P.; Romano, J.; Romano, F.; Ouldridge, T. E.; Tsukanov, R.; Nir, E.; Louis, A. A.; Doye, J. P. K. J. Chem. Phys. 2015, 142, 234901.

45.

Ouldridge, T. E.; Louis, A. A.; Doye, J. P. K. Phys. Rev. Lett. 2010, 104, 178101.

46.

Romano, F.; Chakraborty, D.; Doye, J. P. K.; Ouldridge, T. E.; Louis, A. A. J. Chem. Phys. 2013, 138, 085101.

47.

Schreck, J. S.; Ouldridge, T. E.; Romano, F.; Louis, A. A.; Doye, J. P. K. J. Chem. Phys. 2015, 142, 165101.

48.

Ouldridge, T. E.; Hoare, R. L.; Louis, A. A.; Doye, J. P. K.; Bath, J.; Turberfield, A. J. ACS Nano 2013, 7, 2479-2490.

49.

Ouldridge, T. E.; Louis, A. A.; Doye, J. P. K. J. Chem. Phys. 2011, 134, 085101.

50.

Ermak, D. L.; McCammon, J. A. J. Chem. Phys. 1978, 69, 1352-1360.

51.

Bathe, K. J., Finite Element Procedures, second edition. KJ Bathe, Watertown, MA, 2014.

52.

Fixman, M. Macromolecules 1986, 19, 1204-1207.

53.

Geyer, T.; Winter, U. J. Chem. Phys. 2009, 130, 114905.

54.

Ando, T.; Chow, E.; Saad, Y.; Skolnick, J. J. Chem. Phys. 2012, 137, 064106.

55.

Ando, T.; Chow, E.; Skolnick, J. J. Chem. Phys. 2013, 139, 121922.

56.

Saadat, A.; Khomami, B. J. Chem. Phys. 2014, 140, 184903.

57.

Klenin, K.; Merlitz, H.; Langowski, J. Biophys. J. 1998, 74, 780-788.

58.

Fayad, G. N.; Hadjiconstantinou, N. G. Microfluid. Nanofluid. 2010, 8, 521-529.

59.

Jian, H. M.; Vologodskii, A. V.; Schlick, T. J. Comput. Phys. 1997, 136, 168-179.

60.

Allison, S. A. Macromolecules 1986, 19, 118-124.

ACS Paragon Plus Environment

36

Page 37 of 49

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

61.

Allison, S.; Austin, R.; Hogan, M. J. Chem. Phys. 1989, 90, 3843-3854.

62.

Allison, S. A.; Sorlie, S. S.; Pecora, R. Macromolecules 1990, 23, 1110-1118.

63.

Chirico, G.; Langowski, J. Biophys. J. 1996, 71, 955-971.

64.

Jian, H. M.; Schlick, T.; Vologodskii, A. J. Mol. Biol. 1998, 284, 287-296.

65.

Rotne, J.; Prager, S. J. Chem. Phys. 1969, 50, 4831-4837.

66.

Fayad, G. N.; Hadjiconstantinou, N. G. J. Fluid. Eng.-T. ASME 2013, 135, 024501.

67.

Wocjan, T.; Krieger, J.; Krichevsky, O.; Langowski, J. Phys. Chem. Chem. Phys. 2009, 11, 10671-10681.

68.

Wocjan, T.; Klenin, K.; Langowski, J. J. Phys. Chem. B 2009, 113, 2639-2646.

69.

Vologodskii, A. Biophys. J. 2006, 90, 1594-1597.

70.

Dlugosz, M.; Zielinski, P.; Trylska, J. J. Comput. Chem. 2011, 32, 2734-2744.

71.

Dlugosz, M.; Huber, G. A.; McCammon, J. A.; Trylska, J. Biopolymers 2011, 95, 616627.

72.

Dlugosz, M.; Trylska, J. BMC Biophys. 2011, 4, 3.

73.

Brackley, C. A.; Morozov, A. N.; Marenduzzo, D. J. Chem. Phys. 2014, 140, 135103.

74.

Plimpton, S. J. Comput. Phys. 1995, 117, 1-19.

75.

Jian, H.; Vologodskii, A. V.; Schlick, T. J. Comput. Phys. 1997, 136, 168-179.

76.

Lu, X. J.; Olson, W. K. Nucleic Acids Res. 2003, 31, 5108-5121.

77.

McKinney, S. A.; Declais, A. C.; Lilley, D. M. J.; Ha, T. Nat. Struct. Biol. 2003, 10, 9397.

78.

Fernandes, M. X.; Ortega, A.; Lopez Martinez, M. C.; Garcia de la Torre, J. Nucleic Acids Res. 2002, 30, 1782-1788.

79.

Garcia de la Torre, J.; Huertas, M. L.; Carrasco, B. Biophys. J. 2000, 78, 719-730.

ACS Paragon Plus Environment

37

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 38 of 49

80.

Garcia de la Torre, J.; Huertas, M. L.; Carrasco, B. J. Magn. Reson. 2000, 147, 138-146.

81.

Garcia de la Torre, J.; Perez Sanchez, H. E.; Ortega, A.; Hernandez, J. G.; Fernandes, M. X.; Diaz, F. G.; Lopez Martinez, M. C. Eur. Biophys. J. 2003, 32, 477-486.

82.

Garcia de la Torre, J.; del Rio Echenique, G.; Ortega, A. J. Phys. Chem. B 2007, 111, 955-961.

83.

Sedeh, R. S. Contributions to the analysis of proteins. Massachusetts Institute of Technology, Cambridge, MA, 2011.

84.

Brooks, B. R.; Janezic, D.; Karplus, M. J. Comput. Chem. 1995, 16, 1522-1542.

85.

Ichiye, T.; Karplus, M. Proteins 1991, 11, 205-217.

86.

Kottalam, J.; Case, D. A. Biopolymers 1990, 29, 1409-1421.

87.

Larson, R. G.; Hu, H.; Smith, D. E.; Chu, S. J. Rheol. 1999, 43, 267-304.

88.

Grassia, P. S.; Hinch, E. J.; Nitsche, L. C. J. Fluid Mech. 1995, 282, 373-403.

89.

Zhou, C.; Yang, Z.; Liu, D. J. Am. Chem. Soc. 2012, 134, 1416-1418.

90.

Han, D. R.; Pal, S.; Nangreave, J.; Deng, Z. T.; Liu, Y.; Yan, H. Science 2011, 332, 342346.

91.

Bai, X. C.; Martin, T. G.; Scheres, S. H.; Dietz, H. Proc. Natl. Acad. Sci. USA 2012, 109, 20012-20017.

92.

Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926-935.

93.

Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kale, L.; Schulten, K. J. Comput. Chem. 2005, 26, 1781-1802.

94.

Allner, O.; Nilsson, L.; Villa, A. J. Chem. Theory Comput. 2012, 8, 1493-1502.

95.

Batcho, P. F.; Case, D. A.; Schlick, T. J. Chem. Phys. 2001, 115, 4003-4018.

ACS Paragon Plus Environment

38

Page 39 of 49

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

96.

Martyna, G. J.; Tobias, D. J.; Klein, M. L. J. Chem. Phys. 1994, 101, 4177-4189.

97.

Feller, S. E.; Zhang, Y. H.; Pastor, R. W.; Brooks, B. R. J. Chem. Phys. 1995, 103, 46134621.

98.

Brünger, A. T., X-PLOR Version 3.1: A System for X-ray Crystallography and NMR. Yale University Press: New Haven, Connecticut, 1992.

99.

Watson, J.; Hays, F. A.; Ho, P. S. Nucleic Acids Res. 2004, 32, 3017-3027.

100.

Lavery, R.; Moakher, M.; Maddocks, J. H.; Petkeviciute, D.; Zakrzewska, K. Nucleic Acids Res. 2009, 37, 5917-5929.

101.

Carrasco, B.; Garcia de la Torre, J. Biophys. J. 1999, 76, 3044-3057.

102.

Flyvbjerg, H.; Petersen, H. G. J. Chem. Phys. 1989, 91, 461-466.

103.

Hernandez-Ainsa, S.; Keyser, U. F. Nanoscale 2014, 6, 14121-14132.

104.

Kim, D. N.; Nguyen, C. T.; Bathe, M. J. Struct. Biol. 2011, 173, 261-270.

ACS Paragon Plus Environment

39

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

Structural Model Beam Model

Hydrodynamical Model

Journal of Chemical Theory and Computation

Page 40 of 49

Bead Model

DNA Nanostructure

Langevin Dynamics

 + Zx + Kx = f (t ) Mx Brownian Dynamics Finite Element Method

Stiffness Matrix (K) Mass Matrix (M)

Rotne Prager Tensor (T)

3N × 3N Friction Tensor  = T−1 Z

6N × 6N Friction Matrix (Z)

ACS Paragon Plus Environment

Zx + Kx = f (t )

x : Basepair Displacements f(t) : Stochastic Loading

B E

1 2 3 4 5 6 7 8 9 10

J3 D

J1 G J4 J5

J2 J6

BM b c BD Fit 1 Journal of Chemical Theory and Computation A H C I

BD

0 10–2

100 t [ns]

90°

102

ACF

Page 41 F of 49

ACF

a

BM BD Fit BD

1 0 10–2

100 t [ns]

90°

90°

ACS Paragon Plus Environment 90°

90°

102

90°

BM

102

BDof Fit49 Journal of Chemical Theory and Page Computation 42 τ [ns]

a

101

100

5 Mode number

Mode comparison

10

c

RMSF [nm]

1

BD

BM

1 2 3 4 5 6 7 b8 910 10 11 125 13 14 15 16

NM EM MD Fit

1

ACS Paragon Plus Environment 5 NM

10

0

0

0

1 NMA

BD b 10 Page Journal 43 of of Chemical 49 Theory and ComputationMD d

d [nm]

a

5

p(d) [a.u.]

6

10

σ(d) [nm]

ACF

d [nm]

1 2 0 4 0 2 3 t [μs] 4 d c 510 6 7 85 9 10 100 300 315 330 0 5 11 t [ns] d [nm] 12 f e 13 BD 0.5 14 BD Fit 151 MD 16 MD Fit 170 ACS Paragon Plus Environment 0 18 –2 0 2 1910 10 10 10 0 n t [ns] 20

20

b

90°

Jtwist

c BD MD

J6

BD44 of 49 Page 1

J1 J5

J4 p(Jtwist) [a.u.]

1 2 3 4 5 6 7 8 9 10

J1 of Chemical Theory J2 J3 Journal and Computation

p(Jtwist) [a.u.]

a

–1 MD

J6

1

J6

ACS Paragon Plus Environment –40 –20 0 Jtwist [°]

20

–40 –20 0 Jtwist [°]

20

–40 –20 0 Jtwist [°]

20

J1

J1

J6

–1

b

c d Journal of Chemical Theory and Computation ACF

Page 45 of 49

1 2 3 4 5 6 7 8 9

90°

90°

90°

90°

90°

90°

ACS Paragon Plus Environment

BM BD Fit BD

1 0 100

102 t [ns]

104

102 t [ns]

104

BM BD Fit BD

1 ACF

a

0 100

a

BM

104

BDofFit49 Journal of Chemical Theory and Page Computation 46

102

100

5

10

Mode number Mode comparison

c

RMSF [nm]

1 4 BD

BM

1 2 3 4 5 6 7 b8 910 10 11 125 13 14 15 16

τ [ns]

τ [ns]

NM

2

ACS Paragon Plus Environment 5 NM

10

0

2 4 NMA

b Page Journal 47 of of Chemical 49 Theory 25 and Computation D

D [nm]

a

20

p(D) [a.u.]

270

30

σ(D) [nm]

ACF

D [nm]

1 2 15 0 90 180 3 t [μs] 4 d c5 25 6 7 820 9 10 15 15 20 25 270 11 269.5 t [μs] D [nm] 12 f e 13 BD 0.05 14 BD Fit 151 16 170 ACS Paragon Plus Environment 0 18 0 2 4 1910 10 10 0 10 n t [ns] 20

20

b Journal of Chemical Theory and Computation 90°

100 80

90°

90°

ACS Paragon Plus Environment

10

Page 48 of 49

60 40

BM

1 2 3 4 5 6 7 8

90°

c

τ [ns]

a

5 10 Mode number

Mode comparison

1

5 5 10 NM

5 10 BM

0

Page Journal 49 ofof 49Chemical Theory and Computation

1 2 3 4

D

D [nm]

25 20

ACS Paragon Plus Environment 15 0

90 180 t [μs]

270

p(D) [a.u.]