A full lifespan model of vertebrate lens growth

The mathematical determinants of vertebrate organ growth have yet to be elucidated fully. Here, we utilized empirical measurements and a dynamic branching process-based model to examine the growth of a simple organ system, the mouse lens, from E14.5 until the end of life. Our stochastic model used difference equations to model immigration and emigration between zones of the lens epithelium and included some deterministic elements, such as cellular footprint area. We found that the epithelial cell cycle was shortened significantly in the embryo, facilitating the rapid growth that marks early lens development. As development progressed, epithelial cell division becomes non-uniform and four zones, each with a characteristic proliferation rate, could be discerned. Adjustment of two model parameters, proliferation rate and rate of change in cellular footprint area, was sufficient to specify all growth trajectories. Modelling suggested that the direction of cellular migration across zonal boundaries was sensitive to footprint area, a phenomenon that may isolate specific cell populations. Model runs consisted of more than 1000 iterations, in each of which the stochastic behaviour of thousands of cells was followed. Nevertheless, sequential runs were almost superimposable. This remarkable degree of precision was attributed, in part, to the presence of non-mitotic flanking regions, which constituted a path by which epithelial cells could escape the growth process. Spatial modelling suggested that clonal clusters of about 50 cells are produced during migration and that transit times lengthen significantly at later stages, findings with implications for the formation of certain types of cataract.


Introduction
An important and still largely unresolved issue in developmental biology is the nature of the processes that specify the size and shape of organs. development in Drosophila have been particularly informative, helping identify underlying signalling networks [1] and the potential role of mechanical feedback [2] in the growth process. The pathways that regulate growth of the Drosophila imaginal disc are also present in higher organisms, where they presumably play analogous roles. However, modelling organ growth in vertebrates is a daunting prospect because of the size and complexity of the structures involved.
The lens of the vertebrate eye offers an opportunity to model the growth of a simple vertebrate organ across the entire lifespan and, by doing so, identify key mathematical determinants of the growth process. From a modelling standpoint, the lens has several advantages. Its role in image formation requires a smooth ellipsoidal shape. It contains only two cell types: epithelial cells and fibre cells. The lens cell population (10 5 -10 6 cells) is sizeable, but certainly accessible using modern computing tools. The prismatic fibre cells that make up the majority of the lens volume are packed closely together, leaving little or no space between. Importantly, fibre cells do not turnover; all the fibres that differentiate in the course of development are retained in the body of the lens.
We previously quantified the distribution of proliferating cells on the spherical anterior lens surface [3,4] and used those data to formulate a first generation, branching process model of lens growth [5]. Using an expanded, dynamic version of that model, we here report to our knowledge, the first full lifespan growth model for a simple vertebrate organ system, the ocular lens. We were able to follow lens growth through more than 1000 iterative cycles, during which lens volume increases more than 4000fold. Remarkably, for an organ whose development appears to depend on a stochastic growth engine, the variance in the process was much smaller than predicted by the cell power law [6]. Furthermore, we found that modest adjustment of just two parameters, the rate of increase in cellular footprint area and the proliferation rate, was sufficient to capture the entire growth behaviour of the lens, including radial increment, zonal organization and patterns of cellular immigration-emigration.

Methods
2.1. Age-dependent growth parameters 2.1.1. S-phase labelling Mice (C57BL/6 J) were obtained from Jackson Laboratory (Bar Harbor, ME). S-phase cells were identified following incorporation of 5-ethynyl-2'-deoxyuridine (EdU; Invitrogen, Carlsbad, CA, USA), as described [3]. EdU was administered by intraperitoneal injection and mice were killed 1 h later by CO 2 inhalation. Eyes were fixed in 4% paraformaldehyde/phosphate-buffered saline, embedded in paraffin, and sectioned (4 µm) in the midsagittal plane. EdU-positive cells were visualized using Click-iT (Invitrogen) chemistry with Draq5 (Cell Signaling Technology, Danvers, MA, USA) as a nuclear counterstain. Three sections from each of three lenses were used for each time point.

Measurement of radial growth
Radii of intact, fixed embryonic and early postnatal lenses were determined from digital images. Measurements were supplemented by published data collected from adult mouse lenses [7,8].

Determination of fibre cell dimensions in the equatorial plane
Fibre cell width (w) was determined from the spacing of nuclei in the meridional rows that form during the early stages of fibre cell differentiation. These measurements were made on confocal image stacks collected in the course of an earlier study [4]. Fibre cell thickness (ρ) was determined from vibratome sections of fixed lenses stained with ActinGreen (Invitrogen) to label the sub-membrane actin cytoskeleton and thus highlight the membrane profiles. Measurements were made on fibre cells located 10 cell layers below the equatorial lens surface.

Development of a dynamic model
2.2.1. Anatomy of the growth process, proliferative zones and the applicability of a stochastic model The anatomy of the lens growth process is well understood (figure 1). Mitosis in the anterior epithelium results in displacement of cells from the epithelial margin [9], the 'penny pusher' effect [4]. Stimulated by growth factors present in the vitreous humour [10], displaced cells differentiate into fibre cells and are incorporated into the lens body, increasing its volume and surface area. Mitosis is restricted to the epithelium but the rate of cell division varies greatly with latitude. The epithelium can be subdivided into four zones with respect to the local rate of cell proliferation [4]. From the anterior pole to the equator these are: the central zone (CZ), pre-germinative zone (PGZ), germinative zone (GZ) and transition zone (TZ). In adult lenses, the CZ contains mitotically quiescent cells, the PGZ and GZ contain proliferating cells (the rate of cell proliferation is highest in the latter), and the TZ contains cells that have withdrawn from the cell cycle but yet to differentiate. Cells move toward the equator, traversing the various zones before exiting the TZ and differentiating into fibre cells. Of note, apoptotic cells are rarely if ever detected in healthy lenses [4].
The young lens grows quickly, fuelled by rapid cell division in the PGZ and GZ. In the 30 days following its formation on embryonic day 10.5 (E10.5), the mouse lens increases in volume by more than 4000-fold. It transforms from a hollow ball of ≈1200 cells, into a fully formed lens; a solid, ellipsoidal structure containing ≈200 000 cells (≈50 000 epithelial cells and ≈150 000 fibre cells). Initially, all epithelial cells are dividing, but from E14 onwards the proliferation index drops below 100% and probabilistic models are applicable. Here, we modelled lens growth from E14.5 onwards.

