Heterogeneous migration routes of DNA triplet repeat slip-outs

It is unclear how the length of a repetitive DNA tract determines the onset and progression of repeat expansion diseases, but the dynamics of secondary DNA structures formed by repeat sequences are believed to play an important role. It was recently shown that three-way DNA junctions containing slip-out hairpins of CAG or CTG repeats and contiguous triplet repeats in the adjacent duplex displayed single-molecule FRET (smFRET) dynamics that were ascribed to both local conformational motions and longer-range branch migration. Here we explore these so-called "mobile" slip-out structures through a detailed kinetic analysis of smFRET trajectories and coarse-grained modeling. Despite the apparent structural simplicity, with six FRET states resolvable, most smFRET states displayed biexponential dwell-time distributions, attributed to structural heterogeneity and overlapping FRET states. Coarse-grained modeling for a (GAC)10 repeat slip-out included trajectories that corresponded to a complete round of branch migration; the structured free energy landscape between slippage events supports the dynamical complexity observed by smFRET. A hairpin slip-out with 40 CAG repeats, which is above the repeat length required for disease in several triplet repeat disorders, displayed smFRET dwell times that were on average double those of 3WJs with 10 repeats. The rate of secondary-structure rearrangement via branch migration, relative to particular DNA processing pathways, may be an important factor in the expansion of triplet repeat expansion diseases.


