A Motor-Gradient and Clustering Model of the Centripetal Motility of MTOCs in Meiosis I of Mouse Oocytes

Asters nucleated by Microtubule (MT) organizing centers (MTOCs) converge on chromosomes during spindle assembly in mouse oocytes undergoing meiosis I. Time-lapse imaging suggests that this centripetal motion is driven by a biased ‘search-and-capture’ mechanism. Here, we develop a model of a random walk in a drift field to test the nature of the bias and the spatio-temporal dynamics of the search process. The model is used to optimize the spatial field of drift in simulations, by comparison to experimental motility statistics. In a second step, this optimized gradient is used to determine the location of immobilized dynein motors and MT polymerization parameters, since these are hypothesized to generate the gradient of forces needed to move MTOCs. We compare these scenarios to self-organized mechanisms by which asters have been hypothesized to find the cell-center- MT pushing at the cell-boundary and clustering motor complexes. By minimizing the error between simulation outputs and experiments, we find a model of “pulling” by a gradient of dynein motors alone can drive the centripetal motility. Interestingly, models of passive MT based “pushing” at the cortex, clustering by cross-linking motors and MT-dynamic instability gradients alone, by themselves do not result in the observed motility. The model predicts the sensitivity of the results to motor density and stall force, but not MTs per aster. A hybrid model combining a chromatin-centered immobilized dynein gradient, diffusible minus-end directed clustering motors and pushing at the cell cortex, is required to comprehensively explain the available data. The model makes experimentally testable predictions of a spatial bias and self-organized mechanisms by which MT asters can find the center of a large cell.


Introduction
Spindle assembly in higher eukaryotic cells involves the self-organization of microtubules (MT) into a bipolar structure. During mitosis in animal cells, spindle poles are defined by a pair of centrosomes. However bipolar structures emerge even in the absence of centrosomes during meiosis in vertebrates as well as mitosis in plants. In such acentrosomal spindles, the poles self-organize by the dynamic interactions of MTs with molecular motors, regulatory factors and chromatin. While multiple components of this cellular-scale pattern forming system have been identified, the precise nature of the interactions between the components are still not completely understood.
The meiotic maturation of mouse oocytes is a well studied example of such an acentrosomal spindle assembly system. The first meiotic division is characterized by germinal vesicle breakdown (GVBD) [1], before and after which small aster-like fibrillar structures or microtubule organizing centers (MTOCs) are observed [2]. MTOCs which are nucleated both in the cytoplasmic and peri-nuclear spaces, both aggregate at the center to form a spindle by prometaphase I [3]. Such a convergence of radial MT arrays or asters was reported previously in Xenopus meiosis II oocytes [4]. Using cell-free Xenopus oocyte extracts, this convergence was shown to result from asymmetric centrosomal MT growth due to a gradient of RanGTP [5][6][7]referred to as biased 'search-and-capture' . However during meiosis I in mouse oocytes, experimental perturbation of RanGTP levels does not significantly affect spindle assembly [8,9]. If RanGTP does not act as a guidance cue as reported previously [10], the nature of the directional cue and force generation remains to be understood.
The force required for MTOC convergence to the nuclear region is thought to originate from a combination of MTs, motors and anchorage points. Multiple mechanisms have been reported in the past to drive radial MT array transport in cells-(a) polymerization dependent pushing forces as seen during the centering of asters in vitro [11,12], (b) cortical force-generator based pulling [13], (c) cortical motors which both depolymerize and pull [14], (d) cytoplasmic minusended motors which pull asters in a length-dependent manner [15], (e) cytoplasmic streaming by cargo transport driving aster movement [16,17] and (f) acto-myosin contractility as seen in starfish oocytes [18]. Contact with the cell cortex can move asters when the relative MT lengths is comparable to the cell radius [19]. Both active and passive mechanisms drive the movement of centrosome nucleated asters. However most of the cortical pushing and pulling models are unlikely to affect long-range movement of MTOCs which have MT lengths *3 μm as compared to the cell-radius of *40 μm. Transport of asters by cytoplasmic streaming based on cargo transport by one large aster [17], is also unlikely to drive mouse meiosis I oocyte MTOCs due to their size and number (*80 to 100), which will prevent a coherent and directed flow. Inhibition of acto-myosin contractility has also been shown to have no effect on the centripetal movement of mouse MTOCs [9]. While centrally-anchored MTOCs and cross-linking motors have also been proposed by Schuh et al. [9] to drive the MTOC motility, the movement continues even after nuclear envelope breakdown (NEBD). Thus for a complete theoretical understanding of the mechanism by which MTOCs converge in spindle assembly a mathematical model of the process is necessary to test multiple hypotheses that have been proposed.
Theoretical models have been used to probe the interactions of microtubules and motor complexes and are capable of reproducing in vitro self-organized patterns [20][21][22]. These simulations have been extended to understand the role of multiple components in spindle assembly such as antiparallel interactions [23], pole focussing by minus-end directed motors [24], gradients of stabilization [25] and intra-spindle nucleation and dynamic instability regulation [26]. In recent work, we have demonstrated the centripetal movement of centrosomal MT asters towards surface immobilized chromatin in Xenopus egg extracts can be modeled by a gradient of polymerization dynamics and uniform motor distribution [27]. This is comparable to a model of length-dependent pulling by motors to translocate MT asters during C. elegans embryogenesis [15]. However, neither of these models take into account the relatively shorter MTs seen in MTOC asters, and lack details specific to meiosis I. In search of common design principles in spindle assembly, theoretical modeling of the centripetal motility of MTOC arrays can be used to test the generality of previous results.
Here, we quantify the spatial trends in MTOC motility and find the random and directional components of motility depend on how far the MTOCs are from the cell center. The detailed quantitative analysis allows us to develop and test theoretical models of random walk with drift. Only a spatial gradient of drift can reproduce the experimental data. Such an optimized gradient is further used to model MT dynamic instability and motor distributions, to test the combination of mechanisms that can reproduce the experimental statistics of centripetal MTOC motility.