Model derivation
We previously developed a branching process-based model to examine how the geometry of the lens constrains its growth [5]. We developed formulae relating proliferation in the various zones to the movement of cells across zonal borders and, ultimately, their deposition as differentiated fibres in the body of the lens. The model was tested over the period from one to three months of age, under the simplifying assumption that age-dependent variables remained unchanged over this relatively short interval. In actuality, almost all parameters vary significantly with age and an adequate full lifespan model must account for this. In this study, therefore, we incorporated newly acquired data on age-dependent changes in fibre cell cross-section, cellular compaction, proliferation rate, cell-cycle period, cellular footprint area (a) and the proportional contribution of each zone to the epithelial surface (η). Model simulations were benchmarked against data on lens size, cell proliferation rates and published measurements of epithelial population dynamics [4]. Derivation of the static model has been described [5]. The same axioms were employed here in formulating the dynamic full lifespan model (appendix A). Cell dynamics in the proliferative zones are defined by the following conditions: where X i (t) is the number of cells at time t in zone i, i = 1 ( CZ), 2(PGZ), 3(GZ), 4(TZ), u = the part of the number of cells that move from TZ into the fibre compartment; which consists of elongated cells with hexagonal equatorial plane intersections given by the length ρ(t) and width w(t) of the hexagons. These equations lead to a linear system with a unique set of solutions (appendix A), providing us with the number of cells that move between zones.

Formation of proliferative zones
The fraction of the anterior lens surface covered by each zone (η) was determined previously: η CZ = 0.32, η PGZ = 0.33, η GZ = 0.26, η TZ = 0.09 [5]. In adult lenses, many epithelial cells are in G 0 [11] but this is not the case in the embryo, where Ki67 labelling (a marker for cycling cells) approaches 100% [12]. S-phase cells are distributed uniformly across the embryonic lens epithelium [12] which, as a result, can be considered to represent a single GZ (i.e. η GZ = 1). To delineate the zones in an embryonic lens, the latitudinal distribution of S-phase cells was visualized in lens sections from EdU-treated embryos (figure 2). S-phase cells were numerous throughout the epithelium at E14.5 (figure 2a,b) and the labelling index correspondingly high, with values approaching 50% across much of the epithelium (figure 2d). By postnatal day 14 (P14), however, the labelling index at all latitudes had declined, and S-phase cells were concentrated in the periphery (figure 2c). The labelling index profile (figure 2d) indicated that three zones (TZ, GZ and PGZ) were present. The fourth zone, the CZ, formed between P14 and P21 (figure 2e). The cross-sectional area of newly formed fibre cells is an important determinant of radial growth. Previously, we assumed that width and thickness of the prismatic fibre cells were constant over time. To more accurately follow the growth process across the lifespan, we used confocal microscopy to measure the size and shape of fibres in the equatorial plane ( figure 3). These experiments revealed that the width (w) of the fibres changed significantly during embryonic and early postnatal lens development before stabilizing by ≈6 months of age. By contrast, fibre cell thickness (ρ) was relatively constant across the lifespan.  For the examples shown (E14.5 and P14), the proliferation rate is higher throughout the epithelium at the earlier stage. By P14, three zones (TZ, GZ and PGZ) are discernible, the fourth (CZ) appears a few days later. The fractional contribution (η) of each zone to the surface area varies initially but, by P28, values have stabilized (e). Scale bar (a-c) = 100 µm.

Epithelial cells
Unlike other ectodermally derived structures (e.g. epidermis), differentiated lens cells are not shed from the tissue surface. Instead, they are incorporated into the body of the lens, increasing its volume and surface area. Thus, ongoing cell division in the proliferative zones of the epithelium causes the area of the zones to increase (as the surface area of the lens increases). In such a 'self-inflating' system, the number of epithelial cells must be regulated tightly if positive feedback loops are to be avoided. We modelled epithelial cell dynamics using input data on zonal proliferation rates and cellular footprint areas (a). Data were collected from E14.5 until the end of the lifespan. Proliferation rates and derived η values are presented in table 1 and figure 2d,e. Corrections for fibre cell compaction were applied from six months, as described below. Simulations were based solely on these restricted input data and a random number generator used to mimic branching mechanisms for individual cells. Variations in cross-sectional area (the product of width (w, blue diamond) and thickness (ρ, orange square)) were incorporated into the lens growth model. Note that fibre cell width increases significantly over the first few months of age, whereas cell thickness is relatively constant across the lifespan.
Modelling suggests that the epithelial cell population rises sharply from E14.5-P14 (figure 4), reaching a maximum at ≈4 weeks (w) of age. The rapid increase is a consequence of shortened cell cycles at embryonic stages (appendix B), a uniformly high proliferation rate (figure 2) and relatively small a at all latitudes (which serves to maximize zonal populations). The epithelial cell population overshoots at 4 w, when, for a brief period, it exceeds 50 000. The subsequent population decrease (to ≈43 000 cells by 12 w) reflects a reduction in zonal proliferation rates and an increase in a for all zones. Surprisingly, the growth trajectory of the epithelial population depended simply on the relative rates of change in proliferation and footprint area (appendix C). From 12 w onward, the epithelial population remained constant, despite continuing production of cells in the GZ and PGZ and associated macroscopic growth of the lens. Of note, model runs were nearly superimposable, the stochastic nature of the underlying growth mechanism notwithstanding.
Model simulations were compared with experimental determinations of postnatal lens epithelial populations [4] and newly acquired data on embryonic lens epithelia (table 2). Despite some scatter, the empirical data were consistent with model predictions. Both simulations and measurements suggested a transient overshoot in cell population at ≈4 w declining, thereafter, to a smaller but more stable population.

Fibre cells
The lens growth engine (cell division in the proliferative zones) is particularly active during embryonic and early postnatal life leading, at peak, to the deposition of more than 15 000 fibre cells per day (figure 5a). However, by three months, modelling indicates that the rate of fibre production declines to a few 100 cells per day (figure 5a, inset). The reduction reflects a decline in PGZ and GZ proliferation rates, coupled with an increase in a at all latitudes that serves to deplete the proliferative zones of cells. The slow but steady addition of fibres over the life of the mouse (almost 4 years in our laboratory) results in a substantial accumulation of cells. By the end of the lifespan, modelling suggests that the lens contains 400 000-500 000 fibres (figure 5b).

