# ChromoShake: a chromosome dynamics simulator reveals that chromatin loops stiffen centromeric chromatin

## Abstract

ChromoShake is a three-dimensional simulator designed to find the thermodynamically favored states for given chromosome geometries. The simulator has been applied to a geometric model based on experimentally determined positions and fluctuations of DNA and the distribution of cohesin and condensin in the budding yeast centromere. Simulations of chromatin in differing initial configurations reveal novel principles for understanding the structure and function of a eukaryotic centromere. The entropic position of DNA loops mirrors their experimental position, consistent with their radial displacement from the spindle axis. The barrel-like distribution of cohesin complexes surrounding the central spindle in metaphase is a consequence of the size of the DNA loops within the pericentromere to which cohesin is bound. Linkage between DNA loops of different centromeres is requisite to recapitulate experimentally determined correlations in DNA motion. The consequences of radial loops and cohesin and condensin binding are to stiffen the DNA along the spindle axis, imparting an active function to the centromere in mitosis.

## INTRODUCTION

We have developed a computational polymer dynamics simulator named ChromoShake to aid in understanding the physical and behavioral properties of chromatin in the presence of thermal motion. Unlike molecular-scale simulators, ChromoShake can simulate the dynamics of megabases of chromatin on a timescale sufficient to reach an equilibrium state in a matter of days. Chromatin dynamics can be directly compared with time-lapse microscopy in which chromatin is visualized through DNA-binding green fluorescent protein (GFP) fusion proteins and their target DNA sequences. Moreover, ChromoShake was designed to study the thermodynamics of the class of proteins known as structural maintenance of chromosomes (SMC). These are ring complexes that are involved in chromosome condensation (condensin) and sister chromatid cohesion (cohesin). Using ChromoShake, the three-dimensional (3D) distribution of these proteins can be quantitatively compared with their in vivo distribution (Yeh *et al.*, 2008; Stephens *et al.*, 2011). ChromoShake is designed 1) to increase our intuition on how thermal motion affects the dynamics and overall organization of chromatin of a given geometry and 2) to serve as a physically accurate model for testing chromatin structure. Although ChromoShake can be used to generate a variety of 3D chromatin structures (see *Materials and Methods*), we limit our use of ChromoShake in this study to the budding yeast pericentromere.

We constructed a 3D model of the budding yeast centromere using experimentally determined geometric distributions of fluorescently labeled cohesin (SMC3-GFP), condensin (SMC4-GFP), and lactose operon/tetracycline operon (lacO/tetO) arrays within the pericentric region. Cohesin loading is regulated and largely coupled to replication, thus ensuring that cohesion is established between sister chromatids. Cohesin and condensin are threefold enriched in a 30- to 50-kb of DNA around each centromere, (Megee *et al.*, 1999; Tanaka *et al.*, 1999; Glynn *et al.*, 2004; Weber *et al.*, 2004; D’Ambrosio *et al.*, 2008; Hu *et al.*, 2011). Cohesin function in the pericentromere remains enigmatic. Because sister centromeres are separated by 800 nm, cohesin does not tether sister centromere DNAs to one another in mitosis (Pearson *et al.*, 2001). It has been proposed that individual cohesin complexes link different chromatin strands and different loops of the same strand (Stephens *et al.*, 2013). Condensin is found along the spindle axis and appears variously as spots or linear arrays between the spindle poles in mitosis (Stephens *et al.*, 2011). The position of condensin is dependent on tRNA binding factors that function to cross-link pericentromeres from different chromosomes (Snider *et al.*, 2014). Pericentric chromatin labeled with lacO/LacI-GFP at multiple locations has been found to be radially displaced from the spindle axis when GFP signals appear as foci, but falls closer to the spindle axis when signals appear stretched (Stephens *et al.*, 2011; Haase *et al.*, 2012).

There are several models for the chromatin organization within the centromere. A model in which radial loops emanate from a primary axis between the spindle poles of the mitotic apparatus is most consistent with experimental findings (Yeh *et al.*, 2008; Stephens *et al.*, 2011). Recent work has suggested that crowded radial loops stemming from a primary axis build tension along the axis (Lawrimore *et al.*, 2015). This is a severe departure from the established view that the centromere acts as a passive spring, simply stretched and compressed by microtubule dynamics over which it exerts little control. To explore this hypothesis, we have used molecular dynamics based on in vivo structures to determine the mechanical properties of the centromere. ChromoShake provides the opportunity to quantitatively assess how well a geometrical model that fluctuates to the favored entropic state compares with experimental findings. Moreover, this approach provides the opportunity to differentiate active, directed events in the pericentromere from passive, entropic mechanisms that together contribute to the function and morphology of the centromere in live cells.

## RESULTS

### Thermodynamics collapses radial subloops

The pericentromere is the region of DNA surrounding the centromere on each of 16 chromosomes. The 16 centromeres are replicated before mitosis, become clustered, and biorient relative to the spindle poles (Figure 1A, red disks) of the mitotic spindle apparatus. In our model, the centromere DNA lies at the apex of the primary loop, where the kinetochore is built and promotes microtubule attachment (Figure 1A, kinetochores are white end beads, green rods are microtubules). The primary centromeric axis is parallel to the microtubule spindle axis. To mimic the radial displacement of lacO/tetO arrays observed in vivo (Stephens *et al.*, 2011; Haase *et al.*, 2012), each chromatid is arranged as a loop with four subloops radiating from the primary loop (proximal to the spindle axis; Figure 1A, colored strands). Chromosome arms (Figure 1A, elongated strands) beyond the pericentromere are not modeled in simulations (Figure 1, B–E). Each pericentric chromatid in the model is 50 kb (Figure 1, B–E). The 32 chromatid strands represent the replicated 16 chromosomes. The total amount of simulated DNA in the ensemble of 32 chromatids is ∼1.6 Mb (Figure 1, B–E). The simulated DNA occupies a cylindrical geometry ∼778 nm in diameter and 800 nm in length (Figure 1, B–E).

