Models based in differential expansion of elastic material, axonal constraints, directed growth, or multi-phasic combinations have been proposed to explain brain folding. However, the cellular and physical processes present during folding have not been defined. We used the murine cerebellum to challenge folding models with in vivo data. We show that at folding initiation differential expansion is created by the outer layer of proliferating progenitors expanding faster than the core. However, the stiffness differential, compressive forces, and emergent thickness variations required by elastic material models are not present. We find that folding occurs without an obvious cellular pre-pattern, that the outer layer expansion is uniform and fluid-like, and that the cerebellum is under radial and circumferential constraints. Lastly, we find that a multi-phase model incorporating differential expansion of a fluid outer layer and radial and circumferential constraints approximates the in vivo shape evolution observed during initiation of cerebellar folding.https://doi.org/10.7554/eLife.45019.001
The human brain has a characteristic pattern of ridges and grooves that make up its folded shape. Folds in the outer layer of the brain, known as the cortex, increase the surface area and make more space for cells to connect and form complex circuitries. Different models have been put forward to explain how these folds form during development. Examples include tension from cells pulling areas of the cortex together, or layers of the cortex growing at different rates, causing the cortex to buckle and create folds. Discriminating between these different models requires biological information about the cells and tissue of the brain at the start of the folding process. However, it has been difficult to extract this information when considering the development of the human brain in three dimensions.
Lawton et al. have overcome these difficulties by using a part of the mouse brain called the cerebellum as a simpler system. As in humans, the mouse cerebellum is a densely folded structure, sitting underneath the brain, that plays a major role in regulating movements, as well as cognition. The symmetrical structure of the mouse cerebellum means it can be analyzed in two dimensions, making it easier to track the mechanics of folding.
By applying the extracted biological data onto a mathematical model, Lawton et al. showed folding was driven by a combination of previously unknown features. For instance, that cells in the outer layer of the cerebellum grow faster than cells in the center, with cells growing uniformly across the outer layer. Other features include the fluid-like composition of the outer layer, which allows cells to move freely and regularly change position, and tensions surrounding the cerebellum mechanically straining its growth. Notably, the pattern of cells and tissue fibers in the cortex had no influence over these mechanical properties and provided no pre-indication of where the sites of folding would occur. The data collected deviates from other models, and has led Lawton et al. to propose a new explanation for how the brain folds, incorporating these newly found features.
Problems with brain folding during human development can lead to debilitating conditions. Applying this new model to folding disorders of the human brain could help scientists to understand how these folding defects arise.https://doi.org/10.7554/eLife.45019.002
Recent work to elucidate the mechanics of brain folding has primarily focused on the human cerebral cortex and involved models of directed growth, axonal tension, or differential expansion of elastic materials that generate compressive forces to drive mechanical instabilities leading to folding (Tallinen et al., 2014; Ronan et al., 2014; Bayly et al., 2013; Xu et al., 2010; Hohlfeld and Mahadevan, 2011; Bayly et al., 2014; Lejeune et al., 2016; Karzbrun et al., 2018). Current elastic material models are able to create three-dimensional shapes strikingly similar to the final folds seen in the adult human cortex (Tallinen et al., 2016). A recent multi-phase model (Engstrom et al., 2018) that includes elastic and fluid-like layers, differential expansion and radial constraints takes into consideration that multiple factors could lead to folding in the developing brain. However, the cell and tissue level mechanics actually present at the initiation of folding have not been considered or defined, as technological limitations are significant in animals with a folded cerebrum.
The murine cerebellum has a simple alignment of 8–10 stereotypical folds along the anterior-posterior axis. Combined with the genetic tools available in mouse this allows for precise developmental interrogation to identify and analyze the in vivo cellular and tissue level behaviors driving growth and folding. The developing cerebellum is distinct from the cerebral cortex, as it has a temporary external granule cell layer (EGL) of proliferating granule cell precursors that cover the surface and generate growth primarily in the anterior-posterior (AP) direction (Leto et al., 2016; Legué et al., 2015; Legué et al., 2016). During development a thickening occurs in the EGL at the base of each forming fissure, termed anchoring center (AC) (Sudarov and Joyner, 2007), whereas in the adult cerebellum the inner granule cell layer (IGL), generated by the EGL during the first two weeks of life, is thinnest at the ACs. Previous work on cerebellar folding utilized a tri-layer elastic model incorporating the EGL, the adjacent molecular layer, rich in axons and dendrites, and the IGL (Lejeune et al., 2016). However, neither the molecular layer nor the IGL are present when folding is initiated in the embryo. Therefore we argue that a bilayer system consisting of the EGL and underlying core, is a more appropriate approximation for cerebellar folding.
Here we show that cerebellar folding emerges from differential expansion between an un-patterned, rapidly expanding EGL and an underlying core. Additionally, we demonstrate that the measured stiffness differential, compressive forces, and the thickness variation in the EGL are all inconsistent with traditional elastic wrinkling models driven by differential growth. Furthermore, we demonstrate that the expansion of the EGL is uniform, and fluid-like, and that the cerebellum is under radial and circumferential constraints when folding initiates. Lastly, we constrain the recent multi-phase model with our in vivo data and find we can capture the temporal shape evolution seen during mouse cerebellum folding initiation. The implications of our findings for human cerebral cortex folding are discussed.
It is well known that differentially expanding bilayer systems can wrinkle to relax building stress (Richman et al., 1975; Nelson, 2016; Hannezo et al., 2012; Shyer et al., 2013; Wiggs et al., 1997). We reasoned that in the cerebellum the EGL could behave as a quickly expanding outer layer and its attachment to a more slowly growing core could generate forces that result in a wrinkling-like phenotype. To test whether the cerebellum has differential expansion between the two layers, we measured the expansion of the EGL and the core during the time of initiation of folding from midline sagittal sections (Figure 1a–d). Unlike the cerebral cortex, the unfolded murine cerebellum is a simple cylinder-like structure elongated in the medio-lateral axis (Figure 1e) (Szulc et al., 2015). All folds in the medial cerebellum (vermis) are aligned in the same axis allowing 2-D measurements to estimate expansion in the anterior-posterior axis of the vermis. Therefore the length of the surface of the EGL was used as a measure of the cerebellum surface area and the area of the core as an approximation of cerebellum volume (Figure 1d), and measurements were made each day from embryonic day 16.5 (E16.5) through postnatal day 0 (P0). In cross-section the unfolded cerebellum approximates a semicircle, therefore we reasoned that if the cerebellum were to remain unfolded then the ratio of expansion between the length of the EGL and the area of the core should approximate the ratio of the circumference of a semi-circle to its area. Of significance, we found that at E16.5 and E17.5 the ratios of growth between the EGL and core closely approximated the expansion of a semi-circle. However, at E18.5 and P0 the expansion rate of the EGL was greater than the rate of core expansion (Figure 1f). Thus we uncovered that the cerebellum does indeed go through a phase of differential expansion. We next determined whether differential expansion correlates with when folding occurs by calculating a folding index (the convex curvature of the EGL divided by the length of the EGL) at each stage (Mota and Herculano-Houzel, 2015). Indeed, we found that the cerebellum remains unfolded during the initial proportional expansion between the EGL and core and only folds when the differential expansion is initiated (Figure 1g). These results provide quantitative evidence that cerebellar folding involves tissue level mechanical forces arising from differential expansion.
Since there is differential expansion between the EGL and the core and as this type of expansion is the driver of elastic bilayer models we tested whether the properties of cerebellar tissue are consistent with the requirements and predictions of such models. Briefly, the initial resulting wrinkling instability defines the distance between folds as the initial sinusoidal undulations increase in amplitude to ultimately turn into lobules. The folding wavelength depends on the thickness of the external layer (EGL) and the ratio of the stiffness of the two layers (EGL/core). In particular, for a planar geometry, with the stiffness of the external layer defined as Eo, the stiffness of the core as Ei, and the thickness of the external layer denoted as t, the folding wavelength λ is given by Allen (1969)
If the length of the system is , then the number of folds is inversely proportional to the thickness of the EGL
We explored a standard elastic bilayer model in a circular geometry using the observed ratio of thickness of the EGL to radius of the cerebellum near the onset of shape change (E16.5) and invoked a neo-Hookean elastic solid for both layers (Zhao and Zhao, 2017). The resulting shape change was studied as a function of the ratio of the layer stiffness values (Figure 2a). We found that to produce the observed number of folds (three in the semi-circular cerebellum and six in the circular model) at initiation of folding through wrinkling based models constrained by our measurements of the embryonic cerebellum, a large stiffness ratio was required of around 50. To map the stiffness contrast in the cerebellum we used scanning acoustic microscopy (SAM) to measure the bulk modulus of the cerebellum daily from E16.5 to P18.5 (Figure 2b–c, Figure 2—figure supplement 1) using established methods (Rohrbach et al., 2015; Rohrbach et al., 2018; Rohrbach and Mamou, 2018). For small deformations, the instantaneous bulk modulus should linearly relate to the stiffness and, therefore, the ratio of the instantaneous bulk moduli should scale similarly to the ratio of stiffnesses (assuming the same Poisson’s ratio for the EGL and for the core, neither of which have been directly measured). While this qualitative approach and SAM tissue preparation protocols may not be able to produce the absolute values of the elastic properties of the tissues, it can give a reasonable indication of the relative stiffnesses of different parts in the cerebellum. Using this estimation, we found that the EGL has a slightly higher instantaneous bulk modulus than the core at all stages measured. Unsurprisingly, the ratio (~1.05:1) was not close to being sufficient to produce a folding wavelength similar to that in the cerebellum (Figure 2d). Consistent with our finding, small modulus contrasts have been reported for other brain regions with multiple loading modes, such as shear, compression, and tension (Xu et al., 2010; Lejeune et al., 2016; Budday et al., 2017). Elastic material models with graded growth profiles have been developed that predict folding of cerebral cortex without a large stiffness differential (Tallinen et al., 2014). However, these models are still bound by other measurable requirements as discussed below.
Elastic bi-layer wrinkling models predict compressive forces in the outer layer. Simulations performed of cuts through the outer layer and into the inner layer predict that upon relaxation the outer layer should not open (Figure 2e). We tested whether this prediction reflects the biology using surgical dissection blades to make radial cuts across the meninges, EGL, and into the core of live E16.5 tissue slices. Time-lapse imaging revealed that, in contrast to the prediction, the EGL opens as well as part of the underlying cut in the core (Figure 2f–h, Figure 2—figure supplement 2a–c, and Video 1). This result indicates there is circumferential tension within the outer layers of the cerebellum. This finding also rules out the elastic models with graded growth profiles as they predict compressive forces in the outer region as well.
The elastic bi-layer model requires the EGL to be thinnest at the base of each AC, which are the lowest parts of the cerebellar surface. Thus, the EGL should have an ‘in-phase’ thickness variation. Without this feature, a purely elastic model – bi-layer based or even graded growth profile based – cannot be in mechanical equilibrium (in the quasistatic limit) (Engstrom et al., 2018). However, we previously reported that the embryonic EGL is thickest in the ACs when folding initiates, that is it has an ‘out-of-phase’ thickness variation (Sudarov and Joyner, 2007). To validate this observation, we quantified the thickness variations in the EGL centered at the ACs present at E16.5–18.5. Not all cerebella have visible AC at E16.5. However in the subset that do and in the three ACs present at E17.5, the EGL was found to be 1.2–1.4 times thicker in the ACs than in the surrounding EGL (Figure 2i–l and Figure 2—figure supplement 3). Moreover, the thickness ratio increased to 1.7 times at E18.5 (Figure 2l). As described above, the final thickness variations of the IGL (as well as the molecular layer) of the cerebellar cortex are in-phase, just as the layers of the adult cerebral cortex. These results further show that traditional elastic wrinkling models cannot capture the initiation of cerebellum folding, and highlight the importance of making biological measurements at the time of folding rather than when it is complete.
As elastic bi-layer models do not align with the biology of cerebellar folding, we looked for other drivers of morphometric changes. Since the EGL drives the majority of cerebellar growth (Leto et al., 2016; Legué et al., 2015; Legué et al., 2016), we first tested whether regional differences in EGL proliferation rates are present that could influence the folding pattern of the cerebellum. Proliferation rates (S phase index) were measured in the EGL during folding initiation (E16.5 and E17.5) in the inbred FVB/N strain to reduce variation between samples. First we asked if the regions that will give rise to distinct sets of lobules have different rates of proliferation that could contribute to the larger and smaller sizes that the lobules ultimately attain. We focused on the anterior cerebellum that divides into a larger region with lobules 1–3 (L123) and smaller region (L45), as well as the central area that comprises lobules 6–8 (L678) of the cerebellum (Figure 3a–b). The more posterior cerebellum does not consistently fold at this stage, thus measurements were not included. Interestingly, we found that the proliferation rates were similar in the three regions at E16.5 (Figure 3c). The EGL proliferation rate at E17.5 in L678 was slightly reduced compared to the L123 region, but no other differences were found (Figure 3d). Thus proliferation is uniform just before initiation of folding and the small difference found during folding does not correlate with lobule size. This result indicates that lobule size is not determined by modulating the levels of proliferation at the onset of folding. Rather, lobule size could be set by both the timing of invagination, and the distance between ACs as granule cell precursors in one lobule do not cross the surrounding ACs to contribute to an adjacent lobule (Legué et al., 2015).
Each AC is first detected as a regional inward thickening of the EGL (Sudarov and Joyner, 2007) (Figure 2i–l and Figure 2—figure supplement 3). We measured the proliferation of the EGL specifically within the forming AC regions to test whether altered proliferation rates could explain the thickenings and therefore the initiation of an AC. We found the rate of proliferation within each forming AC region at E16.5 and E17.5 was the same as in the surrounding EGL (Figure 3e,f), thus proliferation within all regions of the EGL at the initiation of folding is uniform. Furthermore, regional modulation of proliferation does not form or position the ACs.
At E18.5, after the initiation of folding, we found that the rate of proliferation was significantly lower in the L678 region compared with the L123 and L45 regions (Figure 3—figure supplement 1a). However, proliferation within the ACs at E18.5 remained uniform with the surrounding regions (Figure 3—figure supplement 1b). Since ACs compartmentalize the EGL, our results show that regional differences in proliferation rates arise in lobule regions after initiation of folding, which thus could be important for determining the ultimate size of the folds.
Changes in cell size and shape have been shown to induce morphological changes (Mammoto and Ingber, 2010; Harding et al., 2014; Stemple, 2005; He et al., 2014). To test if regionally specific regulation of cell shape or size directs folding, we fluorescently labeled cell membranes of scattered granule cell precursors (GCPs) in the EGL using genetics (Atoh1-CreER/+; R26MTMG/+ mice injected with tamoxifen two days prior to analysis). We then segmented the cells in 3D and quantified their sphericity (Figure 3g). We discovered that GCPs in the EGL take on a large variation of shapes and sizes at E16.5 and E18.5. However, we found no difference in cell shape in the different lobule regions of the EGL or between the AC areas and the surrounding EGL at each age (Figure 3h,i). Cell size was uniform at both stages except for a slight reduction in L678 at E16.5 when compared with L123 and the AC regions. However, the size of cells is reduced at E18.5 compared to E16.5 (Figure 3j,k and Figure 3—figure supplement 1c). Thus, the proliferating GCPs that drive expansion of the EGL have both uniform proliferation rates and similar shapes and sizes across the lobule regions defined by the first three ACs at folding initiation.
The EGL is traversed by fibers of Bergmann glial and radial glial cells (Leung and Li, 2018; Yuasa, 1996; Yamada and Watanabe, 2002). We tested whether the fibers are distributed in patterns that could locally change the physical properties of the EGL and induce invaginations. Genetics was used to fluorescently label cell membranes of scattered glial cells (nestin-creER/+;R26MTMG/+ mice injected with tamoxifen at E14.5) (Figure 4a). Fibers crossing the EGL at E16.5 were counted in sagittal slices and aligned relative to the ACs (Figure 4b). This analysis showed that the Bergmann glial and radial glial fibers are distributed evenly along the AP axis of the EGL, and therefore are not directing the positions where folding initiates based on an uneven regional distribution.
Tension based folding models suggest constraints from axons and other fibers could direct folding (Xu et al., 2010; Van Essen, 1997). Since the cerebellum is under circumferential tension, as demonstrated above, we examined evidence of radial tension between the EGL and the ventricular zone (VZ) at the initiation of folding. Cuts were made in live E16.5 tissue slices between the EGL and VZ running approximately parallel to them so that they cut across radial fibers in the anterior cerebellum (Figure 4c). As predicted, after cutting the tissue relaxed revealing tension directed radially within the cerebellum (Figure 4d,e and Figure 2—figure supplement 2 and Video 2). Interestingly, quantification of how the radial and horizontal cuts open revealed that only the horizontal cuts opened along the full length of the cut although they opened more slowly than radial cuts (Figure 2—figure supplement 2g–j), indicating different stress profiles in the two orientations.
Taken together, at the time of folding initiation the EGL, which is driving the differential expansion, is itself growing uniformly and the cerebellum is under both radial and circumferential constraints. Finally, there is no evidence of any pre-patterning in the EGL in either cellular behaviors or fiber distribution.
As the granule cells within the EGL have such varied shapes as shown above, we looked to see if the cells within the EGL were undergoing any rearrangement movements that may indicate fluid properties. A small, scattered fraction of nuclei in the EGL were fluorescently labeled (Atoh1-CreER/+; R26ntdTom/+ injected with tamoxifen two days prior to imaging) and ex vivo slice-culture time-lapse imaging was performed for up to five hours. Tracking the cell positions through time revealed that granule cells within the EGL are highly motile within the EGL. Furthermore, there was no obvious directionality or collectivity to the movement. However, the dynamic motility resulted in the constant exchanging of nearest neighbors over the course of tens to hundreds of minutes and shows that at the timescale of folding the EGL is more fluid-like than a solid epithelial layer (Figure 5 and Videos 3 and 4).
We recently developed a model for folding from differentially expanding bi-layer tissues that takes into account the out-of-phase thickness of the outer layer of several systems and possible contribution of radial mechanical constraints present in neurological tissue (Engstrom et al., 2018). We applied the model here to the initiation of cerebellar folding based on five primary assumptions. First, the core is an incompressible material (μ) as indicated by the bulk modulus measurements. Second the outer layer, that is the EGL, expands uniformly (kt) as shown by the proliferation rate. Third, the EGL is assumed to be a fluid-like material as demonstrated by the live-imaging of neighbor exchanges. Fourth, there is an elastic component radially to the entire cerebellum (kr), seen in the cutting and relaxation experiment and possibly mediated by radial glia. Fifth, the EGL is constrained towards a uniform thickness (β), possibly by Bergmann glia fibers spanning the EGL. Given the interplay between incompressible material, compressible fibrous material, and a proliferating non-elastic EGL, this model is multi-phase.
An energy functional parameterized by both the inner and outer boundary of the EGL and incorporating the above five assumptions into three dimensionless parameters (μ/kr, kr/kt, kt/β) is minimized to yield an equation for a driven harmonic oscillator resulting in sinusoidal shapes for both the inner and outer boundary of the EGL given an initial elliptical shape. In contrast with the elastic bilayer wrinkling model, EGL thickness oscillations are found to be out-of-phase with the surface height (radius) oscillations when 0 < μ/kr <1. Additionally, the model predicts that the ratio of the measured surface height amplitude (Ar) and the EGL thickness amplitude (At) is given by
which need not be as is typical of elastic bilayer wrinkling, and the number of initial folds at E16.5 is determined by
Note that in contrast with elastic wrinkling, the number of initial folds does not depend on the thickness (a length scale) of the EGL, but only on material properties.
Given that our tissue cutting and relaxation experiment revealed circumferential tension in the cerebellum at folding initiation (Figure 2f–h, Video 1), we returned to the mathematics and found a previously unrealized geometric relationship in the circular limit of the model that in fact assumes circumferential tension in addition to the previously discussed radial tension given that the perimeter of a circle is determined by its radius.
To rigorously test the shape prediction of the model, we first constrained 3 of the 5 parameters for a circular version of the model by using both the thickness amplitude, and average thickness of the EGL, as measured at E16.5, and the number of initial folds. Secondly, the parameter μ/kr (denoted as ε) was assumed to scale linearly with time. Together, this allowed for the generation of shape predictions at later developmental stages (E17.5 and E18.5) from the E16.5 starting approximation. Solving the analytical model as constrained by our measured embryonic data we found that it closely approximates the phase and amplitude behavior of EGL thickness and radius oscillations from E16.5 through E18.5. (Figure 6a–c). However, the model is not able to produce self-contacting folds or hierarchical folding, both of which are seen in the cerebellum at later stages.
The cerebellum has hierarchical folding in which the initial folds become subdivided. Given that ACs hold their position during development and compartmentalize granule cells within lobules of the EGL (Legué et al., 2015) we reasoned that the ACs could be acting as mechanical boundaries enabling similar mechanics to drive the secondary folding. To test this possibility we measured the expansion of the EGL and the core of the individual lobule regions from E18.5 to P3. We found that indeed in the lobule regions that undergo folding there is a temporal correlation between when the onset of sub-folding and differential expansion occur (Figure 7a–d). In contrast, the region (L45) that does not fold during the same time period has a different, more rectangular shape, and the ratio of EGL growth to core growth is proportional for a rectangle during the time measured (Figure 7).
Here we have provided experimental evidence that cerebellar folding emerges without obvious pre-patterning. Additionally, the outer layer has fluid-like properties and expands uniformly, and the growth creates a differential expansion between the outer layer and the core. Thus, traditional morphometric cellular behaviors such as changes in cell shape, size and proliferation do not direct where cerebellar folding initiates. Furthermore, our developmental interrogation revealed tissue moduli, mechanical constraints, and emergent thickness variations in the EGL that are fundamentally inconsistent with traditional elastic bilayer wrinkling models. Therefore our results call for a new understanding of brain folding.
By applying a multi-phase model constrained by our measured data we were able to capture the correct shape variations and number of folds at the onset of folding. Our new framework accounts for: the rapidly expanding fluid EGL, whose thickness is proposed to be regulated by Bergmann glial fibers, the slower growing incompressible core, and fibrous material in the form of glial fibers and possibly axons as well as the meninges that potentially provide radial and circumferential tension (Figure 8). This multi-phase model of folding makes many new predictions. One such prediction is that adjusting the amount of tension spanning the cerebellum will change the degree of folding. Indeed, alterations of the cells that likely create tension-based forces could explain the dramatically disrupted folding seen in mouse mutants in which radial glia do not produce Bergmann glia (Li et al., 2014). Without Bergmann glia, the EGL would be expected to not form a layer with regular thickness and it should be more sensitive to variations in radial glial tension. Consistent with this prediction, mutants without Bergmann glia have more localized and less regular folds (Li et al., 2014) Our combination of experimental studies and modeling thus provide new insights into cerebellar folding, including an underappreciated role for tension.
Under the new framework revealed by our measurements made in the developing mouse cerebellum, to approximate the observed shape changes in the murine cerebellum from E17.5 to E18.5 the ratio of the core stiffness over the radial tension must increase. Yet, the measured bulk modulus of the core shows no increase during development. Therefore a second prediction is that radial tension must decrease during development. While the cerebellum is crossed by many fibers at folding initiation, radial glial fibers are an attractive candidate to mediate this change in radial tension (Sillitoe and Joyner, 2007; Rahimi-Balaei et al., 2015). First, they span from the VZ to the surface of the cerebellum at E16.5. Additionally, during folding initiation the radial glia undergo a transition into Bergmann glia where they release their basal connection to the VZ and the cell body migrates towards the surface (Mota and Herculano-Houzel, 2015). This transition could lead to a reduction in the global radial tension and thus would be consistent with our model prediction.
The mechanics underlying hierarchical folding remain an open challenge. However, our developmental data may provide a way forward. As ACs maintain their spatial positions, and as they compartmentalize granule cells within the EGL into the lobule regions (Legué et al., 2015), we propose that they create fixed mechanical boundaries that divide the cerebellum into self-similar domains. These domains, with their similar physical properties to the initial unfolded cerebellum, can then undergo additional folding. Furthermore since ACs compartmentalize granule cells within the lobule regions, once separated the lobule regions can develop distinct characteristics, like the observed differential proliferation rates at E18.5. We speculate, therefore, that the folding patterns seen across cerebella in different species evolved by adjustment of global as well as regional levels of differential expansion and tension which ultimately mold the functionality of the cerebellum.
Finally it is interesting to note the similarities and differences between the developing cerebellum and the cerebral cortex. Radial glia span the entire cerebral cortex just as in the cerebellum (Götz et al., 2002). Furthermore, species with folded cerebrums have evolved outer radial glial cells for which the cell body leaves the ventricular zone to become positioned near the surface while retaining fibers anchored on the surface, similar to Bergmann glia in the cerebellum (Leung and Li, 2018; Reillo et al., 2011). While we have emphasized the notion of tension via glial fibers in the developing cerebellum, axonal tension has been discussed in the context of shaping the developing cerebrum (Van Essen, 1997). Tissue cutting in the cerebral cortex of ferrets has revealed a similar tension pattern during folding as we found in the cerebellum (Xu et al., 2010). We therefore submit that our work calls for a revival of the notion of how tension affects the shape of the developing cerebrum.
Unlike the cerebellum, the cerebral cortex is not divided into a simple bilayer system. However, outer radial glial cells proliferate, much like the GCPs of the EGL, to drive the expansion of the outer regions of the cerebral cortex around the time of initiation of folding (Hansen et al., 2010; Heng et al., 2017; Nowakowski et al., 2016). Moving the zone of proliferation out from the VZ gives more space for the increased proliferation required in folding systems. The cerebellum, housing 80% of the neurons in the human brain may be an extreme example requiring the region of proliferation to be completely on the outer surface (Andersen et al., 1992). Constraining models of folding of different brain regions with developmental data will bring about a more accurate quantitative understanding of the shaping of the developing brain.
All experiments were performed following protocols approved by Memorial Sloan Kettering Cancer Center’s Institutional Animal Care and Use Committee. The inbred FVB/N stain was used for all proliferation rate, area, length, and expansion rate measurements. Atoh1-CreER (Machold and Fishell, 2005), Nestin-CreER (Imayoshi et al., 2006), Rosa26MTMG (Muzumdar et al., 2007), Rosa26Ai75 (Daigle et al., 2018) were used to quantify cell shape and size as well as fiber distribution and were maintained on the outbred Swiss Webster background. The Swiss Webster strain was used for scanning acoustic microscopy. Both sexes were used for the analysis. Animals were kept on a 12 hr light/dark cycle and food and water were supplied ad libitum. All experiments were performed following protocols approved by Memorial Sloan Kettering Cancer Center’s Institutional Animal Care and Use Committee.
The appearance of a vaginal plug set noon of the day as Embryonic day 0.5 (E0.5). All animals were collected within two hours of noon on the day of collection. Tamoxifen (Tm, Sigma-Aldrich) was dissolved in corn oil (Sigma-Aldrich) at 20 mg/mL. Pregnant females carrying litters with Atoh1-CreER/+;R26MTMG/MTMG or NestinCER/+;R26MTMG/MTMG embryos were given one 20 μg/g dose of TM via subcutaneous injection two days prior to analysis. 25 μg/g of 5-ethynyl-2-deoxyruidine (EDU; Invitrogen) was administered via subcutaneous injection one hour prior to collection.
For embryonic stages heads were fixed in 4% paraformaldehyde overnight at 4°C. For postnatal animals, the brain was dissected out first before fixation. Tissues were stored in 30% sucrose. For all proliferation, area, length, and thickness measurements brains were embedded in optimal cutting temperature (OCT) compound. Parasagittal sections were collect with a Leica cryostat (CM3050s) at 10 μm.
Prior to IHC, EdU was detected using a commercial kit (Invitrogen, C10340). Following EdU reaction the following primary antibodies were used either overnight at 4°C or 4 hr at room temperature: mouse anti-P27 (BD Pharmingen, 610241), rabbit anti-GFP (Life Technologies, A11122), rat anti-GFP (Nacalai Tesque, 04404–84). All antibodies were diluted to 1:500 in 2% milk (American Bioanalytical) and 0.2% Triton X-100 (Fisher Scientific). Alexa Fluor secondary antibodies (1:500; Invitrogen) were used: Alexa Fluor 488 donkey anti-rabbit, A21206, Alexa Fluor 488 donkey anti-rat, A21208, Alexa Fluor 488 donkey anti-mouse, A21202, Alexa Fluor 647 donkey anti-mouse, A31571. EdU was detected using a commercial kit (Invitrogen, C10340).
For cell size, shape and fiber density analysis 60 μm parasagittal sections were collected on a Leica vibratome (VT100S). Primary and secondary antibodies were diluted 1:500 in 2% milk and incubated overnight at 4°C.
For scanning acoustic microscopy brains were processed for paraffin embedding and parasagittal sections of 10 μm thick were collected on a microtome (Leica RM2255). Structured illumination and confocal Imaging was done with Zeiss Observer Z.1 with Apatome or Zeiss LSM 880 respectively.
Measurements for all analysis were taken from the three most midline sagittal sections and averaged. The most midline section was determined by dividing the distance in half between the lateral edges where the third ventrical and the mesencephalic vesicle are no longer connected. Quantifications were made using Imaris (Bitplane) and MATLAB (Mathworks) software.
EGL Proliferation rate was calculated as EDU+/(Dapi+;P27-) cells. All cells were counted within the lobule region to the midpoint of the Anchoring Centers. For proliferation measurements through the ACs and the surrounding EGL at E16.5 and E17.5 a 50 μm window measured from the outer surface of the EGL was centered at the AC. The measuring window was centered at every 25 μm anterior and posterior to the EGL for a total distance of 250 μm anterior and posterior to the AC. At E18.5, when the AC is fully formed, everything proximal to the centroid of the cerebellum under the midpoint of the AC was counted as the AC. Non-overlapping regions of 50 μm also were measured in either direction for a total of 200 μM anterior and posterior to the AC. Proliferation was measured in three cerebella at E16.5 and E17.5 and in four cerebella at E18.5.
EGL length was measured from the outer surface of the EGL following the curvature of the EGL. Cerebellar area was calculated as the area within the outer surface of the EGL and the ventricular zone. A short strait edge was made perpendicular to the ventricular zone to close the area back upon to the anterior end of the EGL. The convex curvature of the cerebellum was measured by following only the positive curvature of the EGL. The folding index was determined as FI = 1 - (Positive curvature/EGL length). Data collected for E16.5, E17.5, E18.5 and P0 came from 6,8,7, and 9 cerebella respectively.
EGL thickness was measured by defining the outer and inner curvature of the EGL. The shortest distance lines were drawn to the outer curvature from discrete points distributed at every 12.5 μm along the inner curvature of the EGL. Nine ACs and surrounding regions from five cerebella were quantified at E16.5 and 13 ACs from five cerebellar were analyzed at E17.5. At E18.5 six ACs from two cerebellar were quantified.
Midline sections were imaged with a Zeiss LSM 880. Serial images were taken to cover the entire EGL of lobule regions L123, L45, and L678 and the ACs. Manual cell masks were created with Imaris software defining the curvature at every z-slice. Every cell that was completely included in the imaging window and that was distinguishable from surrounding cells was counted to reduce sampling bias. Cells from three brains were measured at each stage for a total of 131 at E16.5 and 201 at E18.5. Shape was defined via sphericity, which is the surface area of a sphere having the same volume as the cell of interest divided by the surface area of the cell of interest.
Midline sections were imaged with a Zeiss LSM 880. Image tiling was used to cover the EGL. Using Python, a fourth or fifth order polynomial was fitted to the outer edge of the EGL in each image, and five scan lines were positioned at 12.2 μm intervals beneath the surface, and parallel to it. A bin width of 50 μm as measured along the polynomial contour was centered at the AC. Bins of equal distance were extended both anteriorly and posteriorly. Staining intensity was counted along each scan line at every z-slice of the confocal stack. Each image was normalized to the mean intensity and smoothed with a Gaussian filter. Peak counting was done using minimum and maximum filters, keeping neighborhood size and threshold parameters constant for all images. The results from the five scan lines were averaged.
Live cerebella of E16.5 FVB/N mice were collected in dissection buffer as previously described (Wojcinski et al., 2017) and embedded in low-melting point agarose (Invitrogen). Sagittal slices at a thickness of 250 μm were collected. Slices were removed from the agarose and place in petri-dishes coated with Poly(2-hydroxyethyl methacrylate)(Sigma-Aldrich). Tissue cuts (eight horizontal, 10 radial) were made with a 30° Premier Edge stab knife (Oasis Medical). Slices were allowed to relax for 10 min. Time-lapse images were acquired on a Leica MZ75 dissection scope.
Live cerebella of E16.5 Atoh1-CreER/+; R26Ai75/+mice were collected and slices of a thickness of 250 μm were cultured on Millicell cell culture inserts (Millipore) in glass bottom plates (Matek) as previously described (Wojcinski et al., 2017). Image stacks were acquired on a Zeiss LSM 880 at intervals of around 3.5 min for up to 5 hr. Cell positions were tracked using Imaris (Bitplane) software. Three time-lapses were analzyed.
Mechanical tissue properties were analyzed using a 250 MHz Scanning Acoustic Microscope (SAM), described previously (Rohrbach et al., 2015; Rohrbach et al., 2018; Rohrbach and Mamou, 2018). Briefly, 12 μm paraffin sections of mouse embryonic brains were de-parafinized, hydrated in de-ionized water and raster scanned (2 µm steps in both direction) on the SAM to acquire radio-frequency (RF) ultrasound data. At each scan location, signal processing was performed to compute the amplitude, sample thickness, speed of sound, acoustic impedance, attenuation, bulk modulus, and mass density (Rohrbach and Mamou, 2018). Two-dimensional maps of tissues properties were formed using the values obtained at each scan location. Bulk modulus was computed from the product of the acoustic impedance and the speed of sound. Co-registered histology and SAM amplitude images were used to identify regions-of-interest (ROIs) corresponding to the EGL layer and underlying core of the cerebellum in each sample. Bulk modulus was analyzed as a measure of tissue stiffness: ROI measurements were acquired from 3 sections from three embryos at each developmental stage.
The wrinkle of a circular bilayer structure in Figure 3a was simulated with commercial software ABAQUS. Both film and substrate were modeled as incompressible neo-Hookean materials. The ratio between shear moduli of the film and substrate was 50 and the initial radius of the simulated structure was 16 times that of the film thickness. The differential growth of the EGL and core was modeled by an isotropic expansion of the film in the bilayer structure.
To test the elastic wrinkling model, we conducted finite element (FE) simulations for bilayer structures with a film bonded on a substrate, which represents the EGL layer and core structure, respectively. The structures were assumed to be under 2D plane strain deformation to mimic the quasi-2D nature of cerebellum wrinkles. Neo-Hookean model was adopted to describe the elastic properties of both film and substrate, whose strain energy can be expressed as
where is the shear modulus and represents the first invariant of the right Cauchy-Green strain tensor. The Poisson’s ratios for the film and substrate were set to be 0.5, based on experimental observations that the bulk modulus of EGL and core are in the order of GPa, much larger than the shear modulus of soft tissues (~ kPa).
We carried out FE simulations through commercial software ABAQUS. A second order 6 node hybrid element (CPE6MH) was utilized to discretize the film and substrate. Very fine FE meshes were used to make sure the results independent of mesh size. To incorporate differential growth in real EGL layer and core, an isotropic growth deformation tension was applied to the modeled film by decoupling the deformation tenor into elastic deformation part and growth part .
For simplicity, we assume the growth part is isotropic and controlled by a scalar variable
where represents a faster growth in EGL than the core. To trigger instabilities in numerical simulations, random perturbations (e.g., White Gaussian noise with 0.001t mean magnitude) were applied to the nodal positions at the top surface of the film and the interface between the film and substrate.
To qualitatively understand the cut experiments we ran a FE simulation of a pre-cut circular bilayer structure and then assigned swelling strain to the film. This neglected the dynamical process in the real cut experiments and only focused on the final equilibrium of the cerebellum after long time relaxation. All the simulation parameters were the same as those in the wrinkling simulation. The initial cut length a is equal to 8 t. The minimum in-plane principal stress corresponds to the hoop stress in the film.
For a full treatment of the mathematics please see Engstrom et al. (2018).
We, formulated a two-dimensional model based on the parameters of a midsagittal section of the cerebellum. The distance of the outer edge of the EGL and, hence, the outer edge of the cerebellum from the center of the cerebellum was defined as r(θ) with θ as the angular coordinate. We assumed that r(θ) was single-valued. The thickness of the EGL was defined as t(θ). See model schematic in.
Taking into account the four assumptions discussed in the main text, we constructed the following energy functional to be minimized
with kr as the stiffness modulus (a spring constant in one-dimension) of the radial glial fibers and the pial surface contained in the meninges surrounding the cerebellum since the cerebellar radius is proportional to its perimeter, r0 as the preferred radius of the cerebellum, kt denoting a growth potential due to cell proliferation, t0 as thickness of the EGL (cortex), and, β quantified the mechanical resistance to changing the thickness of the EGL. Given our first assumption of an incompressible cerebellar core, we imposed the constraint
with A0 as a preferred cerebellar area. We applied the variational principle to minimize the energy functional subject to the core constraint, that is
where is a Lagrange multiplier. Assuming the preferred radius of the cerebellum is constant and the thickness of the EGL/cortex is also constant, then the preferred cerebellar shape was a circle and the EGL an annulus.
The variational analysis yielded the following equation of shape for t(θ);
with The solution to the equation of shape was
with C1 independent of θ and such that There was an additional equation of shape for r(θ) from the variational principle that depended on t(θ) and so was determined
We used the measured data at E16.5 to set the parameters to make predictions for the shape of both the EGL and core (and so the relationship between the two) at later times. Because we are primarily interested in shape changes, rather than size changes, a nondimensionalized model solution was used, that is we chose units where r0 = 1. This reduces the total number of parameters specifying the model to five dimensionless parameters. Plots assumed a circular preferred shape, and with other parameters as follows: is shown in Figure 6b,c, , , , and . Note that for , these parameters are numerically consistent with our E16.5 measurements: and , as well as the observed number of invaginations in the half circle: . All of these parameters are either constant or depend on the time-like parameter . One of these dependencies has a functional form that is physically justifiable (), another has a form that is biologically justifiable (), owing to the decrease in the number of radial glia over time.
We defined a dimensionless 'shape factor' as half of the perimeter divided by the square root of half of the area as appropriate for a semi-circle. To compare the model’s predictive deviation of this quantity form the semi-circular value we assumed a linear relationship between and time T measured in embryonic days: (T) = 0.3(T-15.5).
Statistical analyses were performed using Matlab software. Significance was determined at p<0.05. Two-way ANOVA was used for proliferation analysis as two variables were tracked, mouse and region. Cell shape, volume, fiber distribution, EGL thickness and bulk modulus were run under a standard ANOVA. After ANOVA analysis a multiple comparison was run with Tukey’s honestly significant difference criterion. F-test for variance and two-tailed student’s paired t-test were used for slice cutting and relaxation quantifications. The degrees of freedom, where appropriate, and P values are given in the figure legends. All error bars are standard deviations. No statistical methods were used to predetermine the sample sizes. We used sample sizes aligned with the standard in the field. No randomization was used nor was data collection or analysis performed blind.
Analysis and Design of Structural Sandwich Panels (1st ed)Oxford, New York: Pergamon Press.
A quantitative study of the human cerebellum with unbiased stereological techniquesThe Journal of Comparative Neurology 326:549–560.https://doi.org/10.1002/cne.903260405
Mechanical forces in cerebral cortical folding: a review of measurements and modelsJournal of the Mechanical Behavior of Biomedical Materials 29:568–581.https://doi.org/10.1016/j.jmbbm.2013.02.018
Rheological characterization of human brain tissueActa Biomaterialia 60:315–329.https://doi.org/10.1016/j.actbio.2017.06.024
Buckling without bending: a new paradigm in morphogenesisPhysical Review X, 8, 10.1103/PhysRevX.8.041053.
Shp2-dependent ERK signaling is essential for induction of Bergmann Glia and foliation of the cerebellumThe Journal of Neuroscience 34:922–931.https://doi.org/10.1523/JNEUROSCI.3476-13.2014
Mechanical control of tissue and organ developmentDevelopment 137:1407–1420.https://doi.org/10.1242/dev.024166
Embryonic stages in cerebellar afferent developmentCerebellum & Ataxias 2:7.https://doi.org/10.1186/s40673-015-0026-y
Fine-resolution maps of acoustic properties at 250 MHz of unstained fixed murine retinal layersThe Journal of the Acoustical Society of America 137:EL381–EL387.https://doi.org/10.1121/1.4916790
Improved High-Frequency ultrasound corneal biometric accuracy by Micrometer-Resolution Acoustic-Property maps of the corneaTranslational Vision Science & Technology 7:21.https://doi.org/10.1167/tvst.7.2.21
Autoregressive signal processing applied to High-Frequency acoustic microscopy of soft tissuesIEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control 65:2054–2072.https://doi.org/10.1109/TUFFC.2018.2869876
Differential tangential expansion as a mechanism for cortical gyrificationCerebral Cortex 24:2219–2228.https://doi.org/10.1093/cercor/bht082
Morphology, molecular codes, and circuitry produce the three-dimensional complexity of the cerebellumAnnual Review of Cell and Developmental Biology 23:549–577.https://doi.org/10.1146/annurev.cellbio.23.090506.123237
On the mechanism of mucosal folding in normal and asthmatic airwaysJournal of Applied Physiology 83:1814–1821.https://doi.org/10.1152/jappl.19220.127.116.114
Cerebellar granule cell replenishment postinjury by adaptive reprogramming of nestin+ progenitorsNature Neuroscience 20:1361–1370.https://doi.org/10.1038/nn.4621
Axons pull on the brain, but tension does not drive cortical foldingJournal of Biomechanical Engineering 132:071013.https://doi.org/10.1115/1.4001683
Cytodifferentiation of Bergmann Glia and its relationship with purkinje cellsAnatomical Science International 77:94–108.https://doi.org/10.1046/j.0022-7722.2002.00021.x
Bergmann glial development in the mouse cerebellum as revealed by tenascin expressionAnatomy and Embryology 194:223–234.https://doi.org/10.1007/BF00187133
Multimodal surface instabilities in curved Film–Substrate StructuresJournal of Applied Mechanics 84:081001–081013.https://doi.org/10.1115/1.4036940
Lilianna Solnica-KrezelReviewing Editor; Washington University School of Medicine, United States
Marianne E BronnerSenior Editor; California Institute of Technology, United States
In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.
Thank you for submitting your article "Cerebellar folding is initiated by mechanical constraints on a fluid-like layer without a cellular pre-pattern" for consideration by eLife. Your article has been reviewed by two peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Marianne Bronner as the Senior Editor. The reviewers have opted to remain anonymous.
The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.
Both reviewers thought that your study of cerebellar folding in developing murine embryo addresses an important outstanding question in the field of developmental neuroscience. They considered your studies of tissue dynamics and cellular behavior analyses, as well as mathematical simulations as comprehensive and compelling. The observations that during the cerebellar folding initiation differential expansion is due to the outer layer of proliferating progenitors expanding faster than the tissue core and that the folding occurs without an obvious cellular pre-pattern, are particularly interesting. These and other observations challenge the existing models and allow putting forward the multi-phase wrinkling model. The manuscript should be of broad interest to developmental biologists, cellular biologist, neuroscientists and applied mathematicians.
One of the reviewers is concerned that in the model the missing two parameters μ/kr are constrained to scale linearly in time. There must be a second constraint to fix two parameters. The motivation behind linear scaling and the full choice of parameters should be more transparent.
Many models were proposed to explain how brain folds. The elastic wrinkling models based on variations in stiffness, compressive forces, proliferation, and thickness has been used to explain the cerebellar folding. In this manuscript, Joyner and her colleagues challenged the folding models by measuring those values and cellar behaviors (e.g. cell shape, motility, and glial fibers) in mouse cerebella when the folding is initiated, and by mathematical simulation. They demonstrated via experimental data: (1) differential expansion and stiffness between the EGL vs. the core, (2) the EGL thickness at the anchoring centers (ACs) increases but not decreases during the folding initiation, (3) uniform proliferation in the EGL, (4) presence of horizontal and radial mechanical constrains, (5) fluid-like property of the EGL. With the mathematical simulation, they propose a novel multi-phase wrinkling model without introducing a cellular pre-pattern. They also explain the subdivisions after the initial folding (they termed hierarchical folding) when they consider ACs as mechanical boundaries. Although I am not familiar with mathematical simulation, I believe that they show enough evidence supporting their hypothesis. The no pre-pattern model was unexpected. The manuscript provides scientists in the research field of developmental biology and neuroscience with a novel idea for brain folding.
The work by Lawton et al. investigates the mechanism behind cerebellar folding by conducting a nicely detailed analysis of the dynamics and mechanical properties of the contributing tissues. What drives brain folding has excited many researchers, developmental biologists, physicists and applied mathematicians alike. Yet, a detailed quantification of the tissue dynamics and tissue properties has so far been lacking. The authors unfold in their very timely contribution how predictions from models guided their quantification of tissue dynamics to step by step rule out any previous models of brain folding and eventually deduct the key ingredients that they consider in a finite element model. They very carefully constrain parameters of their model to then arrive at a very convincing and consistent mechanism of how cerebellar folding unfolds.
My only major concern is that in their model the missing two parameters μ/kr are constrained to scale linearly in time. There must be a second constraint to fix two parameters. The motivation behind linear scaling and the full choice of parameters should be more transparent.https://doi.org/10.7554/eLife.45019.021
One of the reviewers is concerned that in the model the missing two parameters mμ/kr are constrained to scale linearly in time. There must be a second constraint to fix two parameters. The motivation behind linear scaling and the full choice of parameters should be more transparent.
We thank the reviewer for their comment on the parameterization of the model as it gives us a chance to clarify our approach. In the modeling community, it is quickly becoming the norm to reformulate a model in terms of dimensionless quantities, i.e. nondimensionalization. To do so, one works with dimensionless ratios of the parameters, thereby reducing the number of parameters to uncover a smaller set of quantities that the system depends on. In the manuscript, we also work with a nondimensionalized form of the model (i.e., we work in units where r0=1), and so only 5 dimensionless parameters are required to completely specify the model. These may be chosen as μ/kr, kr/kt, At/r0, t0/r0, and q, as discussed in the Materials and methods section entitled “Details of Multi-phase model […]”. Thus, the quantity μ/kr should not be thought of as a ratio of two independent parameters, but rather a single dimensionless parameter. We indeed should have been more transparent about using the nondimensionalized version of the model, and we thank the reviewer for pointing out this oversight. The revised manuscript has language added to the section "Details of multi-phase model […]" to state explicitly that the nondimensionalized model is used to make the plots in Figure 6B-C as follows:
“Because we are primarily interested in shape changes, rather than size changes, a nondimensionalized model solution was used, i.e., we chose units where r0=1. This reduces the total number of parameters specifying the model to five dimensionless parameters.”
In re-reading this section, we also noticed a typo that may have been the source of some of the confusion. The fifth paragraph said "Figure 5F" where it should say "Figure 6B-C". We apologize for this typo, we have fixed it in the revised manuscript, and hope this clarifies how the parameters in Figure 6B-C are chosen.
Regarding time-dependence, the model parameter μ/kr should be increasing over time because kr is decreasing over time, which we suggest in the Discussion can be attributed to the radial glia transitioning to Bergmann glia (see, e.g., Discussion, third paragraph). We assume linear scaling, because that is the simplest possible kind of scaling. That this is an assumption stated in both the legend to Figure 6 and in the section "Details of multi-phase model […]". Having said that, Figure 6C illustrates that this assumption is not too unreasonable, at least concerning the agreement of the experimental and theoretical shape index.https://doi.org/10.7554/eLife.45019.022
- Andrew K Lawton
- Jonathan Mamou
- Teng Zhang
- J M Schwarz
- J M Schwarz
- Alexandra L Joyner
- Alexandra L Joyner
- Alexandra L Joyner
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
We are grateful to Anna-Katerina Hadjantonakis and Songhai Shi for use of Imaris software. We thank Jennifer Zallen and Anna-Katerina Hadjantonakis for experimental advice, Nathanael Kim for help with the acoustic microscopy and Professor Tadashi Yamagushi for his support for MO’s visit to New York. We appreciate the discussions we have had with Alexandre Wojcinski and the entire Joyner Laboratory, and the administrative support from Cara Monaco.
Animal experimentation: All experiments were performed following protocols approved by Memorial Sloan Kettering Cancer Center's Institutional Animal Care and Use Committee. Protocol number: 07-01-001 Mouse Developmental Genetics
- Marianne E Bronner, California Institute of Technology, United States
- Lilianna Solnica-Krezel, Washington University School of Medicine, United States
- Received: January 9, 2019
- Accepted: March 30, 2019
- Version of Record published: April 16, 2019 (version 1)
© 2019, Lawton et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.