Fibre cell compaction
The growth model consistently overestimated the lens radius at later ages ( figure 6). For example, by the end of the lifespan, empirical measurements suggest that the lens radius is ≈1300 µm, whereas, a value        of ≈1750 µm was predicted by the model. The size of newly formed fibres is known ( figure 3) but there is little information on the fate of pre-existing cells located far below the lens surface. In the original model formulation [5], we assumed that mature fibre cells have a stable morphology. However, multiple lines of evidence (reviewed in [16]) indicate that over time cells in the lens core become compacted. If true, this provides a plausible explanation for the increasing discrepancy between model predictions and empirical measurements at later stages. Here, we used linear compaction factors of 0.030% and 0.025% to correct for loss of fibre volume in the periods six months to 2 years, and 2 to 4 years, respectively (figure 6 and appendix A). The rate of compaction is slow, but integrated over the lifespan, the effect is sizable. Expressed volumetrically, the uncompacted tissue would be twofold to threefold larger than actually observed. Thus, the lens would be expected to physically overgrow the eye in the absence of fibre cell compaction.

Cellular flow and flow reversal
The lens growth process is essentially one of cellular displacement. Cells introduced into the epithelium through mitosis in the proliferative zones force the egress of cells at the margin [4]. As a result, cells generally flow from higher to lower latitudes. Localized flow reversals can occur, however, in response to changes in zonal proliferation rate or a. Flow reversals can have the effect of isolating populations of epithelial cells and preventing their progeny from entering the fibre cell compartment.  To see how this operates, consider the relationships between CZ, PGZ and GZ. The proliferative zones, PGZ and GZ, produce a surfeit of cells, sufficient to fuel their own expansion (which occurs passively, in response to lens volume increase) while providing enough cells for export. By contrast, the CZ is mitotically inactive. Thus, expansion of the CZ in response to lens volume increase requires either spreading of constituent CZ cells or immigration of cells from the neighbouring mitotically active PGZ.
We modelled the stochastic flow of cells across the CZ/PGZ border ( figure 7). The flow was positive initially (i.e. cells tended to move from PGZ -> CZ) but, by three months, the net flow was zero, indicating that the rate of CZ cell spreading matched the rate of zonal expansion. Published data (appendix D) show that the footprint areas (a) of cells in the PGZ, GZ and TZ increase monotonically over the lifespan. By contrast, from 3 to 12 months, a CZ decreases. To examine the consequences of such a decrease, lens growth was modelled under conditions where a CZ was reduced in the period 3-12 months (CZ cell contraction case) or allowed to increase in parallel with the other zones (CZ cell expansion case, appendix D). The CZ cell population was modestly decreased in the expansion case because fewer cells could be accommodated within the zone but other population and growth parameters were largely unaffected (figure 7a). Although the zonal populations were not strongly influenced by changes in a CZ , there was a striking effect on the direction of flow at the CZ/PGZ border. If a CZ was allowed to increase in parallel with the other zones, then the net daily flow of cells across the border was zero (figure 7b, lower panel). However, driven by stochastic fluctuations, each day a small number (less than 5) of CZ cells are expected to cross the border into the PGZ (and vice versa). Under conditions where a CZ did not increase (i.e. the normal condition), the flow of cells was always strongly positive ( figure 7b, upper panel)  consequence of the reduction in a CZ during the period from 3 to 12 months was to temporarily corral cells in the central epithelium.
To examine the necessary conditions for cellular isolation, suppose that cells to be corralled are gathered in a zone which covers an area (or a volume) A(t) at time t, and contains X(t) cells of individual size a(t) with offspring branching mechanism Z. At time t, X(t) cells produce a total of O(t) offspring cells. It follows that at time t + 1 there will be no cellular emigration if and only if Considering this condition 'in the mean' we obtain Assuming that the system also grows 'in the mean', i.e. that the right-hand side is positive, one sufficient condition to ensure 'isolation' is to require This condition depends solely on cellular proliferation (distribution of Z) and individual cell sizes (a(t)). If a(t + 1) > a(t), we need to lower proliferation rates sufficiently to satisfy the equation above. If, as in our case, there is no proliferation (i.e. E(Z) = 1), then cells can be held in the quarantine zone by dropping the growth rate of a(t) to zero (i. e. to have a(t + 1) = a(t)), which is what happens in the CZ in the period 3-12 months.

Achieving precision with a stochastic growth engine
Simulations suggest that the ocular lens, a precisely defined biological structure, can be generated using a stochastic engine (figures 4 and 6). Even after more than 1000 iterations, each involving tens of thousands cells acting randomly and independently, the precision of the lens growth process is striking. How might this arise? There appear to be two key elements. First, before E14, the system grows deterministically (but simply, because cells appear to multiply as fast as they can). After E14, we have a fully stochastic mechanism, but by that phase the lens epithelium contains more than 10 000 cells. That number is already sufficiently large for the law of large numbers to operate. This alone would confer a relatively high degree of precision but the effect works in tandem with a second phenomenon. Azevedo & Leroi [6] have shown that the coefficient of variation (CV), where CV = standard deviation mean , is a feature of organ growth [6]. Suppose that an organ is built using a standard branching process (X(t)), where and (Z i ) is a family of iid integer-valued random variables whose distributions are given by the branching mechanism If we start with N cells, it is not difficult to calculate CV at time t; In our model, the proliferative zones (GZ and PGZ) are flanked by mitotically quiescent regions (CZ and the fibre cell compartment) offering a pathway by which supernumerary cells can escape the growth process. This arrangement improves the precision of the process further, by introducing a correction term to the Azevedo-Leroi CV (for calculations and mathematical details see appendix E). The effect of non-mitotic zones is that there exists a factor, say λ, such that emigration can be represented as (1 − λ)(O(t) − X(t)). We consider a constant H := λE(Z) + (1 − λ); observe 1 < H < E(Z) in our situation. The CV for such a process, say CV λ t is given by It follows that In our case, we typically have E(Z) close to 1, implying that the factor E(Z)/E(Z) Hence, a zonal structure improves precision with factor √ λ. For the lens process, λ changes with time and is mostly of the order 1/10. It follows that if a non-zonal arrangement were to achieve (via the law of large numbers, for example) a relative precision of ±5%, then the introduction of a zonal structure would further improve the precision to around ±1%. This calculation applies to any organ-building stochastic system where an active growth region is flanked by non-mitotic zones into which excess cells may escape. will be gradually displaced towards the lens equator, slowly traversing the PGZ and then the GZ. Cells probably accelerate as they draw nearer the equator, owing to the increasing cumulative 'push' imparted by mitosis occurring at higher latitudes [4]. Cells are likely to undergo multiple cell divisions as they cross the PGZ and GZ. Indeed, lineage tracing experiments have identified clusters of clonally related cells near the epithelial margin [17]. But what is the size distribution of emergent clusters, how will this depend on age, and how long does it take for clusters to complete their migration and enter the body of the lens?