To understand how thermal motion affects the structure of the centromere, we built a 3D model based on experimental observations of pericentric chromatin and assessed the range of configurations dictated by thermal forces. The computational model is a way to build intuition about a complex structure that would otherwise be difficult to predict. A robust model of chromatin dynamics can be captured using a simple bead-spring model. In our model, the beads repel each other and are not allowed to pass through one another. This is known as excluded volume and simulates the volume around each chain that is inaccessible to other chains in the system. Hinge-like forces between each bead are parameterized to give the strands the same bending rigidity of DNA, that is, a persistence length of 50 nm (Bloom, 2008). In the model, cohesin is simulated as a 51-nm-diameter ring (Figure 1, white rings) formed by 16 linked beads that encompasses two neighboring chromatid strands. These rings have the same rigidity as the simulated chromatin and can migrate along the chromatin strands. Thermal forces collapse the rings to an average diameter of 44.6, 44.9, and 47.0 nm for cohesin not encompassing chromatin, cohesin encompassing a single chromatin strand, and cohesin encompassing four chromatin strands in the centromere model, respectively (Supplemental Table S1). Condensin is modeled as a spring (not shown in Figure 1) that connects two beads on the same chromatin strand to form radial subloops. Condensin is localized along the spindle axis in fluorescent images (Stephens *et al.*, 2011); therefore we confined our condensin springs proximal to the spindle axis at the base of the radial subloops. The beads at the ends of the primary loops (Figure 1, white beads) are five orders of magnitude more massive than the other beads to “pin” the model in space and mimic attachment of the centromere to microtubules. Microtubule dynamics were not included, as their contribution to fluctuations of pericentromeric DNA is nominal (Lawrimore *et al.*, 2015). We loaded our centromere model into ChromoShake to simulate the thermal motion of beads and bead-ring complexes. The model was allowed to run for ∼0.05 s of simulation time and then rendered (Figure 1, D and E). These images show how chromatin loops crumple toward the spindle axis upon perturbation by thermal forces.

### Radial subloops replicate the position and the motion of chromatin observed in experimental samples

Previous studies have measured chromatin position and fluctuations in vivo using tandem-repeat sequences of the lacO or tetO integrated into chromosomes in cells expressing LacI-GFP or TetR-GFP fusion proteins (Straight *et al.*, 1996; Michaelis *et al.*, 1997). The probability distribution of fluorescently labeled chromatin reveals that the pericentromere DNA is, on average, radially displaced from the spindle axis (Anderson *et al.*, 2009; Stephens *et al.*, 2011). Because ChromoShake only simulates thermal forces, we compared our simulated motion data with experimental data in which ATP was depleted through treatment of cells with sodium azide and deoxyglucose. We performed time-lapse imaging of pericentric chromatin labeled with a 10-kb (256 repeat) lacO/LacI-GFP array centered 6.8 kb from CEN15. We imaged cells for 10 min every 30 s in media containing sodium azide and deoxyglucose. To mimic the 6.8-kb lacO array in the centromere model, we tracked the mean position of a 73-bead array (∼10 kb) centered 49 beads (∼6.8 kb) from each of the centromeres (Figure 1, centromeres are white beads). Given the symmetry of our centromere model (Figure 1), we individually tracked the mean position of all 64 simulated 6.8-kb arrays. We found the 6.8-kb array has a radial displacement (mean ± SD) of 134 ± 74 nm in ATP-depleted cells and 140 ± 51 nm in the simulations (Figure 2, A and B).

The fluctuations of the tetO/lacO arrays have been measured in live cells using single-particle tracking and mean-squared displacement (MSD) analysis, allowing us to compare the fluctuations of our model directly with live-cell measurements. The motion of tetO/lacO arrays is confined, and the extent of confinement is reflected by the plateau value of the MSD curve (Weber *et al.*, 2012; Verdaasdonk *et al.*, 2013; Chacon *et al.*, 2014). The plateau value of a 6.8-kb array in ATP-depleted cells was published as 5600 ± 3200 nm^{2} (Lawrimore *et al.*, 2015), whereas the mean ± SD of the MSD plateau value for the simulated 6.8-kb array was 4000 ± 110 nm^{2} (Figure 2C). Thus our simulation recapitulates the experimental distribution and fluctuations of chromatin in the absence of ATP-dependent processes.

### Cohesin distribution is determined by chromatin subloop size

Cohesin is radially displaced from the spindle axis (Yeh *et al.*, 2008; Stephens *et al.*, 2011). There is no molecular mechanism that accounts for the position and homogeneous appearance of pericentric cohesin. To explore whether the size and position of the chromatin subloops dictate the position of cohesin, we measured the mean radial displacement of the tips (most radially displaced beads at initial configuration) of the chromatin subloops (Figure 3A) and cohesin (Figure 3B) over the initial 0.03 s of simulation time. Mean radial displacement is defined as the mean of all of the distances of tips of the subloops or cohesin from the center of the cylindrical array of the chromosome (Figure 1C). The mean radial displacement of the chromatin subloops (10 kb per subloop) collapsed from 389 to 210 nm, while cohesin decreased from 223 to 170 nm radius (Figure 3 and Table 1). We generated centromere models with various subloop lengths (6–23 kb) and found that the size of the subloop dictates the size of the cohesin barrel (Figure 3 and Table 1).

Sub-loop size (kb) | Initial DNA displacement (nm) | Final DNA displacement (nm) | Initial cohesin displacement (nm) | Final cohesin displacement (nm) |
---|---|---|---|---|

6 | 195 | 149 | 116 | 141 |

10 | 389 | 210 | 223 | 170 |

15 | 584 | 261 | 332 | 199 |

20 | 779 | 325 | 442 | 227 |

23 | 876 | 356 | 497 | 275 |

To determine the relationship between DNA subloop size and position of cohesin, we created centromere models with identical 10-kb DNA subloops (initially 389 nm from the spindle axis) but with various initial cohesin barrel sizes. Despite different initial radial displacements, the mean radial displacement of each cohesin distribution converged to the same mean radial displacement over time (Figure 3C). Cohesin barrels with mean radii < 240 nm linking 10-kb subloops show an initial period of expansion relative to the starting configuration before collapse (Figure 3C). Additionally, the cohesin barrel of the centromere model with 6-kb subloops had a greater final than initial mean radial displacement (Figure 3C). These data are suggestive of a repulsive force away from the base of subloops that is overcome by the collapse of loops greater than 6 kb.

Fluorescently labeled cohesin has a characteristic cylindrical structure around the spindle microtubules during metaphase (Yeh *et al.*, 2008; Stephens *et al.*, 2011). When viewed from the side, fluorescently labeled cohesin has a persistent bilobed appearance. To test whether cohesin maintains this bilobed distribution in ChromoShake, we convolved our cohesin bead rings with a point-spread function from a fluorescence microscope used to image cohesin (Smc3-GFP) using the program Microscope Simulator 2 (Quammen *et al.*, 2008). Centromere models with loops greater than 15–20 kb maintained a bilobed geometry when viewed in sagittal section at thermodynamic equilibrium (Figure 4).

### MSD affected by fluorescent array size and position

