Advance Publication by J-STAGE Journal of Biomechanical Science and Engineering

The cerebellum has a unique morphology characterized by fine folds called folia. During cerebellar morphogenesis, folia formation (foliation) proceeds with granule cell (GC) proliferation in an external granular layer, and subsequent cell migration to an internal granular layer (IGL). GC migration is guided along Bergmann glial (BG) fibers, whose orientation depends on the deformation of cerebellar tissue during folia formation. The aim of this study is to investigate the contribution of the fiber-guided GC migration on folia formation from a mechanical viewpoint. Based on a continuum mechanics model of cerebellar tissue deformation and GC dynamics, we simulated foliation process caused by GC proliferation and migration. By changing migration speeds, we showed that the fiber-guided GC migration caused the non-uniform accumulation of GCs and folia lengthening. Furthermore, the simulation of impaired GC migration under pathological conditions, where GCs did not migrate along BG fibers, revealed that fiber-guided GC migration was necessary for folia lengthening. These simulation results successfully recapitulated the features of physiological and pathological foliation processes and validated the mechanisms that guidance of GC migration by BG fibers causes folia lengthening accompanied by non-uniform IGL. Our computational approach will help us understand biological and physical morphogenesis mechanisms, facilitated by interactions between cellular activities and tissue behaviors.


Introduction
The cerebellum has a unique morphology, comprising fine folds called folia. This morphology is associated with cerebellar function, e.g., motor control and cognition, based on evidence that a dysfunctional cerebellum has morphological defects (Sillitoe and Joyner, 2007;Leto et al., 2016). Therefore, to understand formation mechanisms underpinning the functional cerebellum, it is vital to investigate folia formation during cerebellar morphogenesis.
During cerebellar development, cerebellar morphology is formed via foliation mechanisms including folia formation, elongation, and bifurcation. This foliation process is driven by cellular activities including cell proliferation and migration (Corrales et al., 2006;Legue et al., 2016;Leto et al., 2016;Govek et al., 2018;Lejeune et al., 2019). At initial foliation stages, granule cell (GC) precursors proliferate in the external granular layer (EGL) (Sudarov and Joyner, 2007;Legue et al., 2015;Lawton et al., 2019). These GCs extend parallel fibers along a mediolateral axis, and they then migrate inwards along Bergmann glia (BG) fibers (Edmondson and Hatten, 1987;Yacubova and Komuro, 2002a). After migration Hironori TAKEDA* , **, Yoshitaka KAMEO* , ** , ***, Takahiro YAMAGUCHI* , **, Kazunori NAKAJIMA**** and Taiji ADACHI* , ** , *** through the EGL, molecular layer (ML), and Purkinje cell layer (PCL), GCs halt, and accumulate to generate an internal granular layer (IGL). Consequently, the IGL becomes thicker around a ridge (lobule), when compared to a folium groove (fissure). Additionally, the height of lobule increases relatively as the fissure position is fixed, leading to the formation of characteristic cerebellar morphology (Sudarov and Joyner, 2007). Fiber-guided GC migration is believed to be crucial for foliation mechanisms, as impaired GC migration and disorganized BG fibers cause abnormal foliation, including the decrease in lobule height (Men et al., 2015;Ryan et al., 2017;Hughes et al., 2020). A previous study proposed a mechanism of folia lengthening (increase in lobule height), suggesting that fiber-guided migration caused IGL with non-uniform GC accumulation, and subsequent folia outgrowth (Sudarov and Joyner, 2007). This study suggests that BG fibers radiate inward in a fan shape around fissure depending on the tissue deformation. Since the GCs migrate along BG fibers, GCs spread to larger area in the IGL. Conversely, around lobule, BG fibers radiate outward and GCs concentrate to smaller area in the IGL. Therefore, the oriented BG fibers depending on tissue deformation may cause non-uniform GC accumulation and associated tissue growth in the IGL, which cause folia lengthening.
Based on this foliation mechanism, we investigated the contribution of GC migration along BG fibers toward folia lengthening, using computational simulations. For foliation simulation, a continuum mechanics model of tissue deformation and cell migration was developed with reference to a previous cerebral morphogenesis study (Rooij and Kuhl, 2018). Based on our mathematical model, we simulated foliation caused by GC proliferation in the EGL, and then GC migration from the EGL to the IGL. By altering GC migration speeds, we showed that GC migration caused folia lengthening. Furthermore, a simulation of impaired GC migration under pathological conditions showed that fiber orientation was necessary for folia lengthening.