Cluster size and transit time
In appendix F, we develop a model that traces the path of an epithelial cell as it travels from higher to lower latitudes. We obtain the following basic formula (formula (F 9) from appendix F; where L 0 is the starting layer of the travelling cell, L is the final layer, N ∞ is the travelling time and p k is the proliferation rate during the kth day): After N cycles the original cell, which started on layer level L, would on average reach level L N . Using the Markov property, we obtain Because we know the number of levels, say T, then by solving the equation L N = T, we obtain an estimate for the average number of cycles, N, required to reach the equator. Interestingly, the average cluster size does not change if we change p k values (see formula (F 8) in appendix F); for the mouse lens, we predict clusters consisting of around 50 cells on average. The intuition is that although a higher proliferation rate produces more new cells sooner, it also speeds the cluster in the direction of the equator. Although clone size is not sensitive to proliferation rate, transit time is (see formula (F 10) in appendix F).

Conclusion
Previously, we developed a stochastic branching process-based model of lens growth [5]. The model was formulated from a simple set of axioms and we used it to explore how the unique geometry of the lens influences its growth. We restricted the model to a narrow time frame, so that we could make the simplifying assumption that growth parameters remained constant. In the present study, we modelled growth over the entire lifespan. This necessitated the incorporation of new empirical data on parameters such as age-dependent changes in cellular proliferation rate, fibre cell dimensions and cell-cycle duration. Output of the model was consistent with experimental measurements of cell population dynamics and radial growth [4] although, at later time points, a radial correction factor had to be introduced to compensate for the presumed (but unproven) effects of fibre cell compaction.

Three distinct phases of lens growth
Based on our modelling efforts, we can now distinguish three growth phases: early (from lens formation on E10.5 to eyes opening at P14), middle (from two weeks to six months of age) and late (from six months to the end of life).
Early growth is explosive. In mice, eyes open at two weeks, so the lens must be built quickly. This requires the production of a large number of cells in a short period. Consequently, all cells are cycling and the cycle is short. A temporary capillary network, the tunica vasculosa lentis, envelops the embryonic lens and this may be necessary to support rapid growth. The footprints of individual epithelial cells are small (perhaps because there is insufficient time between cycles for robust cellular growth). As a result, proliferation zones are densely populated, maximizing cell production. Any shortfall in cell production at this stage has a permanent detrimental effect on lens size.
In the middle period, the zonal structure is well developed but there is a risk of overproduction of cells because the growth in epithelial surface area (with a constant η) means that more and more cells are located within the proliferative zones. Exponential overgrowth is avoided by a generalized cell spreading phenomenon, which restricts population increase in the GZ and PGZ and consequently serves as a powerful governor on lens growth. The movement of cells from the CZ into the fibre cell compartment is blocked temporarily by a decrease in a CZ , which reverses the stochastic flow at the CZ/PGZ border, leaving CZ cells marooned in the pupil space. Migrating cells undergo multiple divisions as they traverse the proliferation zones, leading to the formation of clonal clusters. Migration times vary from 50 days at the beginning of this period to more than 400 days by the end.
Late lens growth is characterized by the slow but steady deposition of fibre cells at the lens surface. Cell spreading is completed, but declining proliferation rates in the GZ and PGZ limit the possibility of exponential growth. Late phase growth lasts for several years (at least in laboratory settings) and although the addition of cells is slow, significant radial growth over this extended period is expected. Here, the phenomenon of fibre compaction is important, contributing to the internal refractive properties of the lens [16] and preventing the lens from physically overgrowing the eyeball.

Proliferation rate and a are key determinants of lens growth
An important finding from this work is that lens growth dynamics can be modelled adequately by merely adjusting proliferation rates and footprint areas of individual epithelial cells. Although this effect is not readily apparent in equations wherein all four zones are combined, it can be appreciated in the simplified case, where only one zone exists (this corresponds, for example, to the situation at E14.5, when S-phase cells are distributed uniformly throughout the epithelium). The formulae are provided in appendix C (see equation (C 2)), where we show that the average number of cells in the epithelium will increase, remain the same, or decrease, based solely on whether r · 2ρw > a, or r · 2ρw = a, or r · 2ρw < a, where r := r(t) is the proliferation rate and a := a(t + 1) − a(t) is the footprint area increment. Because the cross-sectional size of the fibre cells (ρw) is constant from six months onward (figure 3), it is noteworthy that such a 'macroscopic effect' (i.e. total number of cells in the epithelium) is determined by the 'microscopic features' r and a. Conceptually, we envisage the contributions of r and a as two pedals in an automobile, one serving as the accelerator, the other as the brake. Thus, the entire growth process is achieved by gently applying pressure on one or other pedal. It is not yet clear whether the 'pedals' have biological identities, but they could be analogous to the Hippo and TOR pathways which, in other organ systems, regulate organ growth by controlling cell number and cell size, respectively [18]. In support of this notion, it was recently shown that inactivation of PI3 K, an upstream regulator of mTOR, results in a lens that is too small but otherwise histologically normal [19].
In our model, cell density and local cell proliferation rates (key determinants of lens growth) are treated as independent variables. In the absence of evidence to the contrary, this seems the most parsimonious approach. The model also assumes that cells flow towards the edge of the epithelium in response to mitotic pressure generated through cell division at higher latitudes. There is some experimental support for this view. For example, in hypophysectomized frogs, lens epithelial cell division is completely blocked. This results in the immediate cessation of fibre cell formation, as expected if cell division were the primary driver of cell migration [9]. We note, however, that other authors [20] have measured cell density distributions across the lens epithelium and related the shape of the density curves to the production and egress of cells (the 'ordered pull-through model' (OPT)). In the OPT model, density and proliferation are not treated as independent parameters. Rather, cell density is proposed to be responsive to both a 'push', provided by epithelial cell production, and a 'pull', resulting from the flow of cells into the fibre compartment. The two lens growth models are not mutually exclusive and, in either case, further experimental data, perhaps on the nature of forces acting on the migrating cells, will be required to evaluate them properly.