A recent study showed that the size of a fluorescent particle in a bacterial cell is inversely correlated with the MSD of the particle (Parry *et al.*, 2014). This parameter is critical in studies deriving time constants to distinguish diffusional from subdiffusional motion. We found an inverse correlation between size and MSD by time-lapse imaging live yeast cells and comparing the MSD curves of a 1.7-kb lacO/LacI-GFP array placed 1.1 kb from CEN3 and a 10-kb array placed 1.8 kb from CEN15 (Supplemental Figure S1). The size and placement of tetO/lacO arrays relative to chromosome landmarks such as loops, insulators, and centromeres can affect the MSD plateau value. As predicted by Parry *et al.* (2014) and our experimental measurements (Supplemental Figure S1), we found that simulated array size is inversely correlated with MSD plateau value (Figure 5A and Table 2). Thus the motions of lacO/tetO arrays of different sizes are not directly comparable. We compared the MSD of our simulated 6.8-kb array (10-kb-long, 256 lacO/LacI-GFP repeat, centered 6.8 kb from CEN15) with a simulated 1.7-kb array (1.7-kb-long, 33 repeat lacO/LacI-GFP repeat, centered 1.1 from CEN3) (Pearson *et al.*, 2001; Chacon *et al.*, 2014). The 6.8-kb array is 73 beads, which comprises most of the first subloop, while the 1.7-kb array is 8 beads, which mostly encompasses chromatin proximal to the spindle axis (Figure 5). Surprisingly, the mean plateau value ± SD of the 6.8- and 1.7-kb arrays are 4000 ± 110 and 230 ± 5 nm, respectively (Figure 5B). Given the smaller array had a smaller MSD plateau value, which was not predicted by our previous results, we queried how the placement of the lacO array within a radial loop affects the MSD plateau value. We found a single bead near the base of a loop has a lower plateau value than a single bead at the tip of a subloop, 357 ± 8.7 nm^{2} versus 15,500 ± 510 nm^{2} (Figure 5B), but was higher than the simulated 1.7-kb array (8 beads) value of 230 ± 5 nm, again demonstrating that array size is inversely related to MSD plateau value when the labels are in a similar location. Thus the position of a lacO/tetO array, either within a radial subloop or along the spindle axis, and the size of the array affect the observed motion of a lacO/tetO fluorescent signal.

Number of beads | % of loop | Average plateau value with SD (nm^{2}) |
---|---|---|

1 | 1 | 15,500 ± 510 |

10 | 14 | 14,400 ± 490 |

20 | 29 | 12,600 ± 400 |

30 | 43 | 10,600 ± 340 |

40 | 57 | 8720 ± 280 |

50 | 71 | 7010 ± 210 |

60 | 86 | 5550 ± 160 |

### Both cohesin linkages and condensin-stabilized subloops are needed for experimentally observed correlated motion of different chromosomes

The motion of lacO/tetO arrays within the pericentromere placed on different chromosomes has been shown to be correlated (cross-correlation = 0.33 ± 0.34; Stephens *et al.*, 2013). Because our centromere model is symmetrical, we can place the simulated 6.8-kb array in any of the outermost (centromere proximal) subloops on either end of model. Therefore a single simulation can contain 64 simulated 6.8-kb arrays, with 32 pairs of cohesin-linked arrays, 32 pairs of neighboring arrays not linked by cohesin, and 32 pairs of “sister” arrays. We performed cross-correlation analysis on our simulated 6.8-kb arrays and found a mean ± SD cross-correlation value of 0.34 ± 0.21 for cohesin-linked arrays (Figure 6A, Table 3, and Supplemental Table S2), in quantitative agreement with the experimental value. The motion of neighboring arrays not linked by cohesin in our centromere model are not correlated (−0.053 ± 0.25; Figure 6B and Table 3). Unlike Stephens *et al.* (2013); our centromere model did not show correlated motion in models depleted of cohesin or between sister arrays with or without cohesin and/or condensin (Figure 6C and Supplemental Table S2), indicative of higher-order features in vivo that remain to be incorporated into the model. Models of the centromere lacking cohesin, condensin, or both cohesin and condensin did not exhibit correlated motion (≤0.12) between any neighboring arrays or between sister arrays (Figure 6 and Table 3). Thus, in the model, the stabilization of chromatin loops by condensin and the cross-linking of those loops by cohesin are sufficient to account for experimental levels of correlated motion between different chromosomes. This provides a physical basis for the enigma of cohesin function in the pericentromere, namely to coordinate motion between adjacent DNA strands in the ensemble of centromeric chromatin.

Comparison | With cohesin, with condensin | No cohesin, with condensin | With cohesin, no condensin | No cohesin, no condensin |
---|---|---|---|---|

Cohesin-linked neighbors | 0.34 ± 0.21 | −0.045 ± 0.31 | 0.040 ± 0.35 | 0.12 ± 0.31 |

Unlinked neighbors | −0.053 ± 0.25 | 0.023 ± 0.32 | −0.013 ± 0.35 | −0.022 ± 0.34 |

Sisters | 0.058 ± 0.23 | −0.0028 ± 0.32 | −0.0043 ± 0.37 | 0.015 ± 0.28 |

### Cohesin and condensin confine centromeric chromatin

Probability distributions of population images of the 6.8-kb lacO/LacI-GFP array revealed that depletion of pericentric cohesin (mcm21Δ) or condensin (brn1-9^{ts}), increases the radial displacement of the array and increases the variance of the positions in the population, that is, it expands the area the 6.8-kb array can occupy (Stephens *et al.*, 2011). To test whether depletion of cohesin and/or condensin would have a similar effect on simulated lacO arrays in the pericentromere, we created centromere models without cohesin, without condensin, and without both. We constructed probability distributions from positions of the simulated 6.8-kb arrays on each outer chromatin subloop (Figure 7). In the absence of cohesin, the radial displacement decreased slightly, 140–134 nm, while the distance from the proximal spindle pole body (SPB) decreased from 419 to 336 nm. In the absence of condensin, radial displacement increased from 140 to 163 nm, while the distance from the proximal SPB decreased from 419 to 395 nm. In the absence of both cohesin and condensin, the radial displacement increased from 140 to 186 nm, while the distance from the proximal SPB decreased from 419 to 345 nm. In all deletions, variance (i.e., area of exploration) increased (Table 4), as has been found experimentally (Stephens *et al.*, 2011). The increased radial displacement and variance was also seen in the tips of the loops and the bases of the loops (Supplemental Tables S3 and S4), although loss of cohesin had little effect on the base beads (Supplemental Table S4). Thus both cohesin and condensin function to confine DNA loops around the spindle in an area more compact than is thermodynamically favorable.

Radial displacement | Distance from SPB | |||
---|---|---|---|---|

Average (nm) | SD (nm) | Average (nm) | SD (nm) | |

ATP-depleted cells, 6.8-kb array | 134 | 74 | 331 | 178 |

With cohesin and condensin | 140 | 51 | 419 | 29 |