Models
We have developed two kinds of models-phenomenological and mechanistic-to address both the general principles and specific mechanisms of the centripetal convergence of small MT asters, MTOCs. The models are: 1. 2D Random walk with drift (RWD) model

2D MT-motor model
The outputs of both models are compared to 2D experimental motility measures of MTOCs from mouse oocyte meiosis I reported previously by Schuh and Ellenberg [9] The RWD model is used to optimize the functional form of the spatial drift field by comparison to experiment, while the MT-motor model is used to test molecular mechanisms which could generate the drift field. Mechanisms that have been previously hypothesized to drive asters to the center of cell involve MT-pushing and pulling [28,29]. A combination of pulling and pushing mechanisms has been experimentally tested in sand-dollar eggs [30] and C. elegans embryos [16,31] and pulling alone in Xenopus egg extracts [5,27], while pushing has been seen in fibroblasts [32]. The MT-motor model is used to test whether the reported self-organized clustering of MTOCs in mouse oocytes [9] is sufficient to result in MTOC convergence, or whether additional mechanisms are necessary. These mechanisms are: 1. Self-organized mechanisms: not requiring explicit spatial localization a. Cortical pushing alone b. Clustering motors and cortical pushing 2. Spatial gradients: requiring explicit spatial localization a. MT dynamic instability gradient b. Dynein motor gradient 3. Hybrid model: combination of self-organized and gradient mechanisms An optimization routine has been developed to compare the outputs of the 'scenarios', i.e. combinations of these mechanisms, to reproduce previous experimental reports of mouse meiotic MTOC motility and distribution [9,33,34].
The choice of a 2D model for a 3D spherical oocyte is determined by the need to compare simulated motility outputs with the only available experimental time-series dataset of MTOC motility in mouse oocyte meiosis I, which is 2D over time [9] (Fig 1A). This compatibility of dimensions is essential since some of the RWD model input parameters are obtained from fits to experiment and MT-motor model mechanisms are optimized based on their ability to reproduce experimental statistics. Additionally, the choice of dimensionality of spatial models is considered to be determined by the balance between the need to capture the behavior of the system sufficiently and the clarity of the model [35].

Random walk in a drift field
The mouse oocyte is modeled in a 2D circular geometry of radius r cell , with concentric circular chromatin of radius r chr . The outer cytoplasmic region has a radius r cyto and r cell = r chr + r cyto (Fig 1B). MTOC asters are modeled as point particles, nucleated uniformly in the cytoplasmic space of the oocyte. The motion of the simulated particles is a mixture of random Brownian and directed centripetal motion, depending on the position of the MTOC in the oocyte ( Fig  1B), based on the previously observed 'stop and go' nature of the motility [9]. In this model, MTOCs are transported to the cell center and 'captured' by the chromatin once they reach the central chromatin mass. The process resembles models of biased 'search-and-capture' used to describe spindle assembly [7,36,37]. Here, we use the model to define the spatial properties of the bias of attraction by comparing the simulation outputs to spatial trends in MTOC motility seen in experiment.
Velocity. The Brownian component of the motion results from radial MTs in an aster interacting with motors in all directions resulting in uniform forces, which fluctuate due to thermal noise and microtubule dynamics [27]. The equation of motion of the particle is given by its velocity ( _ X) and net angle (θ net ). The magnitude of the velocity is described by: where D eff is the effective diffusion coefficient, v eff is the effective velocity and t is time. The values of D eff and v eff are determined by fits to experimental data as described in the section on data analysis and listed in Table 1.