Absence of 'catch-up' growth
Typically, the growth of organisms and organs ceases at maturity. Once adulthood is reached, cell proliferation rates are reduced to balance cell death rates. In that sense, lens growth is unusual because despite the absence of cell death, cells continue to be produced until the very end of life [21]. The eyeball itself stops growing relatively early in postnatal development [22]. Thus, a structure with indeterminate growth (the lens) is housed within a receptacle (the eye ball) with a determinate growth mode. Perhaps because growth is indeterminate, the lens shows no sign of 'catch up growth'. If the growth zones of the lens are not adequately populated early in development, the resulting growth deficit appears to be permanent. For example, inactivation of the gap junction gene Gja8 causes a temporary reduction in epithelial cell proliferation in P2 mice [23]. Although brief, the interruption is sufficient to cause a lifelong growth deficit. A similar effect was noted following conditional deletion of the p110α subunit of phosphoinositide 3-kinase [19]. We replicated this effect with our model ( figure 8) and found that a two-day hiatus in proliferation is particularly consequential if it occurs in the young lens. Cases of microphakia are not uncommon in ophthalmology but there are no reported cases of macrophakia (i.e. an overly large lens), implying that either growth rates are already maximal or that the lens actively compensates for cellular overproduction (for example, by adjusting the fibre cell compaction factor). In mice, lens growth defects often lead to microphthalmia [19,23] suggesting that at least with regard to growth, the eye takes its cue from the lens.

Clones, transit times and cataract
Our calculations indicate that cells are likely to undergo multiple cell divisions as they traverse the proliferation zones. The calculations suppose that the movement of cells is vectorial, and that there are no cellular eddies that might cause a cell to dwell at a particular latitude (details are given in appendix F). The expected maximal clone size of less than 100 cells is consistent with the results of lineage tracing experiments [17]. A counterintuitive finding was that cluster size was not influenced strongly by proliferation rates and that, as a result, clusters of similar size are expected to form in young lenses (where proliferation rates are high) and old lenses (where rates are low). By contrast, proliferation rate had a powerful effect on transit time. For example, progeny of a PGZ cell located close to the CZ/PGZ border at P14 would take 7 weeks to differentiate into fibre cells. For cells at the equivalent location in a six-month-old lens, the expected transit time is more than 13 months.
Deep sequencing studies have recently identified somatic sequence variants in DNA from human lens epithelial cells. The measured allele frequencies imply the presence of millimetre-sized clones [24]. Sunlight-exposed cells are prone to somatic mutations and UV-B exposure is a risk factor for cataract [25]. Cellular flow from the epithelium into the body of the lens provides a route by which mutations generated in central, sun-exposed regions might be amplified and conveyed to the fibre cell mass. If mutations were to occur in genes necessary for fibre cell transparency, they might manifest as spoke-like opacities in the body of the lens (each spoke corresponding to a cluster of mutant fibres derived from a single, sun-exposed epithelial antecedent). In the mouse, cells appear to be temporarily quarantined in the central epithelium ( figure 6) and only released at stage when the projected transit time exceeds the remaining lifespan. It will be interesting to determine if similar mechanisms operate in humans and whether they serve to suppress cataract formation until late in the lifespan. Acknowledgements. The authors thank Dr Alicia De Maria for her help with confocal imaging.

Appendix A. Model description
The axioms underlying the model have been described [5]. In this paper, we allow for cellular dynamics to change with time and also introduce compaction factors which act on the lens late in the lifespan. Consequently, the full lifespan model is more complex, requiring us to adjust the corresponding formulae accordingly.
Consider some time t in the life of a mouse. The lens epithelium has a surface area of A(t) µm 2 and (in principle) is divided into four zones, where index 1 corresponds to CZ, index 2 to PGZ, index 3 to GZ, index 4 to TZ. The surface areas of each of the zones are presented via formulae A i (t) = η i (t)A(t), i = 1, 2, 3, 4. Observe that we allow for some 'i' to have η i (t) = 0(⇔ A i (t) = 0) ; meaning simply that up to time t, the ith zone has not yet developed. The number of cells in zone 'i' is X i (t), and the surface area covered by a single cell is a i (t) Hence, A i (t) = X i (t)a i (t). Let us describe what happens between t and t + 1.
The nature of the TZ is such that cells within TZ do not divide. In principle, we allow for all other cells to divide, remain the same or die in the other three zones, i.e. we consider three probability distributions of the form and then develop the model based on the specifics of actual data on Z i (t). In the case of the mouse lens we have, p 1 1 (t) ≡ 1, p 2 0 (t) ≡ 0 ≡ p 3 0 (t) (see also [4] and [5]). Hence, the information about Z i (t) is in essence q(t) := p 2 2 (t) and r(t) := p 3 2 (t) (observe that p 2 1 (t) = 1 − q(t) and p 3 1 (t) = 1 − r(t)). Assuming that (Z 2,n (t): n ∈ N) are iid copies of Z 2 (t) and (Z 3,n (t): n ∈ N) are iid copies of Z 3 (t), we obtain the total number of offspring cells in the PGZ and GZ between t and t + 1: The interior of the lens consists of elongated fibre cells with hexagonal equatorial plane intersections given by the length ρ(t) and width w(t) of hexagons [5]. In the full lifespan model, we also allow for a compaction effect. Because the compaction phenomenon is not yet well understood, we assume only that it has some global effect on the radius of the lens R(t). The measure of that effect is given by the compaction factor = the percentage of radius that has contracted (0 ≤ < 1, but in general is either 0 or very small ≈10 −4 ). The consequence is that the offspring cells should not only cover the increase in growth but also the effect of the compaction process. Let us denote by v the number of cells required to offset the effect of the compaction. Based on the geometry of the lens we obtain that We denote by x = the number of cells that move from PGZ to CZ, y = the number of cells that move from PGZ to GZ, z = the number of cells that move from GZ to TZ, u = the number of cells that move from TZ to the fibre compartment and are needed to maintain the increase in growth. Observe that the total number of cells that move from the epithelium to the fibre compartment is u + v. The basic immigration-emigration process now follows the one described previously [5]: As v is given in (A 1), we again have four equations with four unknowns, but the dynamics are now different, because we allow for changes in time and have the additional factor . Nevertheless, it is not difficult to solve the equations, and we obtain the following solutions. Observe first that, for every i = 1, 2, 3, 4, we have It follows that, for every i = 1, 2, 3, 4, Using (A 1), we now obtain the complete solution of (A 2): where E i (t) and F i (t) are deterministic and given by and .