With condensin, no cohesin | 134 | 60 | 336 | 45 |

With cohesin, no condensin | 163 | 71 | 395 | 76 |

No cohesin, no condensin | 186 | 78 | 345 | 94 |

### Radial side chains and cohesin can extend chromatin

The structure of the centromere model, radially displaced subloops relative to a central primary loop, resembles a bottlebrush polymer, in which multiple side chains emanating from a single primary chain repel one another, resulting in extension of the primary axis. This is contrary to the behavior of a linear polymer chain (Supplemental Figure S2) that collapses into a random coil. In the bottlebrush polymer, all chains are trying to collapse to a random coil. However, the radial chains are crowded, such that the primary chain is extended and tension is generated along the primary chain (Figure 8A). Bottlebrush polymers with side chains of long lengths and high density (number of side chains divided by length of the primary chain) can generate tension along the primary chain in the nanonewton range, which is sufficient to sever covalent bonds within the primary chain (Panyukov *et al.*, 2009a, b; Lebedeva *et al.*, 2012). To determine the density and length of side chains required for extensional effects, we measured the average radius of gyration (*R*_{g}) of a 1-μm (101 bead) primary chain (Figure 8B) of bottlebrush polymers with increasing densities of 200-nm (20 bead) side chains. *R*_{g} essentially measures the volume a polymer occupies and is defined , where *N* is the number of spherical masses, *r*_{k} is the position of an individual mass, and *r*_{mean} is the mean position of all the spherical masses. The polymers were allowed 0.04 s to collapse into a random coil before calculating the *R*_{g}. The *R*_{g} increased when more than 10 side chains were added. Beyond 40 side chains, the primary axis begins to approach an apparent maximal *R*_{g} (190 nm or 76% of the initial *R*_{g} of 250 nm, Figure 8C). Thus the density of side chains within the range and density predicted for loops in vivo (4 loops × 32 chromatids = 128 over ∼1 micron) is sufficient to extend and stiffen chromatin along the primary axis. To explore the effects of side-chain length, we started with a density of 100 side chains on a 1-μm chain and increased the side-chain lengths in 10-nm increments (Figure 8D). The *R*_{g} gradually increased to ∼180 nm or 72% of the initial *R*_{g} (Figure 8D). The length of loops is within the radial distance of experimental loops (∼100 nm from the central axis).

In addition to the radial subloops that function as side chains in the model, bead rings (e.g., cohesin) may function to extend chromatin along the primary axis. To test this hypothesis, we constructed a chromatin ring 1 μm in diameter (314 beads, equivalent to ∼43 kb) with increasing numbers of cohesin rings (∼51 nm in diameter). We found that increasing the number of cohesin molecules also resulted in an increase in the *R*_{g} (Figure 8, E and F). The *R*_{g} was doubled with 70 cohesin rings, corresponding to one cohesin ring every 620 base pairs. In vivo in budding yeast, cohesin rings are spaced every 9–15 kb and show threefold further enrichment in the pericentric region (Blat and Kleckner, 1999; Laloraya *et al.*, 2000; Glynn *et al.*, 2004; Weber *et al.*, 2004). Thus ring-like proteins at sufficient density can stiffen chromatin, providing a mechanism for shaping chromatin structure in vivo. This points to a novel function for ring-like protein complexes that may have significant biological implications.

## DISCUSSION

### A new model for higher-order chromatin structure and dynamics

Thermal forces are ubiquitous and have a dominant influence on molecular dynamics. Chromatin is a dynamic structure, and the particular geometries (loops, crumples, bodies, etc.) will adopt unique configurations not necessarily intuitive in the absence of molecular modeling. Here we describe a computational tool for exploring how specific chromatin geometries are altered in the presence of thermal motion. Using ChromoShake, researchers can generate models of chromatin and explore the consequences of thermal force on predicted structure. ChromoShake data can be compared with experimental data to determine the structure and forces affecting chromatin in living cells. Our model of the yeast centromere illustrates how chromatin structure can affect protein and DNA distribution and reveals functional attributes not previously reported.

### Size and position of tetO/lacO arrays affect their observed motion

Tandem-repeat arrays of lacO/tetO have been used to label chromatin for decades (Robinett *et al.*, 1996; Straight *et al.*, 1996; Michaelis *et al.*, 1997). We have used these arrays to demonstrate that, on average, the DNA in the pericentromere, like cohesin, is radially displaced from the spindle axis (Anderson *et al.*, 2009; Stephens *et al.*, 2011). A major difference in the motion of DNA arrays in vivo versus in silico is that the in vivo motion is greater than that found in silico (Supplemental Figure S1 vs. Figure 2). While depletion of ATP in vivo largely recapitulates the in silico motion (Figure 2B), the biological condition is much more violent than in simulation. Factors such as remodeling complexes, chaperones, and crowding are constantly bombarding the chromatin. The net effect is increased motion with a random signature, due to the random trajectory of the sum of a number of different inputs. This can be simulated by increasing temperature in silico. One of the outstanding questions is whether there is information in the in vivo trajectories that might reveal hierarchies in the source of this noise.

A second insight from ChromoShake is that size of the tracer (lacO) influences the slope and magnitude of the MSD curve (Figure 5). MSD is a widely used metric for deducing the form of motion (subdiffusional, diffusional, ballistic) and the degree of confinement (for cell work in particular). The size of the array influences the magnitude, due to the averaging that is inherent as the arrays increase in size. The observed motion of one bead will be greater than the motion of the average of 10 beads. A simple analogy would be motion of a crowd in a stadium. The average motion of the entire crowd is static, the stadium remains full. However, the motion of any given individual will vary greatly. This phenomena was recently reported by Parry *et al.* (2014) and is recapitulated herein. Thus, in studies in which time constants are being deduced for models of motion, the size of the arrays must be considered in the resolution of the measurements.

### Homogeneity and function of cohesin in vivo

Cohesin is widely cited for its function in holding sister chromatids together until the appropriate transition from metaphase to anaphase. In the centromere, cohesin is enriched (three times; Blat and Kleckner, 1999; Weber *et al.*, 2004), yet sister centromeres are on average 800 nm apart (Pearson *et al.*, 2001). In addition, cohesin appears as a homogeneous structure radially displaced from the spindle axis (Yeh *et al.*, 2008; Stephens *et al.*, 2011). There is little understanding of the positional determinants or the functional attributes for cohesin enrichment at the centromere. ChromoShake provides critical predictions for these conundrums. First, the distribution of cohesin rings represents the thermodynamically favored position that depends on the size of DNA loops. Both cohesin rings and chromatin loops jiggle toward their entropically favored positions. Cohesin rings proximal to the primary axis initially migrate radially outward (Figure 3C, 180- and 200-nm barrels), away from the crowded primary axis. Conversely, chromatin loops crumple toward the axis, pushing cohesin rings toward the primary axis (Figure 3, A and B). While it is likely that additional biological components may contribute to cohesin’s function in cells, the thermodynamic argument reveals that no other components are required for its unique homogeneity and radial position. Second, we show that cohesin functions to stiffen the DNA to which it is bound (Figure 8, E and F). This is presumably through an excluded volume effect, restricting the DNA from adopting a random coil. Considering the centromere’s function as a spring between sister kinetochores and its role in tension sensing, the role of cohesin in stiffening the spring has important biological consequences.