Model geometry
Oocyte radius (r cell ) 40 μm [9] Chromatin radius (r chr ) 10 μm [9] RWD model Total simulation time (T) 8000 s [8,9] Step time (δt) 0.1 s Optimized for numerical accuracy No. of MTOCs (N p ) 100 [8,9] Effective diffusion coefficient (D eff ) 0.006 μm 2 /s Fit to experiment (S1A and S1B Direction. The net angle (θ net ) of motion, which governs the direction of motion, is determined by: where θ df is the diffusive and θ dr the directed angle. The angle θ df can take random values uniformly distributed between 0 and 2π, while θ dr is determined by the angle of the vector XC ! between the particle (X) and the center (C) of the cell (Fig 1B). Drift field. Sub-cellular transport has been successfully modeled in previous work by combining Brownian motion and directed transport [38]. In our model, the Brownian and directed components are spatially determined. Here, two radially symmetric fields of drift, one repulsive and one attractive, were modeled to represent the effective forces acting on the MTOCs. The attractive field resulting in centrosomal aster convergence in Xenopus oocytes towards chromatin [5], was modeled by a sigmoid gradient [7,27], A similar function was tested for the model of mouse oocytes. Since the attractive (0 a ðXÞ) and repulsive fields (0 r ðXÞ) determine the nature of motion, depending on the position of the particle in space, the net velocity of a particle ( _ X net ðx; yÞ) is given by: where _ X c ðx; yÞ is the velocity at the cell cortex determined by the weighted average of repulsion and Brownian motion. _ X c is determined by: In this model, attraction originates from the central chromatin, since the removal of the nucleus in previous experiments resulted in MTOCs losing their directionality of motion [9,39]. The repulsion from the cell cortex, is based on models of MT-aster pushing at the cell cortex as a mechanism of centering in multiple cell types [28,29]. MTOCs undergo Brownian motion when neither attraction nor repulsion are acting on them. The net angle θ net (X) is the circular weighted average [40] of θ df and y dr ðXC Þ, weighted similarly to the velocity magnitude (Eqs 3 and 4). The fields of 0 a ðXÞ and 0 r ðXÞ are modeled in 2D by radially symmetric functions dependent solely on the distance r from chromatin. A modified sigmoid gradient was tested, based on previous experimental measurements of long-range gradients acting on MTs around chromosomes in Xenopus spindle assembly [5][6][7] and theoretical models testing their form [25,26] and their role in centrosome directional motility [27].
The field (ϕ) of either attraction (ϕ a ) or repulsion (ϕ r ) is determined by: where r is the distance-for ϕ a , r = 0 at the chromatin edge and for ϕ r , r = 0 at the cell boundary. r 1/2 is the distance at which the function takes the half-maximal value and 1/s is the steepness factor. While exponential gradients were also tested, they failed to reproduce the spatial trend in experimental data and are not shown here. Representative gradients of ϕ a and ϕ r are plotted as scaled values between 0 and 1 ( Fig 1C).

Mechano-chemical model of MT-motor interactions
A mechanistically detailed model was developed using Cytosim [41], a C++ Langevin dynamics simulation engine, by building on previously developed models of MT-mechanics [14,42], polymerization kinetics [7,25] and motor interactions [20-24, 26, 27, 43]. To test what minimal components will produce the observed centripetal motility of MTOC asters, a model of an MT stabilization gradient described previously [27] was developed by mapping the gradient shape optimized in the RWD model. This scenario was compared with scenarios where the gradient consisted of immobilized minus-ended molecular motors. These biased 'search-and-capture' scenarios were contrasted with self-organized scenarios which lacked any directional bias, i.e. diffusible motor-complexes and MTOC pushing from the cell boundary. The model is implemented in an oocyte cell geometry with MT polymerization dynamics and mechanics as well as discrete stochastic molecular motors that are either immobilized or diffusible. Cell geometry. As before the oocyte was modeled in a 2D circular geometry with a concentric chromatin region. The chromatin region here is not treated as an absorber in this model, unlike in the RWD model.
MT polymerization dynamics and mechanics. MTs were modeled as discrete polymers undergoing dynamic instability-stochastic switching between growth and shrinkage phasesbased on the four-parameter model [44,45]. Since dynamic instability parameters of meiosis I mouse oocytes are not reported, values were taken from measurements made on centrosomal MTs in mitotic pig kidney cells, which have average MT lengths of 3.2 μm [46], very similar to the average length of MTs in mouse meiotic MTOCs (<L > *3 μm) [2,3,33]. The flexibility of the simulated MTs is defined by a combination of bending modulus (κ), typical cytoplasmic viscosity and thermal energy. The values of these parameters were taken from literature ( Table 2), from the mouse oocytes or related systems, if mouse data was not available.
MTOCs. MTOCs were modeled as hollow circular structures of radius 0.2 μm, initialized randomly throughout the cytoplasmic region. Each MTOC had N MT number of MTs uniformly distributed radially around the centrosome, forming an aster. MTOCs can move throughout the cell interior.
Motor mechanics. Molecular motors were modeled as described previously [27,42] as discrete particles with properties taken from minus-end directed dynein-like motors moving at a velocity v m . Motors bind stochastically to MTs at an attachment rate r attach only if the distance between them is less than a threshold distance of attachment (d attach ). Detachment occurs at a rate r detach . Motors bound to MTs behave like Hookean springs with a spring constant (k mot ), experiencing a extension force (f ex ) parallel to the filament if attached to a motile filament. The rate of detachment increases with extension as: r detach ¼ r 0 detach Á e jf ex j=f 0 (based on Kramers theory [47]), where r 0 detach is the constant basal detachment rate and f 0 stall force. A separate r end detach rate accounts for the rate at which motors detach from the end of a filament-here it is set to a value similar to r 0 detach . As f ex increases to values approaching f 0 , motor step sizes are expected to change, referred to as the 'gear-like' behavior of dynein [48]. This is modeled by a piece-wise approximation, as previously described [27,49]. The parameters of motor mechanics, when available, are taken from experimental measurements of dynein (Table 2). In the absence of reported values estimates were used. Immobilized motors bind to microtubules and generate forces on the MTs causing their movement. Conversely, diffusible motor complexes with diffusion coefficient D c , are modeled as crosslinkers [23], which can bind two different microtubules. Since aster have minus-ends at the center, cross-linking minus-end directed motility of the motors results in coalescence of the MTOCs.
Aster motility. MTOC asters move due to a net force calculated by resolving multiple forces acting on the MTs (using a finite difference scheme to solve the Langevin equation of motion in Cytosim [42]) which are: (a) f bend , arising due to bending of growing MT filaments contact the rigid cell boundary and scaled by the bending modulus of MTs (κ) [11,12,50]; (b) f i ex , the extension force generated when immobilized motors bound to MTs walk on the filament and are stretched by length x, experience a restoring force f i ex ¼ k mot Á x, resulting in filament motion [42]; (c) f c ex , the motor-complex of diffusible tetrameric motors have two binding sites separated by a spring of stiffness k mot . Complexes that are simultaneously bound to two MTs, undergo stretching. The restoring spring force brings the two MTs close and results in clustering [20,22,23]; (d) f diff , the diffusive force is a random normally distributed force [42]; and (e) f drag , the drag force acting on the MTs resulting from the translational and rotational drag forces on the individual points representing the filament [42] based on the cytoplasmic viscosity (η). Parameters are reported in Table 2.

Methods
Data analysis. Experimental 2D trajectories of multiple MTOCs in mouse oocytes meiosis I described previously by Schuh and Ellenberg [9] with Δt % 3 to 4 minutes) ranging between pre-and post-NEBD stages (data kindly shared by M. Schuh) were normalized by origin shift Table 2. MT-motor model parameters. The mecho-chemical parameters of motors were based on reported values for dynein, while motor densities were estimated. MT polymerization dynamics and mechanics is taken from literature, while cell geometry parameters are identical to Table 1.