Appendix B. Cell-cycle estimates
In adult lenses, cell-cycle time is approximately 24 h, with a 12 h S-phase [13]. However, a simple assessment of the increase in cell number during embryonic lens development indicates that the lens cell-cycle duration must be shorter at earlier stages. For example, at E10.5, the lens is a hollow vesicle consisting of some 1200 cells. Of these, only half (those covering the anterior lens surface) are mitotically active after E12.5. By E14.5, we estimate that the epithelium alone contains 11 000 cells (and the interior almost twice that many). Clearly, even if all the original 600 anterior cells were actively cycling throughout the intervening period, it would be impossible to generate the necessary number of cells in the four 24 h cycles between E10.5 and E14.5.
To obtain estimates of cell-cycle lengths in the embryonic lens, we present a simple model relating the observed proliferation rates of cells to their cell cycle. Consider a set of cells within an organism, observed at a particular stage in the life of the organism and at a particular position within the organism. First, we postulate that, at this stage and at this position, the set consists of cells, which are 'cell-cycle homogeneous', in the following sense. Postulate 1. All cells within the observed set have a cell cycle of the same length. Measured in hours, the length of the cell cycle is C. Postulate 2. All cells within the observed set have the same 'probability of proliferation within C hours'; denoted by 0 < p < 1. More precisely, if a cell is selected randomly from the set, then the probability that the cell will start dividing within the next C hours is p. Postulate 3. There exists a fixed number 0 < α < 1, such that the S-phase for each cell within the observed set lasts αC hours. Second, we assume that we can perform an 'S-phase labelling' experiment on such a 'cell cycle homogeneous' set. Thus, for t hours, we expose the set of cells to a label that will be taken up by any cell, which is in the S-phase during any part of the t-hour period. We also assume that our experiment can be performed on a sample of 'cell cycle homogeneous' sets (with the same C, p, α for all elements of the sample).
Third, our experimental procedure should satisfy the following assumptions: (a1) Relative frequencies of the number of labelled cells converge to a fixed number 0 < r < 1.
(a3) The starting times of the S-phases are distributed uniformly within a time interval of C hours.
Remark. (a1) is a standard statistical assumption on the stability of relative frequencies; without it most statistical inference studies would lack meaning. (a2) ensures that the t-interval is short enough that our analysis remains 'within a single cycle'. (a3) is a standard initial assumption if we have no additional information about cycle properties.
Using our postulates and assumptions, we conclude that r is with respect to αC + t as p is with respect to C, i.e. we obtain our basic formula Observe that (B 1) and (a2) imply and Simply stated, (B 3) and (B 4) show that although it is mathematically possible to consider r to be close either to αp or to p, it is not realistic to have such relationships in real-life situations.
Observe that there are five variables in (B 1), so having some knowledge about any of them increases the usefulness of formula (B 1). Typically, we expect to know t, because we control the experiment. In our experiments, we always have t = 1 h. In many cases, the literature on cell-cycle properties provides evidence strong enough to estimate α. In our studies, we have assumed that from E14 on in the life of mouse α = 1/2 (we think that our evidence shows that between E10 and E14 this assumption does not hold, but we intend to devote a separate article to the details of this very early stage of the lens development).
Under assumptions t = 1 and α = 1/2, formula (A 1) leads to the following: and Observe that the estimate for quantity r is given via direct measurements of the results of our experiment.
Let us take the example of the E14-E16 time period to illustrate how to estimate C in actual situations. Consider a realistic estimate of 6 h ≤ C ≤ 24 h (for smaller or larger C, the results are even further from actual measurements). Our measurements give us r ≈ 0.48 at E14.5. Using (B 5), we obtain that p (now as a function of C) gives and we shall go with the 'worst case scenario' and take p = 0.75 for C = 6 h and p = 1 for C = 24 h . Observe that E14-E16 consists of 48 h , so with C = 6 we follow eight cycles, while with C = 24 we have two cycles. We employ our model as developed [5]. Our data show that at this stage the entire epithelium acts as a single GZ with a proliferation rate p, that a single cell covers on average 30 µm 2 , and that equatorial hexagonal intersections have an area of 7 µm 2 . Let us provide a brief sketch of average changes according to our model. Consider a single cycle and denote by X 0 the initial number of epithelial cells and by X 1 the number of epithelial cells at the end of the cycle (we use expected values, i.e. 'averages' for this calculation). Following our model, we have that (epithelial area) = 2 (equatorial area). Hence, if we denote by X ∞ the number of cells that leave the epithelium we obtain which leads to This implies that in a single cycle the radius of the lens will on average change with the factor For p = 0.75, this gives the factor 1.1129 and after eight cycles we obtain the radius estimate at E16 to be around 525 µm (instead of measured value 350 µm). For p = 1, the factor (B 8) is 1.1481 and after two cycles we obtain the E16 radius of 294 µm (instead of measured value 350 µm).
Hence, we end up with C ≈ 12 h (i.e. four cycles within E14-E16 time period). We interpolate linearly for the values of r in between (we have detailed measurements for r at E14.5 and E16.5) and use (B 5) to obtain p.
Appendix C. The relationship between r (proliferation rate) and a (change in the area of an individual cell footprint) determines whether the epithelial population will increase or decrease.
During embryonic development, the lens epithelium acts as a single GZ. During this phase, it is easier to appreciate the relationships between key parameters that regulate the growth (on average) of the number of epithelial cells. The same effect occurs within the four zone structure, but there the formulae are more complicated. We consider a single cycle time period from t to t + 1 and introduce the following notation: Therefore, we can explicitly express m(t + 1) in terms of m(t): Hence, Our data show that during most of the prenatal period ρw ≈ 7μm 2 , which implies that the 'threshold point' between r and a is 14r = a.
This shows how proliferation rate (denoted here by r) and change in the surface area of individual epithelial cells ( a), acting together may regulate the growth process. If 14r ≈ a, then the growth process will maintain a stable number of cells in the epithelium. If 14r > a, then the number of cells in the epithelium will increase. However, if we enlarge a so that it is bigger than 14r, there will be a decrease in the number of epithelial cells.
As an illustration, observe the time period around the birth of the mouse. At E20 we have r = 0.23, i.e. 14r = 3.22. This means that we need to have a < 3.22 µm 2 if we want to keep the growth of X(t). Interestingly, it is after P1 that we observe a more significant change in a (it starts to get bigger).