### Biological consequences of a slip-link network

Ring complexes cross-linking radial loops impart correlated motion across the centromere (Figure 6). In the model, the cohesin rings cross-link different chromatin strands and are free to migrate along the strands. This is a slip-link pulley, reminiscent of “topological gels” that maintain form and flexibility over orders of magnitude of volume changes. Topological gels retain their elastic and tensile moduli properties over several length scales as well (Okumura and Ito, 2001; Granick and Rubinstein, 2004). The cross-linked organization of the centromere has the effect of dissipating changes in motion at the individual kinetochore microtubule attachment sites across the 32 DNA strands. The timescale of DNA motion (Supplemental Figure S2) is much faster than the relatively slow microtubule dynamics (micron/min). Therefore the ensemble of protein rings and individual chromatids behave as a cross-linked network. Cross-linking the network has the attribute of equalizing stochastic dynamics in microtubule growth and shortening across the 16 sister kinetochores. Rather than reading the micromechanics at 32 individual kinetochore–microtubule attachment sites, cross-linking the centromeres averages the mechanics of individual kinetochores across multiple kinetochores. This feature reveals that SMC complexes create local networks or domains within a single or multiple chromosomes for specific regulatory processes, such as segregation. Thus the concentration of SMCs in the pericentromere have the capacity to mold a multichain chromosome domain into a topological network uniquely tuned to integrate the stochastics of microtubule dynamics into a singular output for informing the cell of the state of chromosome attachment.

### Thermodynamic consequences of radial loops

While the radial-loop structure of chromosomes was first noted more than a century ago (Hertwig, 1910) and has been studied in detail (Diaz *et al.*, 1981), the thermodynamic consequences of that structure remain poorly understood. ChromoShake allows researchers to study how the chromatin structure affects chromosome dynamics in detail. We find that a radial subloop structure excludes cohesin molecules from the primary axis as cohesin molecules placed near the primary axis quickly migrate radially outward (Figure 3C). This result is reminiscent of the finding that histones are more quickly turned over in the pericentric region than in the arm region of chromosomes (Verdaasdonk *et al.*, 2012). In addition, the bottlebrush-like shape could also explain the increased turnover of histones, as a bottlebrush polymer has increased tension on the primary axis (Panyukov *et al.*, 2009a, b; Lebedeva *et al.*, 2012). This tension is evident in our bottlebrush chromatin models, as they fail to collapse into a random coil (Figure 8, A–D).

ChromoShake has provided novel insights into centromere function. The centromere is not a passive DNA element, but the unique configuration of loops and rings impart stiffness to the structures that undoubtedly are critical elements in the force balance required for faithful chromosome segregation.

## MATERIALS AND METHODS

### ChromoShake workflow and parameters

The ChromoShake simulator parses a configuration file that contains the coordinates of the spherical masses and the polymer and environmental parameters. These configuration files are generated by C++ programs. Users can specify whether ChromoShake is run on a central processing unit or a graphics processing unit (GPU), the frequency at which output occurs, and the duration of the simulation. ChromoShake outputs a .out text file that can be parsed by various analysis programs and by programs that render the model. An overview of ChromoShake’s workflow is provided in Figure 9. A full list of tunable parameters is given in Table 5. An installer for ChromoShake, the source code, and a user guide are provided in the Supplemental Material.

Parameter | Default value | Parameter specification program |
---|---|---|

Bead separation | 0.01 μm | All configuration programs |

Temperature | 25°C | All configuration programs |

Viscosity | 1 cP | All configuration programs |

Drag/damping radius factor | 0.8 (of bead separation) | All configuration programs |

Time step | 2 ns | All configuration programs |

Modulus of DNA | 2 GPa | All configuration programs |

DNA radius | 0.6 nm | All configuration programs |

Collision radius factor | 0.25 (of bead separation) | All configuration programs |

Random seed | 42 | All configuration programs |

Cohesin distribution radius | 0.22 μm | Centromere configuration program |

Cohesin ring radius | 0.0255 μm | Centromere configuration program |

Collisions | On | ChromoShake simulator |

### Polymer dynamics mathematical model

We model a given polymer’s conformation by approximating the polymer structure as a number of spherical masses. The dynamics of the polymer’s conformation is modeled as the result of deterministic and stochastic forces acting on the masses, yielding a system of stochastic differential equations that can be numerically integrated to simulate the polymer system. We proceed to describe these forces, which model the bending rigidity and tensile stiffness of the polymer, the friction and thermal fluctuations the polymer experiences, and the result of collisions involving the polymer.

#### Description of forces in the model

##### Bending rigidity.

The polymer model possesses the bending rigidity of B-form DNA. In vivo, DNA is compacted into chromatin fibers. The reported bending rigidity of chromatin, typically expressed as a persistence length, is only marginally different from the 50 nm for B-form DNA (Cui and Bustamante, 2000; Dekker *et al.*, 2002), indicating the bending mechanics is dominated by the DNA linkers. We therefore use the DNA bending rigidity values as the default value for our simulation parameters. The bending rigidity parameters for the simulation are determined by first treating the DNA as an elastic beam and applying standard Euler–Bernoulli beam theory (Supplemental Figure S3). The restoring forces as determined through the continuum picture are then translated into a hinge rigidity within the discretized model of the simulation.

We estimate the bending-induced restoring forces for each mass within the bent chain by using the Euler–Bernoulli beam-bending theory. Our estimate is based on the case of a simply supported beam (Supplemental Figure S3) with force applied at its center. These boundary conditions provide well-defined restoring forces at each end and at the center of the beam. For the simply supported beam with central load case shown in Supplemental Figure S3, the relation between the central loading force and the resultant tangent angle at each end is given by (Roark and Young, 1989)

((*1*)) where

*E*and

*I*are the Young’s modulus and moment of area of the beam, respectively (

*I =*π

*r*

^{4}/4 for a beam of circular cross-section), and

*d*is the half-length of the beam and is a variable that can be tuned within the simulator. For our simulations, we use 2 GPa for the Young’s modulus, 0.6 nm as the radius of B-form DNA, and 10 nm for

