Coarse-Grained Modeling and Simulations of Thermoresponsive …
47
other structural parameters as described in References [19, 20, 66, 67]. We direct the
reader to References [19, 20, 66, 67] to understand the variations one would need to
make to the above generic simulation protocol and analyses methods depending on
the system of interest and the question one wishes to answer.
2.4 Key Results
ONA CG Model Parameterization for DNA: In this section, we describe how we
parameterize the above CG model to represent DNA and capture experimentally
observed DNA melting trends. It is known that spacing between nucleobases in a
single-stranded DNA is 0.34 nm and for a DNA in duplex state, it is 0.63 nm [78], thus,
we set the bond length between BB beads, r
BB
0 = 0.84σ (= 0.5 nm) to approximately
represent the average of base spacing distance in single-stranded and double-stranded
DNA. To capture the semi-flexibility of DNA in our CG model, we tune the force
constant of harmonic bond angle potential between three consecutive BB beads (BBBB-BB) to be k
BB
angle = 30ε/rad
2 which leads to a persistence length of 8.4nm for an 8mer single-stranded DNA and 100% G-C content at 307 K. Stacking interactions are
captured phenomenologically via harmonic dihedral potentials with force constants
k
i j
stacking chosen to distinguish the strongest (purine) stacking interactions (GG, GA
and AA) and weakest stacking (pyrimidine) interactions (CC, CT and TT). k
i j
stacking
is set to 9ε for G-G, G-A and A-A and k
i j
stacking is set to 3ε for C–C, C-T and T-T. For
stacking interactions between purine and pyrimidine bases (AC, GT and GC), the
k
i j
stacking is set to 6ε. We set k
i j
stacking = k
ji
stacking as our CG model does not distinguish
between 3
and 5
directions. Lastly, the BB bead of DNA is negatively charged with a
valency of -1 in our CG model to represent the negatively charged phosphate groups.
The non-bonded interactions, i.e., ε
HB
G−C and ε
HB
A−T , are parameterized to ensure they
can reproduce melting curves for DNA for various G-C contents obtained from NN
model [79]. We choose the empirical NN model as compared to experimental data
because the majority of the published experimental literature for short DNA duplexes
are available only at low DNA concentration (< 1μM) and high salt concentration
(> 50mM) both of which are computationally intensive. To optimize the non-bonded
interactions for DNA, we perform simulations at 6μM DNA duplex concentration
and 1mM salt concentration. We systematically vary ε
HB
G−C and ε
HB
A−T to match melting
curves obtained from our simulations to that of empirical NN models. To obtain the
values of ε
HB
G−C and ε
HB
A−T , melting curves from NN model for 100% G-C content and
0% G-C content DNA are used, respectively. In particular, the following sequences
are used: d(GGGGGGGG), d(GCGCGCGC), d(AAAAAAAA) and d(ATATATAT).
Here, d(XYZ) represents a duplex of strand XYZ and its complementary strand, with
“d” representing the duplex. The best fit of melting curves between our simulation
results and NN model is obtained using the values of ε
HB
G−C = 61ε and ε
HB
A−T =
50ε where ε = 0.1kcal/mol. Figure 4 shows the comparison of our simulation
results along with NN melting curves from three sources. This confirms that our CG
47
other structural parameters as described in References [19, 20, 66, 67]. We direct the
reader to References [19, 20, 66, 67] to understand the variations one would need to
make to the above generic simulation protocol and analyses methods depending on
the system of interest and the question one wishes to answer.
2.4 Key Results
ONA CG Model Parameterization for DNA: In this section, we describe how we
parameterize the above CG model to represent DNA and capture experimentally
observed DNA melting trends. It is known that spacing between nucleobases in a
single-stranded DNA is 0.34 nm and for a DNA in duplex state, it is 0.63 nm [78], thus,
we set the bond length between BB beads, r
BB
0 = 0.84σ (= 0.5 nm) to approximately
represent the average of base spacing distance in single-stranded and double-stranded
DNA. To capture the semi-flexibility of DNA in our CG model, we tune the force
constant of harmonic bond angle potential between three consecutive BB beads (BBBB-BB) to be k
BB
angle = 30ε/rad
2 which leads to a persistence length of 8.4nm for an 8mer single-stranded DNA and 100% G-C content at 307 K. Stacking interactions are
captured phenomenologically via harmonic dihedral potentials with force constants
k
i j
stacking chosen to distinguish the strongest (purine) stacking interactions (GG, GA
and AA) and weakest stacking (pyrimidine) interactions (CC, CT and TT). k
i j
stacking
is set to 9ε for G-G, G-A and A-A and k
i j
stacking is set to 3ε for C–C, C-T and T-T. For
stacking interactions between purine and pyrimidine bases (AC, GT and GC), the
k
i j
stacking is set to 6ε. We set k
i j
stacking = k
ji
stacking as our CG model does not distinguish
between 3
and 5
directions. Lastly, the BB bead of DNA is negatively charged with a
valency of -1 in our CG model to represent the negatively charged phosphate groups.
The non-bonded interactions, i.e., ε
HB
G−C and ε
HB
A−T , are parameterized to ensure they
can reproduce melting curves for DNA for various G-C contents obtained from NN
model [79]. We choose the empirical NN model as compared to experimental data
because the majority of the published experimental literature for short DNA duplexes
are available only at low DNA concentration (< 1μM) and high salt concentration
(> 50mM) both of which are computationally intensive. To optimize the non-bonded
interactions for DNA, we perform simulations at 6μM DNA duplex concentration
and 1mM salt concentration. We systematically vary ε
HB
G−C and ε
HB
A−T to match melting
curves obtained from our simulations to that of empirical NN models. To obtain the
values of ε
HB
G−C and ε
HB
A−T , melting curves from NN model for 100% G-C content and
0% G-C content DNA are used, respectively. In particular, the following sequences
are used: d(GGGGGGGG), d(GCGCGCGC), d(AAAAAAAA) and d(ATATATAT).
Here, d(XYZ) represents a duplex of strand XYZ and its complementary strand, with
“d” representing the duplex. The best fit of melting curves between our simulation
results and NN model is obtained using the values of ε
HB
G−C = 61ε and ε
HB
A−T =
50ε where ε = 0.1kcal/mol. Figure 4 shows the comparison of our simulation
results along with NN melting curves from three sources. This confirms that our CG