Consider a random variable O (O stands for offspring) defined by
Using Wald's technique, i.e. splitting into partition sets (X = n), n ∈ N ∪ {0}, one obtains several moment identities and and Suppose now that we have two sums of type (E 2) and that 'counting variables' are heavily dependent. More precisely, consider Using Wald's technique, it is easy to see that Single Zone Epithelium. We treat this case separately for two reasons. First, this is the case which describes the prenatal mouse lens. Second, owing to the simplicity of the formulae involved, it is easier to appreciate the contribution in this case than in the general case.
Suppose that the epithelium consists of X cells with footprint area a and a proliferation law i.e. E(Z) = 1 + p > 1. A single cell cycle is modelled via a family {X; Z i : i ∈ N} of independent random variables with Z i ∼ Z for every i ∈ N. Suppose that u ≥ 0 is the number of epithelial cells that turn into fibre cells with hexagonal intersections given by ρ and w dimensions. If u = 0, then we have a standard branching process dynamic and the corresponding CV is exactly as described (E 6). If u > 0, then we have a process with emigration, which corresponds to our model. In both cases we start with the total number of 'offspring cells', given by If we denote by X the increment in the number of epithelial cells, then the basic equation is Consider the subcase u = 0, first. Using (E 3) and (E 4), we obtain and Var(X + X) = Var(Z) E(X) + [E(Z)] 2 Var(X).
Therefore, the single cycle CV, under u = 0 assumption, is For a process of this kind having t steps (t cell-cycles here) and starting with X 0 ≡ N cells, we obtain and Var This leads to the standard branching process CV at time t formula, i.e.
Let us compare this to the case u > 0, which corresponds to our model. The geometry of the sphere leads to i.e. to Introducing the quantity H := λE(Z) The corresponding (single cycle) CV is In order to compare (E 10) and (E 19), i.e. standard single-cycle CV with 'emigration-corrected' CV, observe the factor λ/H 1/E(Z) in the prenatal period it ranges somewhere between 1/2 and 1/3. Typical behaviour with time is that 2ρw increases moderately, p decreases significantly and a increases rather sharply. For example, for 2ρw = 40, p = 0.1 and a = 200, the factor in (E 20) would be 11/61 ≈ 1/6. For stable parameters through t cycles we would obtain and Hence, we obtain It follows that for λ 'big' (say λ > 1/2 ; which is atypical for our process) the ratio However, for a mouse older than a month, E(Z) goes typically close to 1, while λ would roughly change from 1/5 to 1/7. For E(Z) ≈ 1, we have This principle essentially does not change even within our full model. However, owing to a four zone structure, the formulae are more difficult to follow. In the following section, we provide the details. Four Zone Structure. We follow our complete model description as given in appendix A. For i = 1, 2, 3, 4 we denote by X i the increment at time t in zone i. The equations of our model that provide the basic dynamics between the zones are (see appendix A) where O 2 are offspring cells developed in PGZ with branching distribution and O 3 are offspring cells developed in GZ with branching distribution and u is the number of cells that leave epithelium. This leads to a basic equation of the form where X = X 1 + X 2 + X 3 + X 4 . We consider the case u > 0 here. Denote by G the deterministic function defined by it follows that

X(t) = G(t)A(t). (E 27)
As A = 2ρ(t)w(t) · u, (E 26) leads to Let us introduce deterministic functions λ = λ(t) and µ = µ(t) by and It is now not difficult to check that (E 28) has a form analogous to (E 17), i.e.
As, for i = 2, 3, where p 2 = q and p 3 = r, we obtain In order to calculate E((X+ X) 2 ), we need formulae of the form (E 7). Observe that X i (t) = v i (t)X(t) in practical applications (i.e. there is an integer-valued function h). However, because X is of the order 10 4 , we shall neglect (as we described in [5]) the error smaller than 10 −4 and use the formula that E(X · ν i (t)X(t) ) = E(ν i (t) · X 2 ) = ν i (t)E(X 2 ). Therefore, using (E 7) and analogous formulae, we obtain and This enables us to calculate E((X + X) 2 ) and Var(X + X). Using the notationH := λ(qν 2 + rν 3 + 1) This leads to the formula for CV λ,μ in this case which is analogous to (E 19).