Continuum mechanics modeling for cerebellar tissue deformation
During cerebellar morphogenesis, foliation proceeds with GC dynamics including GC proliferation in the EGL and GC migration to the IGL along BG fibers (Fig. 1A). To investigate foliation mechanisms, we developed a mathematical model of cerebellar tissue deformation, GC proliferation, and migration based on continuum mechanics. The cerebellar tissue was modeled as a continuum body. We considered points in the body, whose position was denoted by at the initial state. In the current (deformed) state at time , the position of the same point was moved to by a deformation map χ ( = χ( , )). Tissue deformation was represented by the deformation gradient tensor defined as = ∂χ( , )/ ∂ . To model volumetric growth, we assumed that the deformation gradient tensor was multiplicatively decomposed into a growth part g and an elastic deformation part e (Rodriguez et al., 1994) as: (1) In this study, we assumed isotropic tissue growth because the anisotropy of GC accumulation in IGL was not observed. The deformation gradient g can be described by using growth stretch as: where was a second order unit tensor. To model tissue growth caused by GC accumulation, we assumed that the volumetric growth depended on cell number density. Thus, the growth stretch was described as a function of cell number density (≥ 0) in the current state as: where and are constants to represent the relationship between tissue growth and cell number density (Rooij and Kuhl, 2018).

Modeling GC proliferation and migration
To model GC proliferation and migration, the balance equation for cell number density in the current state was described as: where (> 0) is the Jacobian of the deformation gradient tensor , the vector is cell migration velocity, is a diffusion coefficient, and is the production rate of cell number density at the current state. The production rate at the initial state is represented using as = .
The second term on the left side of Eq. (4) shows the effect of tissue deformation on cell number density. The first term on the right side of Eq. (4) shows an advection term for GC migration, the second one shows a diffusion term, and the third one shows a production term for GC proliferation.
Under physiological conditions, where GC migration is guided along BG fibers, GC migration velocity is described by a directional unit vector denoting fiber orientation as: By assuming current fiber orientation depends on tissue deformation, thus is described by the deformation gradient tensor as: where 0 is the initial fiber orientation. To model impaired GC migration under pathological conditions, where GCs do not migrate along BG fibers, we randomly set an angle r (− max /2 ≤ r ≤ max /2) between the fiber direction and the migration direction . Through abovementioned modeling, we constructed the mathematical model for GC migration guided by the fiber orientation , which changes from the initial orientation 0 depending on tissue deformation , as formulated in Eqs. (5) and (6). The cell number density of GC temporally changes through GC migration and production according to Eq. (4) and affects tissue growth and subsequent tissue deformation, as per Eq. (1-3).

A simulation model for cerebellar foliation
To investigate foliation mechanisms, we performed numerical simulations based on the above models, using a finite element method. Cerebellar tissue, consisting of EGL, ML, PCL, and IGL, was modeled in a two-dimensional space by assuming plane strain condition (Fig. 1B). Based on experimental observations (Sudarov and Joyner 2007), the tissue had a length = 400 μm and a thickness = 50 μm, at the initial state ( = 0) and each layer had the same thickness. By presupposing the presence of folia at the embryonic day E16.5 (Sudarov and Joyner 2007), the initial tissue shape was slightly curved and followed a sine function having an amplitude of 5 μm. In numerical simulations, the outof-plane displacement on the lateral surfaces was constrained with fixed top points.
During cerebellar foliation, not all the GCs in the EGL proliferate. We used a cell number density c of the migrating GCs, which do not have a proliferative ability, and assumed that a fixed number of proliferative GCs existed in EGL. Therefore, to model GC production in the EGL, the production rate of cell number density (= ) was set to = 0 ( 0 > 0), and = 0 in the other layers. To model GC migration from the EGL to the IGL, the migration speed | | was set to | | = 0 ( 0 ≥ 0 ) in the EGL, ML and PCL, and | | = 0 in the IGL. The migration direction was determined depending on BG fiber orientation as described in Eq. (5). The directional unit vector of fiber orientation 0 at the initial state was set to 0 = 0 and 0 = −1, where 0 and 0 are xand y-component of 0 . These model parameters were defined as listed in Table 1.