Motor parameters
Motor stiffness (k mot ) 0.1 pN/nm [53] Motor (dynein) speed (v m ) 2 μm/s [54][55][56] Attachment rate (r attach ) 12 s −1 [27] Basal detachment rate (r to the respective end-points. Thus, multiple experimental trajectories could be compared. For comparison with experiment, simulated trajectories were down-sampled to a time-interval comparable to experiment (Δt = 3.5 minutes), and the following measures of motility were evaluated from both experimental and simulated data: Motility statistics: The instantaneous velocity v = δL/δt and directionality or tortuosity (χ) was estimated as χ = d net /L. A value χ = 1 indicates directional motion, while χ * 0 indicates a random (tortuous) path. The time taken by MTOCs from their time of nucleation to arrive at the edge of the central chromatin region is used to calculate a 'capture time' from both experiments and simulations. The mean square displacement (msd) of experimental and simulated trajectories was estimated by: where r is the displacement, t is the time point, δt is the time increment. A sliding window of δt from the smallest simulated time step to 3/4th of the trajectory length was used [27]. This empirically determined cut-off reduces artifacts arising from the fact that the number of steps of time-length greater than this threshold are extremely small [58,59].
In order to estimate input parameters for the RWD model (Eq 1), the msd profiles of experimentally measured MTOCs were fit to the model of diffusion and transport: On the other hand, the msd output from the motor-MT interaction model (Models section) were quantified by fitting to an 'anomalous diffusion' model: to estimate the apparent diffusion coefficient (D 0 ) and the anomaly parameter (α) as described previously [27] (S1 Text). All data fitting was performed using the trust region reflective least square fitting algorithm implemented in the Optimization Toolbox of MATLAB (Mathworks Inc, USA). Simulations. The RWD simulation time steps were optimized for numerical stability, by simulating a normal random walk (α = 1) with no boundary conditions and minimizing the error in input and fit diffusion coefficient (D). Simulations with typically 100 particles were run for 8 Á 10 3 seconds and took *20 minutes. The explicit MT-motor model (Cytosim) with 80 MTOCs was run for 1.2 Á 10 3 seconds, typically requiring *8 hours on a 12 core Intel Xeon machine with 16 GB RAM.
Model optimization to experiment. The parameters of the model were optimized by a rank minimization method based on a modified weighted root mean square error () of simulations as compared to experiment. Only one published dataset by Schuh and Ellenberg [9], has reported detailed XY over time trajectories of MTOC centering motility during mouse oocyte maturation, as a result of which quantitative model optimization was performed on this dataset (kindly provided by M. Schuh). The error ((k)) for a parameter set k was estimated between the i th experimentally measured (e i ) and simulated (s i ) values by: The weight (w i ) was determined based on whether the simulated value was within one σ of the average experimental measure or outside that range: where w m = 2, n i is the number of experimental data-points for the i th value and n max is the maximal number of data-points in experiment. This weighting scheme allows us account for the increased uncertainty in experimental measurements with fewer data points. The choice of w m was set empirically to 2 to ensure a high penalty (since n i /n max 1). For each parameter set the of multiple variables (v) were ranked (R v (k)). The directionality (χ) distribution with distance (r) and frequency distribution of capture times (t c ) were used to calculate errors. The sum rank, R s = S v R v (k) was minimized to obtain a parameter set, referred to as optimal.

RWD simulated MTOC motility compared to experiment
The experimental trajectories of MTOCs show a distinct centripetal motion as seen in the time-projected trajectories ( Fig 1A). The input parameters for diffusive and directed motion in the RWD model ( Fig 1B) were obtained from fitting Eq 7 to experimental msd profiles (S1A The outputs of such a gradient result in XY trajectories which are directed inwards at the cell boundary, random in the mid-zone and directed closer to chromatin (Fig 1D), qualitatively comparable to experiment (Fig 1A). Quantitative comparisons between experimental and simulated χ profile reflects this trend in motility-particles at the cell boundary and near chromatin are more directed, than those in the mid-zone ( Fig 1E). The simulated capture time distribution also matches with experiment ( Fig 1F). In order to understand the mechanism underlying the experimental motility observed, trajectories were further analyzed for their time-dependence.

Pushing and pulling profiles in MTOCs transport
The experimentally measured MTOCs have heterogeneous distance-time profiles, with some MTOCs moving rapidly in < 40 min, while others undergo a delayed (> 40 min) inward movement (Fig 2A). The optimized RWD model profiles qualitatively match those from experiment. In previous work, distance-time plots with a sigmoid profile have been interpreted to mean pulling forces are at work, while a parabolic has been interpreted to mean pushing is at play [15]. Here, we use a fit function with three parameters, n: a measure of the shape of the profile (n > 1: sigmoid and n 1: parabolic), T half : the time at which the distance travelled is half-maximal, and d max : the maximal distance travelled, as follows: The experimental and RWD simulation distance-time plots were fit to obtain d max , T half and n.
A value of n > 1 is taken to indicate pulling, while n 1 is taken to indicate pushing. Representative data from experiment ( Fig 2C) (Fig 2E). The n value from simulations is higher in the mid-cell region as compared to near chromatin or at the cell boundary. This can be understood in terms of the sharp transition in the attractive gradient of drift (ϕ a ) at r % 15 μm, resembling a 'pulling' process. The phenomenological model does not allow an interpretation of the origin of pushing or pulling forces. We therefore proceeded to add detailed molecular-motor and MT polymerization dynamics to the model to make experimentally testable predictions about the system.

Mechanistic models of MTOC centering: Gradient and Self-Organized Models
The RWD model predicts two drift fields-attractive from the center and repulsive from the cell boundary-are required to reproduce experimental statistics. Here, we proceed to test molecular mechanisms which combine molecular motors, MT-dynamics and mechanics of MTs and membrane interactions, with the aim of developing a mechano-chemical understanding of the effective drift fields. We categorize these mechanisms into two types: 1. Self-organized

Gradient-based
We then systematically evaluate plausible mechanisms and evaluate them for their ability to result in MTOCs finding the center of the oocyte within an experimentally observed time-scale (*20 min).
(i) Self-organized mechanisms. (a) Cortical pushing. In somatic cells the MT-asters find the center of a cell entirely based on pushing interactions at the cell cortex [32]. We test if pushing from the cell boundary alone, could result in MTOC centering, while maintaining uniform MT polymerization and immobilized dynein-like motors (Fig 3A). MTOCs at the start of the simulation are localized randomly throughout the oocyte (Fig 3B). However even after 20 minutes, asters did not move centripetally (Fig 3C, S1 Video).
(b) Clustering motors and cortical pushing. Based on the self-organized centripetal motility of MTOCs observed experimentally in mouse oocytes, clustering motor complexes had been proposed to be the major driver of MTOC motility to the center of the cell [9]. As a result, we modeled diffusible minus-end directed (dynein-like) motor complexes uniformly in the cytoplasmic region (Fig 3D and 3E). These complexes couple two nearby MT filaments and as a result of minus ended-motility, result in clustering of the asters (since asters have their minusends at the MTOC center). In this scenario, dynamic instability parameters were uniform and no immobilized motors were modeled. The complexes diffuse through the whole cell including the chromatin region within % 30 s, mimicking start of NEBD (S2 Video). However even after 20 min of simulations, only a small fraction of the simulated MTOCs are inside the chromatin region ( Fig 3F, Table 3).
(ii) Gradient-based mechanisms. (a) MT dynamic instability gradient. Based on evidence of a gradient of dynamic instability in multiple cell types [5,7,10,36], we tested a gradient of f cat and f res as a mechanism for centering of asters, while maintaining a uniform surface-immobilized motor distribution ( Fig 3G) and random MTOC nucleation (Fig 3H). However, this too did not result in a perceptible increase in accumulation of MTOCs at the chromatin center at the end of 20 minutes (Fig 3I, Table 3, S3 Video).
(b) Dynein motor gradient. It was only when motors were localized in a chromatin-centered sigmoid gradient and dynamic instability parameters were homogeneous (Fig 3J), did the randomly nucleated MTOCs (Fig 3K) dramatically converge to the center (Fig 3L, S4 Video).
Thus, a minimal model of asymmetric pulling forces resulting from a gradient of immobilized motors can move MTOCs to the center of the cell in a time-scale comparable to experiments (*20 min). However, a quantitative comparison of the simulation statistics with experiment is necessary to understand the critical parameters in this model.

Model sensitivity to stall-force and density of motors and MTs per aster
In the previous section, qualitatively, self-organized mechanisms of MTOC centering could not drive centripetal motility. The spatially binned directionality (χ) of simulated MTOCs in the absence of any gradient further quantifies this (Fig 4A). Neither tetrameric motor complexes, nor uniform surface immobilized motors and pushing at the cell boundary result in a trend in χ comparable to experiment. Cross-linking by motor complexes of different stall forces (f 0 = 2 and 7 pN) and densities (N c m ¼ 10 3 and 10 4 motors/oocyte) were tested and higher stall forces with high densities result in high values of χ throughout the cell (> 0.5). A directional bias in the form of a field of f cat and f res (based on Eq 5) resulting in asymmetric MT lengths, also fail to reproduce the directionality trends (Fig 4B). Confirming our qualitative observations from the simulation visualization, only a gradient of motors (f 0 = 7 pN, N i m ¼ 10 3 ) can reproduce most of trend in directionality as a function of distance (Fig 4C and 4D).
In the absence of experimental estimates of the number of motors and their stall forces from meiotic mouse oocytes, we explore two extreme values of f 0 (2 and 7 pN) reported in literature for dynein and scan N i m over thee orders of magnitude. We find high density (N i m ¼ 10 4 motors/oocyte) of weak motors (f 0 = 2 pN) ( Fig 4C) and a lower density (N i m ¼ 10 3 motors/ oocyte) of strong motors (f 0 = 7 pN) (Fig 4D), can both reproduce experimental profiles of directionality. The proportion of MTOCs captured at the chromatin boundary was evaluated by following the distance of MTOC centers as they entered the chromatin mass. Strikingly the proportion of MTOCs captured in the first 20 minutes of simulation from the motor gradient with immobilized motors of f 0 = 7 pN and N i m ¼ 10 3 most closely matched experiment (Table 3). However, the mean velocity values from simulations were insensitive to either f 0 = 2 or 7 pN and motor densities over a range N m = 10 2 to 10 4 motors/oocyte (Table 4).
Additionally to test how the motility was affected by the total number of MTs per aster in the scenario of a motor gradient, we examined the χ and <v> of simulated MTOCs which were initialized at a distance of 15 μm from the chromatin edge. As a result we expect these asters to experience the maximal force asymmetry, as they are at the lower end of the motor gradient. We find for increasing motor densities (N i m ) directionality χ continues to increase, but increasing N MT per aster appears to rapidly saturate the χ value for any given N i m value for both f 0 = 2 and 7 pN (Fig 5A and 5B). Increasing N i m leads to a marginal increase in the mean velocity (<v>) for a fixed N MT value. However increasing the value of MTs per aster, to our surprise, does not affect <v> (Fig 5C and 5D). We interpret this to be the result of the uniform radial distribution of MTs in the aster and a tug-of-war arising from it.

Model
Motor type Localization f 0 (pN) N m % captured In order to further understand the role of motor-density changes and their effect on MTOC motility, we evaluate the random walk statistics of the motility and compare it to experiment.

A hybrid model of a motor-gradient and self-organized clustering enhances MTOC centering
Simulations of multiple gradient forms, motor types and densities in this work suggest a motor gradient as a minimal model to understand the centering motion of MTOCs in    (Fig 4A), due to MTOC aggregation. We hypothesized that a hybrid mechanism combining an immobilized motor-gradient with diffusible clustering motor-complexes might reproduce the complete experimental distance-dependent directionality profile. Visually the MTOCs appear to find the center more efficiently (Fig 6A and 6B). Here 10 4 clustering motors per oocyte of stall force 2 pN were combined with 10 4 immobilized motors per oocyte and stall force 2 pN. A systematic screen of the effect of increasing tetrameric complex density while keeping the density of immobilized of motors constant was evaluated in terms of the measure of directionality. The radial distance profile of χ increases in the mid-range when clustering complex density is increased from 10 3 to 10 5 motors/oocyte (Fig 6C). Either 'weak' diffusible dynein-like complexes (f 0 = 2 pN) with a high density (N c m 10 5 ) or 'strong' motors  (Fig 6C and 6D). This result suggests that while clustering alone cannot center MTOCs in oocytes, it improves the fit to experiment.
Thus we believe our model reproduces both the qualitative and quantitative nature of the MTOC centering motility and provides novel insights into the sensitivity of this model. It demonstrates how a combination of directional cues and self-organized clustering can center small MT asters in a large cell such as an oocyte.

Discussion
Meiotic spindle assembly in mammalian cells in the absence of centrosomes involves the nucleation of MTOCs in cytoplasm and their coalescence and sorting around chromosomes, resulting in bipolar spindle assembly. The nucleation of MTOCs in cytoplasm and the centripetal motility of small radial MT asters has been previously observed during the first meiotic division in mouse oocytes [3]. Similar convergence was also observed in Drosophila oocytes [60] and quantitative analysis and model calculations were used to infer that directed transport for the MTOCs was essential for the coalescence of MTOCs [61]. In Drosophila the unconventional kinesin Ncd (minus-end directed) was implicated to play a role in this inward motility [62]. In mouse oocyte MTOC coalescence, cytoplasmic dyneins or comparable minus-end directed motor have been implicated in force generation for the centripetal movement [9]. Here, we analyze experimental data from the early stages of the meiotic maturation of mouse oocytes. The experimental analysis is used to constrain a field of spatially inhomogeneous drift in a random walk. The optimized model of drift is used to model gradients of the motor dynein and MT dynamic instability. By comparing the simulation outputs with experiments, we arrive at a minimal model of a gradient of motors, essential to reproduce the experimentally observed statistics.
The experimental statistics of the MTOC motility from mouse oocytes are mostly representative of post-NEBD dynamics, allowing us make the simplifying assumption that the nuclear envelope plays no explicit role in the process. The movement of these radial MT arrays appears visually to have both an effectively diffusive and a transport component (Fig 1A). The frequency distribution of velocity from experiment is long tailed and fit to a lognormal function (S4 Fig), suggesting anomalous super-diffusive transport. While the motility had been previously described as 'stop-and-go' [9], we find little evidence of 'stop' or pause events in the motility. Our quantification of MTOC motility in cells, demonstrates the motility is qualitatively comparable to previous estimates of centrosomal aster movement observed in C. elegans fertilization [15], MTOCs in Drosophila oocyte meiosis I [61] and centrosomal asters in Xenopus meiotic extracts [5,27]. This suggests a common theme underlying the transport of radial MT arrays in meiotic spindle assembly.
The distance travelled or displacement from the start-point plotted over time of MT arrays have been used previously [15] to distinguish between "pulling" and "pushing" modes of motility of centrosomal MT arrays during pronuclear migration. In the case of mouse meiotic MTOCs, these plots also help to distinctly separate trajectories into two sub-populationsthose which show an initial rapid rise followed by capture (< 40 min), while others do not move much for a long time (> 40 min) and after this delay are captured at chromatin (Fig  2A). We further improve on the work of Kimura et al. (2005) [15] by fitting the data with a saturation model with cooperativity (Eq 11), and use it quantitatively to distinguish between pushing (n * 1) and pulling (n > 1) mechanisms (Fig 2C). Every trajectory from experiment was fit to obtain a profile shape (sigmoid or parabolic) term n (S3 Fig) and comparing n from simulation and experiment shows a pushing mode of transport (n * 1) close to chromatin and the cell boundary, while those in the mid-range are pulled (n > 1) (Fig 2E). The distance travelled plots from the mechanistic motor-gradient model with MT asters sorted by nucleation distance demonstrate an increase after a short delay near chromatin (d n = 0 to 10 μm) due to pulling (S5A Fig). Those in the intermediate range (d n = 10 to 20 μm) have a longer delay and then appear to be pulled (S5B Fig). Close to the cell-boundary (d n > 20 μm) they increase rapidly and then saturate, due to MT-pushing at the membrane (S5C Fig). Thus a model of a gradient of chromatin-centered motors and a rigid cortex can reproduce the qualitative differences observed in the distance-time profiles of MTOC transport. However in these experiments, the absence of clear 'pulling' in experiment (n > 1) near chromatin might have been missed, since the MTOC identity is lost as it nears the chromatin. While cell cortexbased pushing mechanisms have been demonstrated to generate centripetal movement [12,50], the nature and localization of minus-end directed pulling motors in the oocyte remains to be determined. Recent evidence of from mouse oocytes during MTOC fragmentation [63] and meiotic maturation [34] suggest dynein anchored at the nuclear envelope might influence both processes. A careful study of dynein localization dynamics in this and related systems, could by used to test our model prediction.
A 'tug-of-war' in the transport of anti-parallel MTs moving on a surface coated with motors arises from the action of the same species of motor acting against each other, with small asymmetries in length, amplifying the velocity of transport [64]. Two-fold length asymmetries (10 to 20 μm) have been previously observed in centrosomal asters [5,7] and simulations of such asymmetric asters on sheets of dynein motors resulted in aster transport towards chromosomes [27]. However, the mouse MTOC radius is in the range of 2 to 3 μm and no appreciable asymmetry of lengths has been reported [9]. In this work, neither homogeneous motor distributions (S1 Video), nor tetrameric motor-complexes (S2 Video), nor a dynamic instability gradient (S3 Video) result in convergence to chromatin in the *20 min time scale seen in experiment. Taken together, a gradient of motors is necessary in a minimal model with a mouse oocyte geometry for MTOCs to converge to the chromatin center (S4 and S5 Videos). This suggests aster motility in meiosis I of mouse oocytes differs from that in meiosis II of Xenopus oocytes. In the latter, length asymmetry can result in directional transport, however in mouse oocytes the MTOC asters are *4 fold shorter in MT length, resulting in smaller forces being generated. The motor numbers would then be insufficient to successfully resolve the tug-of-war in mouse oocytes by a simple *2-fold changes in MT length, as seen in Xenopus oocytes. In contrast, the gradient of immobilized motors result in a comparable spatial distribution of velocity for minus-ended motors with f 0 = 2 pN (S6A Fig) and f 0 = 7 pN (S6B Fig). Spatial gradients of RanGTP [6,65] are thought to direct centrosomal MT asters to chromatin during spindle assembly in meiotic extracts [5][6][7]. However the outcome of our simulations predicts that in meiotic MTOC motility, motor gradients are necessary. Such a gradient could arise from self-organized diffusion and attachment of minus-ended motors on MTs nucleated at the chromatin periphery. In addition, the motors could be immobilized on intracellular organelles, as shown in centrosomal aster centration in C. elegans [16]. It remains to be seen which of these specific mechanisms results in the the chromatin centered gradient of anchored motors. Some evidence from experiments in mouse oocytes serves to support our motor gradient model, where pre-NEBD maturing oocytes showed a high-concentration of dynein motors around the nucleus [66]. More recently, evidence of dynein localized at the nuclear periphery fragmenting MTOCs [63], suggests a validation of our model of a chromatin centered gradient of dynein-like motors. Further testing is however required to examine its dynamics.
The role of chromosomes in the centering activity of MTOCs in maturing mouse oocytes during meiotic I spindle assembly had been tested by experimentally removing the nucleus [9,39]. The MTOCs of such oocytes have a scattered appearance and fail to assemble bipolar spindles. While chromosomes are considered essential for meiotic spindle assembly [67], we hypothesize that they also serve as the primary guidance cue of centripetal MTOC motility, in a manner comparable to other aster cell centering systems [5,27,36]. To test this, the spatial organization of MTOCs in enucleated oocytes was measured using published experimental data reported by Schuh and Ellenberg [9]. The image analysis of MTOCs positions from experiment (S7A and S7B Fig) resulted in a radial density distribution around the cell center (S7C Fig). This distribution is comparable to a simulation of randomly localized MTOCs, suggesting the chromatin serves as an important guidance cue for the centripetal motility of MTOCs, and in it's absence the directional motility is lost. The molecular mechanism which converts the positional information of the chromatin into a gradient of molecular motors, still remains to be understood.
The msd profiles of MTOCs from simulations that were close to chromatin transition from super-diffusive to sub-diffusive motility (S1 Text, S8 Fig), as estimated by α, the measure of anomalous diffusion (S8F Fig). This transition results purely from changing the single dynein motor stall force from that reported for bovine brain cytoplasmic dynein (f 0 = 2 pN) [48] to the yeast cytoplasmic dynein (f 0 = 7 pN) [56] value. The MTOC convergence in the 2 pN stall force calculations is slow (S4 Video) as compared to when the motors had a higher stall force (7 pN) (S5 Video) at the same motor density (N i m ¼ 10 3 motors/oocyte). Once the MTOCs reach the center of the motor gradient, they undergo effectively sub-diffusive movement, unable to generate enough force asymmetry to escape. With an additional centrosomal fluorescence label in mouse MTOC motility, the intra-nuclear motility and reorganization of MTOCs could be studied in future. This will also help better understand the subsequent bipolar spindle formation and its connection with chromosome biorientation in meiosis I of mouse oocytes [68,69].
The sensitivity analysis of the gradient model demonstrates that MT number per aster does not affect velocity and directionality of the asters. On the other hand the model is sensitive to motor density. The motor density in 2D area-density ranging between N i m ¼ 10 3 and 10 4 motors/oocyte corresponds to a physical density of % 0.2 to 2 motors/μm 2 for an oocyte of radius 40 μm. The addition of motor complexes, which cross-link MTs and walk to the minusends of the respective MTs, also change the dynamics of MTOC transport. At high densities of motor complexes, the asters coalesce to a few (*10) clusters and centripetal convergence results. Yet the directionality profiles are qualitatively different from the experimental measures ( Fig 6C). While the densities of both immobilized and diffusible motors tested are mean area densities, over-expression of dynein motor proteins could serve as an experimental test of our model predictions.
A limitation of the models described here is that both the phenomenological and mechanistic models assume a 2D circular geometry, while the mouse oocyte is a 3D sphere. So asters finding the cell center in 3D is expected to take longer due to dimensionality in biological search problems [70], and a further parameter optimization would be required. Our effort here serves to reduce the search for 'scenarios' , i.e. combinations of mechanisms such as motors, MT-dynamic instability, gradients and clustering motors, by optimizing simulations to experimental data. In future, a full 3D model would then only require parameter optimization to a 3D experimental dataset. At the other extreme, a 1D model could have further simplified the geometry based on the radial-symmetry of the system, as has been assumed in the 'slide and cluster' model of linear MT filaments in spindle assembly [71]. However this would ignore the orthogonal interactions between MTs of multiple asters, which are only possible in an explicit 2D geometry. Thus our choice of spatial dimensions is driven by an attempt to capture the important qualitative behavior of the system, while keeping the model simple enough for clarity and calculation speed. While a 3D simulations of cellular processes are more 'complete' , it has been suggested the choice of spatial dimensions should simplify the system to sufficiently capture the spatial dynamics [35,72]. A further limitation of both the phenomenological RWD and MT-motor model is the availability of only one time-series dataset. In future, additional experiments with mouse oocytes could help test the model predictions. While this model has been developed to test the generality of mechanisms on the mouse MTOC motility, it would be useful in future to explore the relevance of this model to other aster-centering systems such as in C. elegans [15], Xenopus [27], Drosophila [73] and sea urchin [74]. Additional mechanisms such as a contractile actin network has been shown to drive spindle assembly in meiosis I of starfish oocytes [18], but are ignored in this model since in experiments, it does not affect the process [9]. Additionally the potential role of MT-dependent MT nucleation [75] remains to be explored as a self-organized mechanism. Testing alternative cellular geometries and additional mechanisms could in future result in a more complete exploration of MT aster centripetal transport and the physical constraints. This in turn might help us better understand the physical basis of the evolutionary diversity of aster-centering mechanisms.
The models of of MTOC centripetal motility explored in this study, provide insights into several aspects of the early stages of self-organized spindle assembly. For instance, the RWD model demonstrates that for a cell of diameter % 80 μm, a simple random-walk strategy of MTOC particles 'search-and-capture' at the chromosomes is insufficient within the time-scale of spindle assembly in meiosis (20-30 min) [9], and a long-range bias in motility is essential. The mechanistic model demonstrates that a bias of motors is minimally capable of reproducing the dynamics. MT mass increase which could arise due to increasing MT lengths (stabilization) [5,6] or number (nucleation) [75], both appear to be insufficient to affect the dynamics of motility. We also find, dynein-like clustering motor complexes, which diffuse and cross-link MTs, at a high density result in an aggregation of MTOC asters to the cell interior. Such a selforganized mechanism however is inefficient in resulting in MTOC capture percentages comparable to experimental values (Table 3) in comparable time-scales (% 20 min). Additionally MT density per aster does not affect MTOC centripetal motion but the density of immobilized motors (N i m ) localized in a gradient does change the dynamics. Such parameters are likely to be adjusted by cells, depending on the specific geometry and time-constraints. The results of this study could be used to predict the nature of MTOC centripetal transport, when the size ratios of asters to cell sizes are also comparable to the mouse oocyte. A more complete picture of the evolutionary constraints on mechanisms driving radial MT arrays to find the center of a cell will however require further quantitative studies from acentrosomal spindle assembly from more organisms, as well as calculations that explore a greater parameter space.
The model presented here predicts the functional form of a drift field in agreement with experimental data of MTOC motility. A mechanistic model of the gradient with immobilized minus-ended motors minimally reproduces experimental dynamics and allows us to test the effect of single molecule characteristics of motors on collective transport statistics. A hybrid model with a gradient of immobilized motors and diffusible clustering dynein-like complexes, combined with cortical pushing (Fig 7), satisfies all experimental measures available. Measuring the localization, density and mobility of dyneins in the mouse oocyte during meiosis I, would be a useful test of the model predictions. While the model is fit to spindle assembly during meiosis I in mouse oocytes, it also predicts the design constraints in terms of MT lengths and the localization, mobility and density of molecular motors when MT radial arrays are required to search space and undergo capture at the cell center. All those parameter sets (k) which fell within the top 15% of the sum of scores ranking are plotted with the component rank for capture time (R t (k)) (grey) and directionality (R χ (k)) (blue). The parameter sets (k) of the top-ten ranks are listed in S1 Table with the values of r 1/2 and s for ϕ a and ϕ r . For a representative subset of the optimization scheme, the error () in (B) χ and (C) t c were evaluated keeping the attractive gradient constant (r a 1=2 ¼ 10 μm and s a = 1) and varying the repulsive gradient parameters r r 1=2 (y-axis) and s r 2 (x-axis). The colorbar indicates the value of . MTOCs (grey) were simulated in cytoplasmic space of the 2D oocyte geometry marked by the outer cell boundary (outer blue circle) and inner chromatin region (blue circle) in the presence of minus-end directed motor complexes (purple dots) initialized in the cytoplasm with density N c m ¼ 10 3 motors/oocyte. These diffusible motors with stall force f 0 = 7 pN, can cross-link MTs and result in clustering by walking towards the minus-ends of neighboring asters that they crosslink. (MP4) S3 Video. MTOC motility in a gradient of MT dynamic instability. 80 MTOCs (grey) were simulated in the 2D oocyte geometry with the outer cell boundary (outer blue circle) and N i m ¼ 10 3 uniformly distributed surface immobilized motors with f 0 = 7 pN. The f cat and f res parameters were distributed in a sigmoid gradient (Fig 3G) originating from the center of the chromatin region (inner blue circle). (MP4) S4 Video. MTOC motility in a gradient of weak motors. 80 MTOCs (grey) were simulated in the 2D oocyte geometry with the outer cell boundary (outer blue circle) and N i m ¼ 10 3 surface immobilized motors with f 0 = 2 pN distributed in a sigmoid gradient (Fig 3J) originating from the center of the chromatin region (inner blue circle), and homogeneous dynamic instability. (MP4) S5 Video. MTOC motility in a gradient of strong motors. 80 MTOCs (grey) were simulated in the 2D oocyte geometry with the outer cell boundary (outer blue circle) and N i m ¼ 10 3 surface immobilized motors with f 0 = 7 pN, distributed in a sigmoid gradient (Fig 3J) originating from the center of the chromatin region (inner blue circle). (MP4) S6 Video. MTOC motility in a gradient of immobilized motors and diffusible motor-complexes. 80 MTOCs (grey) were simulated in the 2D oocyte geometry with the outer cell boundary (outer blue circle) and motors which are immobilized at a density of N i m ¼ 10 3 motors/ oocyte (green dots) in a sigmoid gradient originating from the center of the chromatin region (inner blue circle). The diffusible minus-end directed motor-complexes with N c m ¼ 10 4 motors/oocyte (purple dots) bind to 2 MTs and walk simultaneously on them, generating a clustering force on the MTOC asters. For both kinds of motors f 0 = 7 pN. (MP4) S1 Table. Optimized RWD gradient parameters. The top ten sum ranks are listed in ascending order, based on the ranks from the directionality and capture time error ((k)), with the corresponding parameters of the attractive (r a 1=2 and s a ) and repulsive (r r 1=2 and s r ) gradients (see