Appendix F. Cluster size estimates and transit time
We provide elements of the spatial movement of individual cells and their progeny through the epithelium. The analysis is simplified on several accounts and is a study of averages only. Nevertheless, it gives us the essence of the laws that characterize such movements through the epithelium. We consider the epithelium of the mouse lens as a half-sphere of radius R. During early stages, the entire lens epithelium acts as a single GZ, and is covered by a single layer of cells. On average each cell covers the surface area a (given in square micrometres). Furthermore, we assume that the epithelium is covered by L ∈ N concentric latitudinal layers (lines of latitude on the surface of the Earth are a good analogy), where every layer has the same width, equal to √ a. Hence, As the zonal structure of the epithelium does not have too much of an impact on (F1) until P7, we obtain the following Between P7 and P21, the zonal structure of the lens develops, where PGZ and GZ are zones with active proliferation. By P21, the number of layers in both zones stabilizes at around 50 layers.
The second set of simplifications covers the description of the spatial movement of a cell and its daughter cells. Suppose that we observe a single epithelial cell within a zone with proliferation rate p ∈ 0, 1. We are interested in two quantities, the position of the cell and the total number of its daughter cells, i.e. a stochastic process (Z t , Y t ). The interpretation of Z t ∈ N is the number of layers above (i.e. 'to the north of') the cell at time t ∈ N 0 . Suppose that we have m ∈ N layers above our cell. We assume that the probability of each layer 'pushing' our cell for a single layer 'southward' is p. According to our general model, cells within the epithelium act independently, which implies that layers act independently (because they are disjoint subfamilies of an independent family of random variables). The consequence is that the 'total push' of m layers acts as a binomial random variable B(m, p). Hence, the probability that at time t + 1 the cell will be in a position m + k, 0 ≤ k ≤ m is It follows that the average movement 'southward' in a single time interval, starting from a position m, is mp.
( F 2 ) As the cell started at m, at t + 1 it will be, on average, at m + mp = m(1 + p).

(F 3)
Owing to the independence assumption, the process (Z t ) has a (discrete time and discrete state space) Markov property. It follows then (using (F 3)) that, if the cell was initially positioned at Z 0 ≡ m, and proliferation rates in the next N time intervals were p 1 , p 2 , . . . ,p N , then, on average, the position of the cell as time N is The interpretation of Y t ∈ N ∪ {0} = N 0 is the total number of daughter cells at time t. Hence Y 0 ≡ 1. It follows that (Y t ) is a standard branching process with the branching mechanism distribution at t = k given by The well-known formula for the average value is (1 + p k ).

(F 6)
We are interested in a 'cluster of cells' produced by the original cell. By 'cluster' we mean all daughter cells that remain in 'close proximity' throughout their movement within the epithelium. Therefore, the worst case scenario is and we shall use (F 7) as the 'cluster size' (it is actually an average upper bound on the 'cluster size'). Suppose that our cell starts at some level L 0 and, on average, after N time intervals (for us a single time interval corresponds to a single cell cycle; after E18 this is a period of 24 h) our cluster reaches some level L ∞ , where L 0 < L ∞ . Then the cluster size does not depend on proliferation rates; it depends only on the ratio of the levels, i.e. we have The result follows easily from our formulae This result is perhaps unexpected, but the intuition behind it is that higher proliferation rates will increase the cluster size faster, but at the same time the cluster will travel 'southward' at a higher speed.
Our primary interest is in the cluster size that will enter the fibre compartment, i.e. in the case L ∞ = L. Because for most of a mouse life L is generally around 100 (with a high probability the initial position of interest is going to be L 0 = 2, L 0 = 3).
Another interesting consequence of our computations is that, although cluster size is not affected by the changes in proliferation rate, travel time is. Let us denote by N ∞ the time (number of cycles, later number of days) required for the average position to reach the equator. It follows that N ∞ satisfies the  At a later stage, when the proliferation rate stabilizes, we obtain the formula N ∞ = lnL/L 0 ln( 1 + p) .
(F 10) Let us illustrate the importance of this result on a particular example. Suppose that a CZ cell enters PGZ on the level L 0 = 2, and then travels, say L = 100, to the equator. If that occurs in a three-month-old mouse, average p will be 4%, i.e. N ∞ = ln 50 ln(1.04) ≈ 100 days.
Hence, we can expect a cluster of 50 cells to have entered the fibre compartment by the time the mouse is about seven months old. If, on the other hand, the migration process commenced in a six-month-old mouse, the average p will be 1%, i.e. N ∞ = ln 50 ln(1.01) ≈ 395 days.
Thus, the mouse would be over a year and a half old by the time the cell cluster arrived at the lens equator.
Among the most detailed information, we possess are data on proliferation rates. This enables us to calculate exactly the key factor N k=1 (1 + p k ). This factor appears in most of our formulae. For example, if we start at t = 0 which corresponds to E14.5, then we obtain the following This implies that every cell with L 0 ≥ 3 would have enough average time to leave epithelium by P1 (of course, that does not mean that it will actually do so; the outcomes of this process are stochastic).
If the cell travels within the PGZ only, the formulae above apply directly. However, when the cell (or its cluster) moves into GZ, we have to account for the differing proliferation rates of the two zones. Let us denote the rate in PGZ by q and the rate in GZ by r. We always have 0 < q < r < 1. Suppose that PGZ has the total of :L PGZ layers and that our cluster is positioned so that there are m layers within GZ above it. Then, the probability that at the same time PGZ layers produce a 'push' of size j and GZ layers a push of size k is L PGZ j q j (1 − q) L PGZ −j · m k r k (1 − r) m−k .
Hence the average movement 'southward' is qL PGZ + rm, and the average position at the end of the cycle is L PGZ + m + qL PGZ + rm = (L PGZ + m)(1 + r − ), (F 11) where 'the correction factor' is given by Observe that the factor (1 + r − ) is not going to cancel with a corresponding factor in (E 6) (which is of the form (1 + r)). In other words, formula (E 8) will be corrected via a finite number of factors of the We can estimate based on our data. Using midpoint estimate for m ≈ 1/2 L PGZ , we obtain ≈ 2/3(r − q). The ratio of r and q changes from say, a one-month-old mouse, when q = 1/3r to a 3+ year-old mouse, when q = 1/5r. Hence, for a ≈one-month-old mouse (r = 0.06) the factor in (F 13) is approximately 1 + r 1 + (5/9)r ≈ 1.0258, while for a ≈3-year-old mouse r = 0.01 it is 1 + r 1 + (7/15)r ≈ 1.00531.
The factor 1.00531 is very close to 1 and over a period of a few weeks will not have a significant effect on the cluster size. However, in the unlikely case (even under laboratory conditions) that a 3-year-old mouse lives to be 4 years old, the effect becomes significant. More precisely, 1.00531 365 ≈ 6.91.
This suggests that an expected cluster of, say size 100, may turn out to be of the size close to 700, which would cover most of the lens equator.