*d*.

Within the model, there will be an elastic restoring force at each mass proportional to the relative angle of adjacent straight segments connecting that mass to its neighbors and in a direction of the bisector of the included angle between segments. A force of half this magnitude and in the opposite direction is applied to both neighboring masses (ϕ in Supplemental Figure S3D). Supplemental Figure S3 illustrates how, within the mass chain, each hinge angle is defined and the corresponding restoring forces are associated with each hinge bend. The force magnitude associated with each hinge bend, given by

((*2*))

is applied to the central and adjacent masses as indicated in Figure S3D.

##### Tensile stiffness.

The chain includes a tensile stiffness that resists extension of the segments connecting each mass. Again, we used the literature values for the Young’s modulus of DNA to set the magnitude of the segment stiffness. The standard stress versus strain equation describes tensile stretching of a beam:

((*3*)) ((

*4*)) ((

*5*))

Substituting Eqs. 4 and 5 into Eq. 3 and rearranging yields the restoring force for tensile extension, Δ*L*, of the segments joining each mass:

*6*))

The quantity within the parentheses of Eq. 6 acts as a spring constant for the segments joining the masses. Like a conventional spring, it scales inversely with the length of the segment. For the chosen DNA parameters of *E* = 2 GPa and *r* = 0.6 nm, this yields a “spring” constant for each segment of

*7*))

##### Drag force.

We describe a drag force on a spherical mass moving with velocity, , by , where γ is the drag constant. Owing to the spherical shape of the masses and the small Reynolds number system, we obtain the drag constant as a function of our viscosity (η) and the radius of the sphere (*a*) according to Stokes’ law: . By default *a* is 8 nm (0.8 × node separation of 10 nm).

##### Thermal fluctuations.

The result of thermal fluctuations is inherently stochastic, yielding a Brownian-like motion on individual particles. These fluctuations can be modeled as a normally distributed random force averaged over any chosen time interval, Δ*t*. To estimate the magnitude of the random forces, we first consider the root-mean-squared (RMS) displacement of a freely diffusing particle in one dimension:

*8*)) where is the RMS displacement,

*D*is the diffusion constant of the particle, and Δ

*t*is the time interval. The diffusion constant is given by ((

*9*)) where γ is the drag coefficient of the particle. This is the Stokes–Einstein equation relating diffusion with viscous drag. For a real diffusing particle, the thermal force will fluctuate in magnitude and direction and will have a non-zero time average over a given time interval resulting in finite displacement. This time-averaged thermal force can be estimated by approximating the motion of the particle over the time interval as directed motion through a viscous fluid. The motive force required to balance the drag force of a sphere moving at speed

*v =*Δ

*x*/Δ

*t*is given by ((

*10*))

We set the SD of the random force in one dimension equal to the force in the previous equation. This yields

((*11*))

The drag coefficient for a sphere is

((*12*)) where η is the fluid viscosity and

*R*is the radius of the sphere (set as bead separation). Thus, over chosen intervals Δ

*t*, we model the thermal fluctuations on each sphere as a normally distributed random force with mean centered at zero with an SD of

*F*

_{B}in each dimension.

##### Collision force.

While the spring and hinge forces prevent neighboring masses of the polymer model from crossing through each other, they do not prevent distant sections of the polymer from overlapping. To prevent such behavior, we included an additional collision force in the model. This collision force is modeled as a repulsive Hookean spring force between spheres when they are closer than an equilibrium distance, discouraging spheres from overlapping or crossing through each other. The repulsive force is defined as , where is the collision force, *k*_{c} is the universal collision spring constant, and *x* = 2*R*_{c} − *l*, where *R*_{c} is the collision radius and *l* is the distance between the center of the two masses. When spheres are farther apart than the equilibrium distance, *l* ≥ 2*R*_{c}, the collision force is zero.

### Creating geometric models

When constructing the model, we must specify the locations of all the spherical masses and the respective parameters of the forces holding them together. Several spherical masses are defined in three dimensions in external model configuration files that are parsed by ChromoShake. The model configuration files can be created by C++ programs and provide ChromoShake with the relevant polymer information (mass and location of the spherical masses, location and strength of hinges and springs, etc.) and environmental information (temperature and viscosity) to construct and simulate thermal forces on a polymer model. Programs to generate basic geometric shapes (i.e., linear chains and circles) are provided in the Supplemental Material. The geometric model of the yeast centromere was based on a 3D geometric model we originally produced using the open-source Blender rendering program (Frankel and DePace, 2012). The control points of the splines (subloops) were taken from the Blender model, and a C++ program was created that generated nonuniform rational B-splines (NURBs) to recreate its loops. The code was adapted from programs in *An Introduction to NURBS* (Rogers, 2001). The NURBs program to generate the model of the centromere is provided in the Supplemental Material.

### ChromoShake program

ChromoShake loads a 3D polymer model, introduces thermal forces, and evolves the state of the model with respect to time by solving the differential equations of motion at a given temperature and viscosity. In essence, the solver keeps track of the state (positions and velocities) of a system of multiple spherical masses. Forces are generated to model different types of interactions, which are summed for each mass. Jacobian matrices of the derivative of the forces with respect to the positions and velocities are computed and solved implicitly. The solver describes how the sums of these forces and Jacobians affect the evolution of the system state over discrete time steps. Given the potentially large number of masses (>20,000) in the polymer models, ChromoShake was written in OpenCL to enable parallelization of the code with GPUs. ChromoShake’s solver uses the conjugate gradient method to solve the linearized implicit equation for position and velocity. The solver was constructed using the ViennaCL sparse linear algebra library http://viennacl.sourceforge.net, enabling it to run in parallel using OpenCL. Moreover, ViennaCL enabled the use of compressed sparse-row sparse-matrix data structure, significantly simplifying the differential equation–solver algorithm for faster performance.

### Collision detection

Given the large number of masses in a polymer model, accelerated collision detection can be used for large (mass number >>100) polymer models. Collisions can be detected in two ways using ChromoShake. The simplest but most computationally expensive method is to compare the location of each sphere with the location of every other sphere. For acceleration of collision detection, a spatial subdivision scheme was implemented to divide the space in the simulation into a finite grid of 3D boxes. The sizes of these boxes are set so that a mass whose center is in one box can only intersect with masses whose centers are either in the same box or one of the 26 adjacent boxes. The scheme uses a partially parallelized collision detection algorithm for use with GPUs, based on *GPU Gems 3* (Nguyen and NVIDIA Corporation, 2008; chap. 32). This scheme is much faster and is the default scheme for machines with GPUs.

### ChromoShake output