INTRODUCTION
An increasing number of human disorders are attributed to expansions of short tandem repeats. These repeat expansion diseases (REDs) are monogenic, with toxicity mediated at the level of DNA, RNA, or protein.(1) Although expansion is linked to replication, (2,3) growing evidence also supports a major role for replication-independent somatic instability,(4) mediated by mismatch repair (MMR). (5)(6)(7)(8) Repeat expansion is also associated with polygenic diseases such as autism spectrum disorder. (9,10) In spite of the variation in RED progression and pathogenesis, there are a number of common features, two of which are particularly relevant here. Firstly, the number of repeats correlates with disease severity and with an earlier age of onset. Secondly, the tandem repeat sequences have an ability to form secondary structures such as hairpins, triplexes, and quadruplexes. (11) These secondary structures are key intermediates in several models for repeat expansion involving processes such as replication, transcription, and repair. (1,6,7) The most widely studied secondary structure is the hairpin slip-out, which can form when the repeat has partial self-complementarity. (8,(12)(13)(14)(15)(16)(17)(18) Hu et al. recently reported single-molecule FRET (smFRET) studies of DNA three-way junctions where one arm is a slip-out of trinucleotide repeats (CAG or CTG) of varying length. (19) Expansion of CAG or CTG is the cause of several triplet REDs including Huntington's disease, myotonic dystrophy type 1, and various spinocerebellar ataxias.(1) Hu et al. found that the 3WJs where the triplet repeat sequence extends into the adjacent duplex (so-called mobile 3WJs) can undergo two distinct types of rearrangement: a localized conformational change at the branchpoint (Fig. 1 a) and a longerrange interconversion that was assigned to branchpoint migration (Fig. 1 b). The latter finding was rather surprising, since it had been previously assumed that long slip-outs were unable to undergo branch migration. (1,13) The report from Hu et al. (19) built on earlier work, which recognized the dynamical nature of small repeat slip-outs (20)(21)(22) and isolated hairpins. (23)(24)(25) The significance of such processes is not yet established, but they may play a role in disease progression and could be important as therapeutic targets. (19) It is, therefore, important to fully explore the dynamics of these slip-out DNA structures and to establish the mechanism and properties of branch migration.
In this work, we have performed a detailed analysis of the dynamics of smFRET for 3WJs to find the optimal number of FRET states, and to probe the state-to-state kinetics via the dwell time distributions. We also examined the dependence of slip-out size on the dynamics by extending the range to 40 repeats, which is a repeat number that is within the symptomatic range for a number of CAG-associated triplet REDs. (1) To complement the smFRET experiments, we have also employed coarse-grained modeling based on the oxDNA2 model of DNA. (26) (27) Through a combination of experiment and modeling, we have provided compelling additional evidence for the branchpoint migration model that was proposed recently. (19) The coarse-grained modeling confirms the feasibility of migration and gives insight into the complexity of the conformational landscape, which is supported by both the heterogeneity in the smFRET trajectories and the significantly slower kinetics displayed by the 3WJ with a (CAG) 40 slip-out. We discuss how a competition between branch migration and expansion-prone processing might lead to the characteristic repeat threshold for triplet REDs.

RESULTS smFRET time-trace analysis workflow
The smFRET trajectories of three mobile DNA 3WJs were analyzed, where mobile refers to those 3WJs that are capable of undergoing branch migration between positional isomers ( Fig. 1 b). (19) A 3WJ with either 10 or 40 CAG repeats in the slip-out, henceforth referred to as (CAG) 10 and (CAG) 40 , respectively, and one with 30 CTG repeats in the slip-out, named (CTG) 30 , were studied (see Fig. 2 for all samples studied and Supporting material for details of the oligonucleotides). The (CAG) 40 sample has not been reported previously. The (CAG) 10 40 samples were chosen because this number of repeats is close to the threshold between healthy and symptomatic ranges for several triplet REDs (e.g., 36 CAG for Huntington's disease). (1) The overall workflow of the smFRET analysis is illustrated in Fig. 3. First, the time trajectories were fit across a range of states using the HaMMy smFRET analysis software. (28) It is recommended for HaMMy to fit the time traces to K max þ 2 states, where K max refers to the number of states one expects to observe. (28) Based on the simple model previously proposed, the mobile 3WJs could exist in eight different states (i.e., 2 conformations for each of the four positional isomers). (19) In the previous study, analysis was carried out on the stitched trajectories by constraining HaMMy to K ¼ 6. This value was chosen after visual inspection of the fitted time traces and comparison with confocal single-molecule FRET data for freely diffusing 3WJs. (19) In the analysis detailed below, we extended the state range to K ¼ 6 -10. After processing the traces, the FIGURE 2 Structures studied in this work. See Supporting material for oligonucleotide sequences. The circles represent the donor dye, Alexa488 (green), and the acceptor dye, Cy5 (red). FIGURE 3 smFRET time-trace analysis workflow. smFRET trajectories are first analyzed via HaMMy to obtain discretized FRET traces. The output files are then processed via the MASH-FRET analysis package to produce transition density plots (TDPs) and to calculate the BIC scores to find the optimal number of states K opt . transitions detected by HaMMy were binned and plotted as transition density plots (TDPs) using the MASH-FRET data analysis package, following the procedure described by Hadzic et al. (29) The TDPs were then fit with two-dimensional Gaussians, and the models were optimized by calculating the Bayesian information criterion (BIC). The goal here was to determine whether a consensus could be reached across the models obtained after analysis, allowing us to determine the optimal number of FRET states, K opt , for the trajectories.

Model selection in the stitched and unstitched smFRET trajectories
In the previous study, (19) the FRET trajectories were stitched into one long trace before analysis by HaMMy. Stitching was performed to avoid the need to carry out postprocessing steps to infer an overall model. However, care must be taken with stitching, since artefactual state-to-state transitions at the boundaries between each trace can be introduced. Moreover, molecular heterogeneity may be obscured during analysis as only one trajectory is considered. (30) In order to examine the effects of stitching, both stitched and unstitched trajectories of the (CAG) 10 and (CTG) 30 samples were analyzed according to the workflow above (Fig. 3).
When the stitched trajectories of (CAG) 10 and (CTG) 30 were considered, we found that the optimal number of states, K opt , identified by BIC was highly dependent on the number of states K provided during HaMMy analysis (Figs. 4 a, S1, and S2, and Tables S1 and S2). In general, the value of K opt was observed to be either the same or one state less than K. This can be understood by considering that the algorithm was constrained to analyze only one long trajectory, so that any molecule-to-molecule variation will result in the fitting of additional states. Since these states are well defined in the TDP, the BIC analysis will also return an overfitted model. In short, particular care must be taken when stitched trajectories are employed for maximum likelihood BIC (ML-BIC) estimation. Specifically, a smaller range of possible states may be more appropriate, based on prior knowledge on the dynamics of the system, as was done previously (19).
There was a marked difference between the ML-BIC fitting for the unstitched trajectories (Figs. 4 b, S3, S4, and Tables S3-S6) in comparison to stitched (Fig. 4a). For (CAG) 10 , a consensus across the various results was found for K opt ¼ 6, suggesting that 6 states can be reliably resolved from the data (Fig. 4b). Specifically, we note that K opt was found to be different than 6 only for K ¼ 10; the requirement to fit seven states could be a result of overfitting during HaMMy analysis. However, even for K ¼ 10, the BIC scores for K opt ¼ 6 and 7 were close ( Table S3), suggesting that the data can be reliably  (19).
For the (CTG) 30 experiment, a consensus was also achieved for K opt ¼ 6 (Figs. 4 b). As with (CAG) 10 , a different value was obtained for K ¼ 10, where K opt ¼ 5. This could be ascribed to states with very similar FRET efficiencies that were clustered into one single state after BIC estimation. Similar to the case of (CAG) 10 , the BIC scores for K opt ¼ 5 and 6 were very close (Table S3)  Inspection of the TDP clusters of both samples (Figs. 5 and S1ÀS4) shows that the transition frequencies vary according to the states involved, as supported by the kinetic analysis (see next section). Based on these results, the analysis of individual traces allows assignment of the number of states from complex trajectories in a more unbiased manner, in comparison to those from the stitched trajectories. Specifically, by carrying out model selection via ML-BIC on TDPs, a resistance to overfitting was observed, leading to a consensus across the various discretized state trajectories even for higher values of K.

State-to-state smFRET kinetics
In the previous report from Hu et al., (19) a detailed analysis of the dwell times was not performed. Having established that the data can be described by six FRET states, we next examined the state-to state kinetics of 3WJs with repeat slip-outs. To do this, we use MASH-FRET (29) to generate the normalized cumulative probability distributions (CPDs) of the dwell times for each FRET state. For a Markovian state-to-state process, we would expect these CPDs to fit to a single exponential. To validate the approach, we first analyzed data for a DNA hairpin that was previously shown to obey twostate dynamics. (31) As expected, the TDP (Fig. S5) and dwell times for both states fit to a single exponen-tial ( Fig. S6 and Table S7). We then analyzed data for four static 3WJs with (CAG) 10 slip-outs, a fully complementary 3WJ where the hairpin is composed of (CTG) 5 (CAG) 5 , and mobile 3WJs with (CAG) 10 or (CTG) 30 slip-outs; preliminary data had been reported previously for these samples (see Fig. 2 for schematic representations of the structures). (19) The static 3WJs have almost identical sequences to the mobile (CAG) 10 3WJ but lack the contiguous CAG repeats in the duplex immediately adjacent to the slip-out; therefore, they cannot undergo branch migration. Previously it was shown that the static 3WJs could interconvert between two states, which was attributed to a conformational change at the branchpoint, based on earlier work with fully complementary 3WJs. (32,33) We performed the same TDP-BIC analysis as for the mobile 3WJs, with up to four FRET states. For the static 3WJs, we found two differences in comparison with the previous report. (19) Firstly, although two of the static isomers fitted to two states, there was an improvement in the TDP-BIC for two of the isomers by adding a third state (Fig. S7). Secondly, several of the CPDs required a second component to adequately fit the data (Figs. S8-S11). For the fully complementary 3WJ, the time-trace data could be fitted to two FRET states (Fig. S12), both of which exhibited single exponential dwell times; no improvement was found in fitting to a biexponential ( Fig. S13 and Table S7). The CPDs for the mobile 3WJs FIGURE 5 Transition density plots for mobile 3WJs. Representative TDPs obtained after clustering the stitched (top) and unstitched (bottom) trajectories of (CAG) 10 (left) and (CTG) 30 (right). The TDPs for the trajectories were plotted after running HaMMy with K ¼ 6 and K ¼ 8 for the stitched and unstitched, respectively. The TDP fittings for other values of K are in the Supporting material (Figs. S1-S4). Each colored cluster represents all the transitions assigned from one FRET state (FRET before ) to a different one (FRET after ) centered at the FRET before value. The diagonal data arises from traces with no dynamics or low-amplitude artefactual jumps around the same state. with (CAG) 10 slip-outs ( Fig. S14 and Table S8) and (CTG) 30 slip-outs ( Fig. S15 and Table S8) were qualitatively similar to those of the static samples. The dwell time distributions were best fitted by biexponentials, with the exception of state 5 and state 6, which were monoexponential.
Dependence of FRET dynamics on repeat length: (CAG) 10 versus (CAG) 40 Hu et al. (19) presented evidence for a trend toward longer total dwell time with repeat length. In that study, the maximum repeat length was 30. In order to examine the effect of repeat length using the TDP-BIC approach, we also collected data on a mobile 3WJ with 40 CAG repeats (Fig. 2). Apart from the number of repeats, the hairpin is identical to the sample with 10 CAG repeats. Confocal and total internal reflection fluorescence (TIRF) data for freely diffusing or immobilized (CAG) 40 , respectively, recorded in the same way as in the previous study were qualitatively similar to other mobile CAG 3WJs (Fig. S16). (19) As with the mobile (CAG) 10 and (CTG) 30 samples discussed earlier, the optimal number of FRET states (for K ¼ 6 -10) was found using MASH-FRET to be 6 ( Fig. S17). The six states had mean FRET values of 0.08 (0.05), 0.26 (0.06), 0.42 (0.05), 0.58 (0.05), 0.77 (0.06), and 0.91 (0.05). The only exception was for K ¼ 7, where K opt ¼ 7; however, the BIC score was very similar between K opt 6 and K opt 7 (Table S3).
As above, the dwell time distributions were examined for each of the six FRET states, using the data obtained for K opt ¼ 6. Although some states could be reasonably well fitted to single exponential decays, improvements for all states were found for biexponential fitting ( Fig. S18 and Table S8). In order to compare the effect of repeat length on dwell time, the biexponential fits for mobile (CAG) 10 and (CAG) 40 were used to produce an average dwell time (inset in Fig. 6) In all cases, particularly for state 3 and state 4, the dwell times for the longer repeat are observed to be longer. The average lifetime, taken from the biexponential fit, which represents the area under the dwell time curve, quantifies this (Fig. 6). The average lifetime increased by a factor ranging from 1.2 to 3.38, with a mean increase (SEM) of 2.0 (0.35). Some transitions will be linked to branchpoint folding rather than migration, and it is possible that migration will have a stronger dependence on the size of the slip-out than folding.

Coarse-grained modeling
Coarse-grained modeling was performed using the oxDNA2 coarse-grained model of DNA (26,27) and its implementation in the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) (34) via the CG-DNA package (35). As a model, we used a 3WJ with reverse polarity to the (CAG) 10 that was studied with smFRET ( Fig. 1 b). Fig. 7 shows representative conformations during a simulated slip-out event. Fig. 7 a gives the initial configuration where the nucleotides that hybridize during the slippage event have been marked for clarity in light purple (orange strand 5' -GACGACGAC -3 0 ) and dark purple (blue strand 5' -GTCGTCGTC -3 0 ). We also highlight two nucleotides at the top of the hairpin loop (orange strand 5' -CG -3 0 ) for better visualization of the rearrangements within the hairpin loop. At first, the 3WJ denatures around the branchpoint (Fig. 7 b). This results in two very flexible ssDNA regions (Fig. 7 b), which quickly form stacks of nucleotides (Fig. 7 c), thereby temporarily reducing the free energy again. Unstacking and opening of the hairpin from the stem (Fig. 7 d) is followed by the complete opening of the hairpin loop ( Fig. 7 e) and reclosure from the top (Fig. 7 f), whereby the hairpin undergoes a subtle rearrangement that allows it to take in the slack. Then the hairpin stem closes again (Fig. 7 g), forming a very stable intermediate state. Then it opens again at the stem (Fig. 7 h) and forms again a short ssDNA sequence of stacked nucleotides, visible in Fig. 7 i directly in front of the branchpoint. Unstacking (Fig. 7 j), opening of the hairpin from the stem (Fig. 7 k), complete denaturing of the hairpin loop (Fig. 7 l), closure from the top (Fig. 7 m), and closure at the stem with re-hybridization occur similarly as before (viz., Fig. 7 d-i), completing the branch migration and slippage event in Fig. 7 n. It should be noted that Fig. 7 depicts one typical sequence of events. But although the exact sequence differs slightly in individual runs, the shown conformations are generic and representative for all runs we conducted.
We were able to identify several rate-limiting steps: 1) First of all, basepairs in the hairpin stem close to the branchpoint need to denature in order to create partner nucleotides for the slip-out event.
2) The nucleotides in the arm on the right in Fig. 7, which provide the partner nucleotides in the hairpin loop and receiving arm on the left, need to denature and create a bulge of nucleotides in the process that is ready to move up into the hairpin loop. 3) New basepairs need to form in the receiving arm on the left in Fig. 7 through hydrogen bonding with vacant, complementary nucleotides. 4) Finally, the hairpin loop needs to open up and rearrange to take in the ssDNA bulge and reduce stress along the duplex axes. Particularly the last event can take longer as the ssDNA bulge shows a tendency to form an array of stacked nucleotides, which is energetically favorable but prevents the bulge from moving up into the hairpin loop. When this occurs, the 3WJ appears quite distorted and does not resemble a clear T-or Y-shape. This means that nucleotides that are partners in post-slippage basepairs can be held back at a distance from each other, making the hybridization rather unlikely.
The slippage event seems to occur in two separate, but more or less identical steps: one to two basepairs are formed rapidly in the receiving arm, but further formation is stalled until the aforementioned steps one through four have taken place. When this has happened another two to three basepairs form, bringing the total number of slipped-out nucleotides up to four to five basepairs. Then, steps one through four have to take place again for the remaining nucleotides to undergo slippage, bringing the total number of slipped basepairs up to the final nine in the receiving arm.
It is interesting to analyze the free energy landscape that determines the thermodynamic aspects of the branchpoint migration. Fig. 8 shows the averaged free energy from seven independent runs versus the two collective variables, namely the number of basepairs that need to open up (labeled "denatured bps") and form (labeled "hybridized bps") during the slippage; the letters correspond to the snapshots shown in Fig. 7. Both the initial state (a) in the bottom left corner and the final state (n) in the top right corner are characterized by noticeable minima with values around -32 k B T.
From the heatmap, it becomes obvious there are two regions around (g) and (h), which are relatively easily accessible from the initial state with free energy values around -27 and -23 k B T, respectively. These regions represent partially slipped states that correspond to the conformations shown in Fig. 7 g and h with around two to four slipped basepairs. In these states, some new basepairs have already formed in the receiving arm on the left, but the hairpin loop opens only partly at the stem, and the nucleotides in the ssDNA bulge undergo repeated stacking events that block further slippage. However, the state (g) features a hairpin that has taken in the slack and is therefore very stable and relaxed. The state (h) has progressed slightly further in the two reaction coordinates and represents an intermediate conformation before the big slip event. Similar intrastrand interactions have been observed previously in trinucleotide-repeat DNA hairpins, (23) where it was concluded that the repetitive, high mismatch content and self-complementarity gave rise to misfolding and metastable states.
The two basins formed by the initial state (a), including its extension to the states (g) and (h), and the final state (n) representing a full slip-out are separated by another region with free energy values around -17 k B T, where the intermediate states (k) and (l) are situated. This leads to an effective barrier of around 15 k B T between the initial state (a) and the final state (n), a value that is comparable with that found in another single-molecule fluorescence study on DNA hairpin dynamics.(36) Using Kramers' reaction rate theory for diffusive barrier crossing, the transition rate between the two positional isomers can be expressed with k f and k b as forward and backward rate, respectively, DG z as the height of the free energy barrier between the two positional isomers, and k 0 as preexponential factor. Note we assume here that forward and backward rate are the same, which is justified as FIGURE 6 Dependence of dwell time on repeat length. Biexponential dwell-time fits for (CAG) 10 (red) and (CAG) 40 (black). Panels (a)-(f) correspond to FRET states 1-6, respectively. The average dwell times from the biexponential fits are shown (inset). the free energy values in the respective minima related to the two positional isomers in Fig. 8 are virtually identical. Although the exact pre-exponential factor is a piori not known, it is plausible to assume values of the order k 0 zOð10 6 Þs À 1 as these have been typically found in studies of hairpin dynamics. (36,37) Together with our value of the free energy barrier, this leads to characteristic transition rates of Oð1Þs À 1 , which are very compatible with our experimental observations. It should be noted that we opted here for the simplest bias that allowed us to observe full slippage events between the initial and final state. A more detailed study and statistical analysis of the slippage events with different initial configurations and longer hairpin loops and different triplet repeat sequences will be conducted in the future.