Under physiological conditions GC migration causes non-uniform cell accumulation in the IGL, and folia lengthening
To understand the contribution of GC migration along the BG fibers (oriented based on tissue deformation) to cerebellar morphogenesis, we performed numerical analyses of cerebellar tissue morphogenesis, by altering the cell migration velocity 0 as 10 μm/h with reference to the experimental measurement (Yacubova and Komuro, 2002b) and as 5 and 0 μm/h for comparison. Our simulation data showed that GC proliferation in the EGL and their subsequent migration from the EGL to the IGL, caused folia formation (Fig. 2 A-C). The distribution of GC cell number density was temporally changed ( Fig. 2A-C; right section). When the migration speed 0 was 10 μm/h, GC increased in the IGL when compared to other layers, over time = 20 h ( Fig. 2A). Cell numbers density in the IGL had a non-uniform distribution, and the density were higher around the lobule than the fissure. As the migration speed decreased ( 0 = 5 μm/h), cell number density decreased in the IGL, but increased in other layers, including the EGL, ML, and PCL, over time = 20 h (Fig. 2B). When GCs did not migrate ( 0 = 0 μm/h), cell number density became uniformly distributed in the EGL, over time = 20 h (Fig. 2C). These data suggested that higher migration speeds caused a non-uniform GC distribution.
To quantitatively analyze the non-uniformity of IGL thickness, temporal changes in difference in IGL thickness around lobule l and fissure f were plotted (Fig. 2D). Differences in the IGL thickness were larger when GCs migrate ( 0 = 10, 5 μm/h) than when GCs do not migrate ( 0 = 0 μm/h); this result suggested that GC accumulation mediated by GC migration contributed to the formation of an IGL of non-uniform thickness. To analyze the effects of GC migration on folia lengthening, temporal changes in the lobule height ℎ was plotted (Fig. 2E). We observed that lobule height increased as migration speed increased. These data suggested that fiber-guided GC migration caused the formation of IGL of non-uniform thickness and folia lengthening, which are characteristic morphogenesis features of the cerebellar tissue.

Impaired GC migration causes defects in folia formation under pathological conditions
As proposed in a previous study, the orientation of BG fibers, depending on tissue deformation, may be key to the formation of IGL of non-uniform thickness and folia elongation (Sudarov and Joyner, 2007). To investigate the contribution of fiber-guided GC migration, we performed numerical simulations on impaired GC migration, under pathological conditions. Impaired GC migration was modeled, as the migration direction was disturbed by a random angle , ranging from − max /2 to max /2. We fixed the migration speed at 0 = 10 μm/h.
The direction field in the EGL, ML and PCL at time = 20 h was distributed as shown (Fig 3 A-C; grey arrows in the enlarged view). When the migration direction generated no disturbance ( max = 0) corresponding to a physiological condition as previously simulated (Section 3.1), the migration direction spread radially toward the IGL around the fissure, and the migration direction converged toward the IGL around the lobule (Fig. 3A). When the random angle range max was set to π and 3π/2 under pathological conditions, the migration direction was disoriented ( Fig. 3B and C; grey arrows in the enlarged view). Cell number density decreased in the IGL, but increased in other layers as random angle range increased (Fig. 3A-C; lower section). To investigate the effects of fiber-guided migration on folia morphology, we plotted temporal changes in difference in IGL thickness between the fissure and lobule (= l − f ), and lobule height ℎ ( Fig. 3D and E). These data indicated that the difference in IGL thickness and the lobule height were smaller as migration randomness increased. Thus, we suggested that an elongated lobule and non-uniform IGL thickness during cerebellar morphogenesis were caused by deformation-dependent GC migration.