The output file from ChromoShake contains three parts. The first is the configuration file that was used to generate and run the model. The second is a series of integers that specify the color of each mass when rendered (used to visually segment the model). The last section contains the *x*-, *y*-, and *z*-coordinates of all the simulated masses at each time point.

### Rendering ChromoShake output

ChromoShake output files can be rendered in two ways. The first is using the ChromoView program that is packaged in the included Windows installer. The second uses a custom Python script to convert ChromoShakes’s output file into a .blend file for the open-source rendering software Blender (www.blender.org).

### Diffusion validation

An automated C++ program runs the ChromoShake simulator on 100 spherical masses with collision detection turned off. The lack of collisions therefore makes diffusion from thermal forces the only external force acting on the masses and provides 100 independent replicates per run. The squared displacement of each bead is calculated over a specified time step. The default time step is 9 ns, with each bead taking 10 steps before an MSD calculation. Then the distribution of each bead’s MSD is compared with the expected values based on , where *k*_{B} is Boltzmann’s constant, *T* is the temperature, η is the viscosity, *R* is the radius of the spherical bead, and Δ*t* is the elapsed time (time step × step number, 90 ns; Supplemental Figure S4). Supplemental Figure S4B shows that the beads are freely diffusing, since plots of the have a slope of ∼1 (Saxton, 1994; Rubinstein and Colby, 2003). Because diffusion is a random process, given enough iterations, approximately half of the values should be below the expected value and half should be above. We define a diffusion simulation as successful if more than 10% of calculated MSDs are greater than and more than 10% are less than the expected value, that is, the values are distributed around the expected value (Supplemental Figure S4C).

### Rouse relaxation time validation

The linear chain models are generated in an extended state and then collapse into a random coil; therefore we can compare the time it takes for a chain to collapse against the predicted Rouse relaxation time (Rouse, 1953; Rubinstein and Colby, 2003). The Rouse relaxation time is, , where η is viscosity, *b* the Kuhn length, *k*_{B} the Boltzmann constant, *T* temperature, and *N* the number of Kuhn length segments in the chain. The time to collapse, determined by viewing a plot of the radii of gyration and end-to-end lengths over time, were consistent with the calculated Rouse relations times (Supplemental Figure S2). The good agreement between measured and predicted values is a strong indicator that our simulation accurately models polymer dynamics.

### Radius of gyration and end-to-end distance validation

Linear polymer models were generated and simulated using ChromoShake. The radius of gyration, , where *N* is the number of spherical masses, *r*_{k} is the position of an individual mass, and *r*_{mean} is the mean position of all the spherical masses, and the end-to-end length, distance between the end masses, was calculated for 1-, 2-, and 3-micron chains. Radii of gyration and end-to-end lengths were measured after twice the Rouse time and compared with expected values given and , where *b* is Kuhn length and *N* is the number of segments (Rubinstein and Colby, 2003). Average radii of gyration and end-to-end lengths were similar to expected values (Supplemental Table S5 and Supplemental Figure S2).

### Strain information and growth conditions

The strains KBY8065 (MATa CEN(15)(1.8)-GFP[10kb] ade2-1, his3-11, trp1-1, ura3-1, leu2-3112, can1-100, LacINLSGFP:HIS3, lacO::URA3, Spc29-RFP:Hyg) and KBY9471 (MATa, YEF473a, trp1-63 leu2-1 ura3-52 his3-200 lys2-801, Spc29-RFP:Hyg, Smc3-GFP:URA(pLF639)), were grown in rich yeast–peptone–dextrose media with 0.5 mg additional adenine per milliliter of media to logarithmic phase at 24°C. For depletion of ATP from KBY8065, 20 min before imaging, cells were washed and resuspended in yeast casamino acids (YC)-complete media containing 0.5 mg additional adenine per milliliter of media, 0.02% sodium azide, and 1 μM of deoxyglucose, and incubated at 24°C.

### Microscopy

Cells containing lacO/LacI-GFP (KBY8065) were imaged in YC-complete media containing 0.02% sodium azide, and 1 μM of deoxyglucose for 10 min every 30 s using a Nikon Eclipse Ti wide-field inverted microscope with a 100× Apo TIRF 1.49 NA objective (Nikon, Melville, NY) and Andor Clara CCD camera (Andor, South Windsor, CT) using Nikon NIS Elements imaging software (Nikon) at room temperature (25°C). At each interval, a 7 Z-plane section image stack with a 400-nm step size was taken.

Cells containing fluorescently labeled cohesin (KBY9471) were imaged in YC-complete media with 2% filter-sterile glucose using a an inverted, wide-field microscope (Eclipse TE2000-U; Nikon) with a 100× Plan Apo 1.4 NA digital interference contrast oil-immersion lens with an Orca ER camera (Hamamatsu Photonics, Bridgewater, NJ) with MetaMorph 6.1 software at room temperature (25°C). At each interval, a 5 Z-plane section image stack with a 300-nm step size was taken.

### Image analysis