DISCUSSION
We recently showed that mobile 3WJs with hairpins formed from CTG or CAG repeats exhibit two distinct types of dynamics, involving two-state rearrangement at the local branchpoint and interconversion of positional isomers. (19) In this work we have examined these 3WJ dynamics in detail using smFRET and coarse-grained modeling.
A drawback of smFRET is that the signals are inherently noisy, which has led to the development of specialized methods of analysis.(30) One of the most popular approaches is hidden Markov Modeling (HMM), which is well-suited to analyzing noisy signal trajectories. (38) HMM formulations seek to find the hidden variables within observed time series and identify the most probable path through them.(39) As a result, they allow identification of FRET states that may be hidden by noise in the smFRET signal.(30) Dedicated HMM analysis software for smFRET data has been freely available for over 15 years, (28) and this continues to be an active area of research. (40)(41)(42)(43)(44)(45)(46)(47)(48) Previously, Hu et al. (19) employed the software HaMMy developed by McKinney et al. (28) to identify the number of FRET states for DNA 3WJs containing trinucleotide repeats. (19) To estimate the best fitting path for a given signal trajectory, HaMMy uses a maximum likelihood (ML) approach, which aims to find the parameters that maximize the probability of observing a given set of data. McKinney et al. (28) also showed how to use the BIC(49) from the discretized HaMMy traces to infer the most optimal number of states from the HaMMy data by selecting the model that returns the lowest value of the BIC. BIC estimation further improves on direct ML inference for model selection, as it introduces a penalty term for complexity. By avoiding model overfitting, BIC inference can be employed as a powerful tool to obtain the simplest model that best fits the data. (50) We analyzed individual unstitched trajectories via ML-BIC on TDPs using the MASH-FRET package, which demonstrated a resistance to overfitting and provided a consensus about the number of FRET states. This was not the case for the stitched trajectory, suggesting that particular care needs to be taken when analyzing such complex datasets by this method. Nevertheless, the overall conclusion of the previous work in terms of the conformational and migrational dynamics holds. Rather than finding additional states with the improved time-trace analysis approach employed here, we have resolved the same number of states.
By comparing the FRET states for the static and mobile (CAG) 10 samples, we can see that all of the static FRET states are present in the time-trace data for the mobile 3WJs. In other words, the mobile 3WJs can interconvert between all four positional isomers. Our simple model of two-state conformational dynamics for each of four positional isomers would produce up to eight FRET states. We therefore might expect one or two of the FRET states to be triply or doubly degenerate, respectively. However, heterogeneity seems to be the norm, rather than the exception, and it is also evident in the static 3WJ samples. Although, hidden smFRET states have been documented previously, our system (with six degenerate FRET states) is more complex than most systems that have been reported to date. (44)(45)(46)(47) As a result of the heterogeneity in FRET states, we are unable to unambiguously assign the individual FRET states to particular 3WJ structures. However, it is clear that almost all transitions are possible. This means that, within the timescale of the experiment, there must be transitions from one positional isomer to its neighbor, but also as far away as the next nearest neighbor.
To support the smFRET experiments, we also adopted a coarse-grained method to study 3WJ dynamics. We studied a 3WJ with a (GAC) 10 slip-out, using oxDNA2 in combination with the well-tempered metadynamics method.(51) The oxDNA2 model applies a top-down approach to coarse-grained modeling and represents each nucleotide as a single rigid body with effective bonded and pair interactions between the nucleotides. oxDNA2 is able to describe the thermodynamics of duplex formation very accurately and provides a good average representation of the structure of DNA with major and minor groove and the mechanics of both single-and double-stranded DNA and its assemblies. These thermodynamic and structural aspects are very important for investigating structures such as 3WJs, where the conformations deviate significantly from the canonical B-DNA form.
The coarse-grained modeling indicates that a simple bulge-loop model is not appropriate for the hairpins studied here, which are long enough to form stable hairpins. Instead, more complex, metastable, secondary structures prevail, which are formed by intermediate states between the two positional isomers and introduce the requirement of coordinated structural rearrangements for slippage events to occur. The source of both the branch migration and the dynamic heterogeneity is likely to be the local branchpoint and hairpin structure. It is worth noting that dynamic heterogeneity has been reported for the related Holliday junction, (52) where it was attributed to Mg 2þ -induced folding; the various topologies revealed by modeling here (Fig. 7) are consistent with this. Importantly, this dynamic heterogeneity is absent from both a simple hairpin made from ssDNA and a fully complementary 3WJ hairpin. In other words, the key to the heterogeneity appears to be the mismatches in the hairpin and near the branchpoint. These mismatches presumably allow more conformational flexibility, in agreement with ensemble NMR experiments of static 3WJs. (53) This flexibility could facilitate the branch migration in the mobile 3WJs. (19) For longer slip-outs, additional heterogeneity could result from various multiloop configurations or even tertiary structures. (54) The possibility of repeat hairpins migrating was discussed previously in the context of spontaneously formed slip-out hairpins, though there was no experimental evidence at that time. (55,56) For very short repeat slip-outs, which are essentially unpaired loops, (19) it has also been proposed that branch migration occurs via soliton-like movement of small bulge loops (Fig. 9 a). (20)(21)(22) Related to this, isolated hairpins formed from longer repeat sequences have been shown to undergo dynamic rearrangements, which have been attributed to small loop movement. (24,25) These isolated hairpins could form via a nick in a single strand of duplex DNA, a potential source of instability. (25) Although the dynamics observed in the 3WJs is clearly related to those of the small bulge loops and isolated hairpins above, placing a longer repeat hairpin in the context of duplex DNA is rather different, as was noted previously. (55) Based on earlier work showing that a single base mismatch is a major barrier to Holliday junction branch migration, (57) it was assumed that 3WJ branch migration does not to occur.(1) However, up to four contiguous internal mismatches were required to significantly slow strand displacement by an invading single strand in dsDNA.(57) Indeed, a single strand of DNA with a run of three mismatched bases was capable of displacing an otherwise complementary strand from a 3WJ, albeit at a reduced rate.(58) Therefore, mismatches per se are not a major obstacle. We have also presented evidence that the 3WJ dynamics increase with repeat length. Since stable hairpins can migrate (Fig. 9 b), this opens up the possibility that the likelihood of harmful expansion might relate to the relative rate of branch migration and expansion-prone processing. In support of this, it is known that the length of the pure CAG repeat length in Huntington's disease is a better predictor for the age of disease onset than the length of the translated polyglutamine (from CAG with CAA interrupts).(59) The model of branch migration presented earlier (19) and supported here relies on the self-complementarity of the repeat, and it would be impeded by such interrupts.
Migration could involve movement of the slip-out to a less or more error-prone site, or it could result in merging of hairpin loops on opposite strands to recover the duplex (Fig. 9 b). It is possible that similar migrational dynamics also occur for other secondary structures (e.g., triplexes or quadruplexes). At some length of repeat, there could be a tipping point where secondary structure relocation or reversal slows and expansion-prone processing takes over. With this model, the repeat length threshold for a specific RED would be set by the relative rate of secondarystructure dynamics and the specific downstream cause of pathogenesis.
Although this is plausible, evidence from replicationindependent somatic expansion suggests that DNA expands and contracts by only a few repeats at a time, which points toward smaller bulges as the key secondary structure that are incorrectly processed. (4,5,(60)(61)(62) However, if it is a smaller loop that is processed, why is there such a strong dependence on repeat length, and why does there appear to be such an important role for longer slip-outs? Small loops are essentially fast moving bulges and would be expected to migrate relatively freely. (20)(21)(22) Although longer repeats would require small loops on opposite strands to traverse greater distances before merging, this does not seem to explain the rather dramatic dependence on repeat length. Furthermore, CAG and CTG repeats certainly form stable hairpins, and these have been implicated in various expansion mechanisms. (2,3,7,63,64) Hairpins have also been observed in vivo, where somatic expansion correlated with levels of slipped DNA. (65) It is possible that there is more than one mechanism for expansion, depending on the length of the slip-out, the repeat sequence, the sequence context, the tissue, etc. However, it is also feasible that there is a single model that might reconcile the different viewpoints. Kaplan et al. postulated that a universal model for trinucleotide expansion diseases might exist, which would involve an increasing probability of somatic expansion as a function of repeat length.(66) More recently, a dependence of expansion on repeat length has been demonstrated. (67) We propose here a mechanism of somatic expansion that involves both large and small loops (Fig. 9 c). Although a short repeat sequence may preferentially form small loops, these might be expected to migrate quickly and recover the stable duplex ( Fig. 9 a). As the number of repeats increases, it might be expected that larger and larger hairpins can form (Fig. 9 b). Although not processed by MMR, longer hairpins are recognized by the mismatch repair machinery. (68) We speculate that the large slip-outs may have a role in blocking the migration and merging of the small loops. As the larger hairpins grow, they might become increasingly stable and less prone to migrating and re-annealing, giving time for the cellular machinery to process the small loops, ultimately leading to expansion.
Annealing of samples for branched DNA was carried out in buffer (20 mM Tris (Sigma-Aldrich) and 50 mM NaCl (Fluka), pH 7.5) with the ratio of donor strand to acceptor strands at 1:2. Samples were heated to 90 C in a water bath and left to cool down slowly overnight. For TIRF measurements (see below), the buffer contained 20 mM Tris-HCl, 10 mM NaCl (pH 7.8) with 6% glucose (w/w), 2 mg/mL glucose oxidase, 0.08 mg/mL glucose catalase, and 1 mM Trolox to reduce the rate of blinking and photobleaching of the dyes. The measurement buffer contained 1 mM MgCl 2 (Fluka). All measurements were recorded at 20 C 5 1 C.
Objective-type total internal reflection fluorescence microscopy data TIRF experiments were performed, and text files containing the TIRF trajectories were generated as reported previously. (19) Briefly, the TIRF signals were recorded with an exposure time of 50 ms and were processed using TwoTone (v3.1).(69) Photobleaching removal and generation of donor and acceptor intensity time traces was achieved using a custom MATLAB script.

Evaluation of the number of states
Several hundred single-molecule FRET trajectories for each sample were analyzed to evaluate the optimal number of states, K opt . In the previous study these traces were stitched. Here we used the same traces either stitched or unstitched. Trace processing was performed using the HMM software HaMMy (v4.0). (28) In order to test the robustness of the analysis method, the FRET efficiency trajectories were fit across a range of states, specifically from 6 to 10. HaMMy required a reasonably long time to analyze the stitched trajectories, especially when the value of K exceeded 7 or 8. For example, when K was chosen to be 10, HaMMy took up to 3 hours to provide the results on a laptop with operating system Windows 10 Pro and an Intel Core i5-6200 processor. This is attributed to the size of the stitched trajectories (up to 45,000 data points), which is close to the limit of datapoints that can be analyzed via HaMMy. (28) Overall, when the unstitched trajectories were analyzed, the computational time reduced significantly across the whole dataset, typically taking half of the time required for the stitched ones.
After running, the program outputs a file containing the idealized FRET trajectories with the intensity values of each transition and the idealized FRET state. The output files were imported in MASH-FRET 1.3.1 in the "Transition analysis" tab to build TDPs and obtain the optimal number of states K opt . (29) In the TDPs, the transitions from one FRET state to another were counted and evenly spaced out in 100 bins between 0 and 1. The TDP was further convoluted to obtain a smooth plot (s 2 ¼ 0.0005). (28,70) To perform model selection, the TDPs were fitted to a mixture of K Â (K -1) isotropic 2D Gaussians. FRET values were derived from the Gaussian means and the associated errors (in parentheses) from the average Gaussian sample standard deviations. The value of K opt was inferred via ML-BIC optimization of models comprising K ¼ 1 to 10 states with 10 initializations. Overall, a total of 10 TDPs (five for each system) were fitted, and the results were compared to determine the most likely model.
Finally, the analyses were carried out using MATLAB R2020b on a laptop with operating system Windows 10 Pro and an Intel Core i5-6200 processor. Trace analysis was performed using HaMMy v4.0(28) and MASH-FRET 1.3.1. (29) Dwell time cumulative probability distributions CPDs were plotted using the cumulative dwell times histograms provided by MASH-FRET after TDP clustering. (29) The cumulative counts were imported to Origin v2020b and plotted as histograms. The rate coefficients were then recovered by fitting the cumulative counts to either single exponential decay (ExpDec1) or biexponential decay (ExpDec2) using nonlinear curve fitting (Levenberg-Marquardt).
Average dwell times (t av. ) for biexponential fits were calculated by weighting the dwell times (t) by the amplitudes (A): t av. ¼

Coarse-grained modeling
We made use of the well-tempered metadynamics method for enhanced sampling of rare events through the PLUMED library, (71,72) which is interfaced with the LAMMPS code through its PLUMED package. In particular, we used the LAMMPS version from July 2, 2021, and PLUMED version 2.7.1.
A bias in two independent collective variables measuring the approximate number of basepairs in specific locations was applied. One collective variable measured the number of basepairs that were to form in the hairpin loop and in the arm receiving the slipped-out sequence during a slippage event. The other collective variable consisted of the number of basepairs that would denature in the hairpin loop and in the other arm. Both collective variables were defined through a PLUMED contact map function, which measures the proximity of specific entities through smooth, but sharply decaying step functions. This entails that the number of basepairs is not integer valued and only approximate to the number of actual hydrogen bonds. Nevertheless, it gives a well-defined and clear control handle for the biasing algorithm and characterization of the individual states.
In the initial state the value of both collective variables is minimal and increases during the slippage event. The metadynamics algorithm biases the free energy landscape by adding gradually small Gaussian bias potentials, therefore forcing the system to move away from the energetically more favorable states and towards other, less favorable states that would normally not be accessible during realistic simulation times. Examples of these states include the intermediate states in Fig. 7 k and l between the initial and fully slipped-out final state. By recording the applied bias over time, running the metadynamics algorithm until convergence and tallying the recorded bias afterward, it is possible to extract an accurate picture of the free energy landscape of the unbiased system as the inverse sum of all FIGURE 9 Branch migration of small unpaired loops and stable hairpins. (a) Small loops can move freely across a duplex composed of fullypaired repeats. (b) Large hairpins can migrate more slowly than small loops. (c) Large hairpins might act as roadblocks for small loops, preventing loops on opposite strands from merging and increasing the time available for processing of small loops (e.g., by mismatch repair). applied bias potentials. The fact that our free energy values are fully compatible with the experimentally observed transition rate of Oð1Þs À 1 provides strong evidence that the detected values are indeed realistic. All simulations were run at a temperature T ¼ 300 K using a metadynamics bias factor g ¼ 8.0 and different random seeds in a Langevin thermostat to sample independent trajectories. The electrostatic interaction is modelled by a Debye-H€ uckel approximation with implicit ions. Although the salt concentration of the presented runs is 0.2 M, the oxDNA2 model has been parameterized to minimize the difference between the observed melting temperatures T m in oxDNA2 and the SantaLucia nearest-neighbor model for a range of salt concentrations between 0.1 and 0.5 M. This leads to a difference DT m of around 0.25 K. (27) The metadynamics algorithm was deemed as having converged when frequent transitions between the two basins of the initial and final state were observed and the height of the applied Gaussian bias potentials (as controlled through the algorithm) had significantly tapered off. For more information on the metadynamics method, we refer to a recent review. (73)

Data availability
The data that support the findings of this study are available from the University of Glasgow data repository Enlighten at http://dx.doi.org/ 10.5525/gla.researchdata.1337.

DECLARATION OF INTERESTS
The authors declare no competing interests.