Discussion
During cerebellar morphogenesis, the characteristic morphology of the cerebellar tissue is formed via GC migration along BG fibers (Corrales et al., 2006;Legue et al., 2016;Leto et al., 2016;Govek et al., 2018;Lejeune et al., 2019). To investigate the GC mediated foliation, we performed numerical simulations of the cerebellar tissue, driven by GC proliferation and migration. Our simulation results under physiological conditions showed that GC migration from the EGL to the IGL caused elongated folia and a non-uniform IGL which was thicker at the lobule, than the fissure (Fig. 2). Moreover, impaired GC migration under pathological conditions, modeled as random migration directions, resulted in reduced folia lengthening and formation of IGL with uniform thickness (Fig. 3). These results suggested that GC migration along BG fibers was implicated in the formation of characteristic cerebellar morphology, including nonuniform GC accumulation in the IGL, and lengthening folia.
This study validated foliation mechanisms via computer simulation, i.e., folia lengthening accompanied by nonuniform IGL appeared to be caused by GC proliferation and migration along BG fibers (Sudarov and Joyner, 2007). The direction of GC migration corresponding to BG fiber orientation (Fig. 3A) agreed with a previous study (Govek et al., 2018). With impaired and random GC migration, we observed that abnormal foliation, including reduced folia lengthening and formation of IGL with uniform thickness, was caused ( Fig. 3B and C). Importantly, these results agree with experimental data; impaired GC migration leads to widely-distributed GCs and a thin IGL (Govek et al., 2018). Therefore, our numerical simulation data successfully recapitulated the features of physiological/pathological foliation processes previously observed in experimental studies (Sudarov and Joyner, 2007;Govek et al., 2018).
We simulated and identified the physical mechanisms that underpin deformation-dependent cell migration and cause non-uniform layer formation during foliation. Such layer formation has been previously explained as the result of biochemically-regulated neuronal migration (D'Arcangelo and Curran, 1998;Sekine et al., 2012;Caffrey et al., 2014;Matsunaga et al., 2017). In addition, the identification of mechanical perspectives is also important in understanding layer formation mechanisms (Lejeune et al., 2016;Engstrom et al., 2018;Lawton et al., 2019). It was previously shown that tissue elastic deformations, such as stretch and compression, during wrinkling caused differences in thickness of multi-layered tissue between lobule and fissure (Tallinen et al., 2014;Holland et al., 2018). In addition to this mechanical mechanism, we also proposed a mechanism based on plastic deformation caused by changes in cell number via fiberguided migration. In physiological cerebellar tissue, GC number around lobule is clearly more than around fissure (Sudarov and Joyner, 2007). Because fiber-guided migration is commonly observed during brain morphogenesis such as neuronal migration along radial glial fibers in cerebral development (Borrell and Gotz, 2014;Del Toro et al., 2017;Rahimi-Balaei et al., 2018), our physical mechanism will be important for our increased understanding of layer formation during this process. In our simulations, we assumed the thickness of each layer at the initial state to be equal (Fig. 1B). The thickness of each layer is considered to affect the number of produced GCs and extent of GC migration. In addition, we assumed the initial folia height as a specific value, which can influence the magnitude of deformation. However, because these assumptions affected only the quantitative results, the setting of the initial tissue state does not affect the qualitative results for the proposed foliation mechanism, where the folia elongate via non-uniform GC accumulation caused by fiberguided GC migration. For the simulation of a specific condition, a further analysis is required. To apply this model to a three-dimensional analysis for the entire cerebellum, the effects of a cerebellar tissue curvature on foliation must be investigated in the future.
The foliation mechanism investigated in this research presupposed the existence of the initial fissure, even though initial foliation was important for the determination of cerebellar morphology. To comprehensively understand the foliation, it is important to investigate the formation of an initial fissure by considering other mechanisms, such as fluidlike tissue behavior resulting from cell migration in the EGL (Engstrom et al., 2018;Lawton et al., 2019). Moreover, folia bifurcation in parallel with folia lengthening is important for the study of cerebellar morphogenesis. To investigate folia bifurcation, other mechanical factors that cause bifurcation must be considered, such as lobule-lobule contact. In addition to mechanisms implicated in the formation of characteristic tissue morphology, further developments may reveal mechanisms for establishment of the physiological function by considering the tissue morphology associates with its function (Sillitoe and Joyner, 2007;White and Sillitoe, 2013;Leto et al., 2016).

Conclusion
In this study, we constructed the computational model for cerebellar foliation caused by fiber-guided cell migration. Through finite element analysis, we successfully recapitulated the features of foliation processes and validated foliation mechanisms, i.e., folia lengthening accompanied by non-uniform IGL appeared to be caused by GC proliferation and migration along BG fibers. Our computational approach has provided insights on the biological and physical mechanisms that characterize cerebellar morphogenesis, and may facilitate a greater understanding of interactions between cellular activities and tissue behaviors.