The ND2 files were converted to TIFF file format using either Nikon NIS Elements or ImageJ (National Institutes of Health, Bethesda, MD; http://imagej.nih.gov/ij). The in-focus planes for each GFP and RFP foci for each Z-stack were determined for each time point and then combined into separate image stacks using MetaMorph 7.7 imaging software (Molecular Devices, Sunnyvale, CA). The in-focus GFP and RFP signals were tracked using the MATLAB program Speckle Tracker (Mathworks, Natick, MA; Wan *et al.*, 2009, 2012) which is a graphical user interface that detects RFP and GFP foci in an image, tracks the foci through the time lapse, performs 2D Gaussian fitting on a 5 × 5 pixel region using MATLAB’s nonlinear curve-fitting methods (*lsqcurvfit*), and exports the subpixel *x-* and *y*-coordinates to an Excel spreadsheet (Microsoft, Redmond, WA). The in-focus plane of the cohesin barrel (Smc3-GFP) was defined as the plane with the highest SPB signal intensity using MetaMorph 7.7.

### Statistical analysis

Calculation of radial displacement was performed using the C++ program cohesin_analysis.cpp. This program calculates all of the distances from the origin (center of the model in Figure 1C) for the beads specified by an input file. The program then calculates the mean distance, the maximum distance, and the minimum distance for each time point and then prints the information. Calculating the average position for multiple bead arrays was done using the C++ coord_summary.cpp. MSD of bead position and the average positions of arrays were calculated using the PERL script MSD_3D.pl. The correlated motion of the average position of the simulated 6.8-kb array was performed using a custom MATLAB script with the *corrcoef* function. Heat maps of experimental and simulated data were generated using two custom MATLAB scripts. Heat maps of experimental data were performed on the subpixel coordinates generated from the Speckle Tracker MATLAB program. Heat maps of the simulation were generated using a custom MATLAB program that parsed the coordinates from a ChromoShake output file, calculated the distance of each *x*-coordinate from the proximal SPB, and calculated the distance of each *y*-coordinate from the spindle axis. The spindle length was set as the average spindle length, 1.27 μm, from population of the azide and deoxyglucose-treated cells. Both experimental and simulated heat maps were mirrored about the Y-axis. Only the center 200 nm of the Z-axis was sampled for simulated heat maps. Radius of gyration analysis, , of bottlebrush polymers and the circular chromatin polymer was performed using a custom MATLAB script. Radius of gyration analysis for polymer chains was performed using the C++ program radius_of_gyration.cpp. Tukey’s honest significant difference test was perform using the functions anova1 and multcompare(stats,‘alpha’,0.05,‘ctype’,‘hsd’) in MATLAB. The mean diameters of individual cohesin rings were calculated by multiplying the mean *R*_{g} that was calculated using custom MATLAB programs by two. For cohesin not encompassing chromatin, three independent simulations were run for 0.05 s of simulation time. For cohesin encompassing a single chromatin strand, one simulation with 10 cohesin rings on the same circular chromatin strand was run for 0.08 s of simulation time. For cohesin encompassing two chromatin strands, the 48 cohesin rings of the centromere model that was run for 0.10653 s of simulation time were measured. The mean diameter for the final 2000 time points of each simulation was calculated. The ensemble mean and SD between the time-averaged diameters of each cohesin ring are reported in Table S1.

## FOOTNOTES

This article was published online ahead of print in MBoC in Press (http://www.molbiolcell.org/cgi/doi/10.1091/mbc.E15-08-0575) on November 4, 2015.

**Abbreviations used:**

3D | three-dimensional |

GFP | green fluorescent protein |

GPU | graphics processing unit |

lacO | lactose operon |

MSD | mean-squared displacement |

NURBs | nonuniform rational B-splines |

R_{g} | radius of gyration |

RMS | root-mean-squared |

SMC | structural maintenance of chromosomes |

SPB | spindle pole body |

tetO | tetracycline operon. |

## ACKNOWLEDGMENTS

We thank Belinda Johnson for her initial contributions to a previous version of ChromoShake. This work was funded by National Institutes of Health (NIH) R37 grant GM32238 (to K.B.), NIH T32 grant GM 007092-39 (to J.L.), and NIH P41 grant EB002025 (to R.M.T. in the Center for Computer Integrated Systems for Microscopy and Manipulation).

## REFERENCES

**Mol Biol Cell**

*20*, 4131–4139. Link, Google Scholar (

**Cell**

*98*, 249–259. Crossref, Medline, Google Scholar (

**Chromosoma**

*117*, 103–110. Crossref, Medline, Google Scholar (

**J Cell Biol**

*205*, 313–324. Crossref, Medline, Google Scholar (

**Proc Natl Acad Sci USA**

*97*, 127–132. Crossref, Medline, Google Scholar (

*cis*-acting sites for condensin loading onto budding yeast chromosomes.

**Genes Dev**

*22*, 2215–2227. Crossref, Medline, Google Scholar (

**Science**

*295*, 1306–1311. Crossref, Medline, Google Scholar (

*Notophthalmus*.

**Cell**

*24*, 649–659. Crossref, Medline, Google Scholar (

**Visual Strategies: A Practical Guide to Graphics for Scientists & Engineers**, New Haven, CT: Yale University Press Google Scholar (

*Saccharomyces cerevisiae*.

**PLoS Biol**

*2*, E259 Crossref, Medline, Google Scholar (

**Nat Mater**

*3*, 586–587. Crossref, Medline, Google Scholar (

**Curr Biol**

*22*, 471–481. Crossref, Medline, Google Scholar (

**Lehrbuch der Entwicklungsgeschichte des Menschen und der Wirbeltiere, Jena**, Germany: G. Fischer Google Scholar (

**Curr Biol**

*21*, 12–24. Crossref, Medline, Google Scholar (

**J Cell Biol**

*151*, 1047–1056. Crossref, Medline, Google Scholar (

**J Cell Biol**

*210*, 553–564. Crossref, Medline, Google Scholar (

**Proc Natl Acad Sci USA**

*109*, 9276–9280. Crossref, Medline, Google Scholar (

**Mol Cell**

*4*, 445–450. Crossref, Medline, Google Scholar (

**Cell**

*91*, 35–45. Crossref, Medline, Google Scholar (

**GPU Gems 3**, Upper Saddle River, NJ: Addison-Wesley Google Scholar (

**Adv Mater**

*13*, 485–487. Crossref, Google Scholar (

**J Phys Chem B**

*113*, 3750–3768. Crossref, Medline, Google Scholar (

**Phys Rev Lett**

*102*, 148301 Crossref, Medline, Google Scholar (

**Cell**

*156*, 183–194. Crossref, Medline, Google Scholar (

**J Cell Biol**

*152*, 1255–1266. Crossref, Medline, Google Scholar (

**Eurographics Workshop Vis Comput Biomed**

*2008*, 151–158. Medline, Google Scholar (

**Roark’s Formulas for Stress and Strain**, New York: McGraw-Hill Google Scholar (

**J Cell Biol**

*135*, 1685–1700. Crossref, Medline, Google Scholar (

**An Introduction to NURBS: With Historical Perspective**, San Francisco: Morgan Kaufmann Google Scholar (

**J Chem Phys**

*21*, 1272–1280. Crossref, Google Scholar (

**Polymer Physics**, Oxford, UK: Oxford University Press Google Scholar (

**Biophys J**

*66*, 394–401. Crossref, Medline, Google Scholar (

**J Cell Biol**

*207*, 189–198. Google Scholar (

**J Cell Biol**

*193*, 1167–1180. Crossref, Medline, Google Scholar (

**J Cell Biol**

*203*, 407–416. Crossref, Medline, Google Scholar (

**Curr Biol**

*6*, 1599–1608. Crossref, Medline, Google Scholar (

**Cell**

*98*, 847–858. Crossref, Medline, Google Scholar (

**Mol Biol Cell**

*23*, 2560–2570. Link, Google Scholar (

**Mol Cell**

*52*, 819–831. Crossref, Medline, Google Scholar (

**Mol Biol Cell**

*23*, 1035–1046. Link, Google Scholar (

*et al.*(2009). Protein architecture of the human kinetochore microtubule attachment site.

**Cell**

*137*, 672–684. Crossref, Medline, Google Scholar

**PLoS Biol**

*2*, E260 Crossref, Medline, Google Scholar (

**Proc Natl Acad Sci USA**

*109*, 7338–7343. Crossref, Medline, Google Scholar (

**Curr Biol**

*18*, 81–90. Crossref, Medline, Google Scholar (