Neuroevolutionary evidence for a universal fractal primate brain shape
Abstract
The cerebral cortex displays a bewildering diversity of shapes and sizes across and within species. Despite this diversity, we present a universal multiscale description of primate cortices. We show that all cortical shapes can be described as a set of nested folds of different sizes. As neighbouring folds are gradually merged, the cortices of 11 primate species follow a common scalefree morphometric trajectory, that also overlaps with over 70 other mammalian species. Our results indicate that all cerebral cortices are approximations of the same archetypal fractal shape with a fractal dimension of d_{f} = 2.5. Importantly, this new understanding enables a more precise quantification of brain morphology as a function of scale. To demonstrate the importance of this new understanding, we show a scaledependent effect of ageing on brain morphology. We observe a more than fourfold increase in effect size (from two standard deviations to eight standard deviations) at a spatial scale of approximately 2 mm compared to standard morphological analyses. Our new understanding may, therefore, generate superior biomarkers for a range of conditions in the future.
eLife assessment
This study presents valuable framework and findings to our understanding of the brain cortex as a fractal object. Based on detailed methodology, the evidence provided on the stability of its shape property within 11 primate species is convincing, as well as the scalespecific effects of ageing on the human brain. This study will be of interest to neuroscientists interested in brain morphology, and to physicists and mathematicians interested in modeling the shapes of complex objects.
https://doi.org/10.7554/eLife.92080.4.sa0eLife digest
Many of the brain’s essential functions – from decisionmaking to movement – take place in its outer layer known as the cerebral cortex. The shape of the cerebral cortex varies significantly between species. For instance, in humans, it is folded in to grooves and ridges, whereas in other animals, including mice, it is completely smooth. The structure of the cortex can also differ within a species, and be altered by aging and certain diseases.
This vast variation can make it difficult it to characterize and compare the structure of the cortex between different species, ages and diseases. To address this, Wang et al. developed a new mathematical model for describing the shape of the cortex.
The model uses a method known as coarse graining to erase, or ‘melt away’, any cortical folds or structures smaller than a given threshold size. As this threshold increases, the cortex becomes progressively smoother. The relationship between surface areas and threshold sizes indicates the fractal dimension – that is, how fragmented the cortex is across different scales.
Wang et al. applied their model to the brain scans of eleven primates, including humans, and found the fractal dimension of the cortex was almost exactly 2.5 for all eleven species. This suggests that the cortices of the different primates follow a single fractal shape, which means the folds of each cortex have a similar branching pattern. Although there were distinctions between the species, they were mainly due to the different ranges of fold sizes in each cortex. The model revealed that the broader the range of fold sizes, the more folded the brain – but the fractal pattern remains the same.
The brain melting method created by Wang et al. provides a new way to characterise cortical shape. Besides revealing a hitherto hidden regularity of nature, they hope that in the future their new method will be useful in assessing brain changes during human development and ageing, and in diseases like Alzheimer’s and epilepsy.
Introduction
The morphological complexity of the mammalian cerebral cortex has fascinated scientists for generations, with cortices across and within species, exhibiting a large diversity of shapes and sizes. Such diversity is not arbitrary, however. The mammalian brain folds into stereotypical, hierarchicallyorganised structures such as lobes and major gyri. In fact, qualitative and quantitative regularities in cortical scaling have often been suggested and observed (Zhang and Sejnowski, 2000; Francis et al., 2009; Karbowski, 2011; Mota and HerculanoHouzel, 2014). More specifically, through modelling the mechanism of cortical folding from a statistical physics approach, we have previously derived a theoretical scaling law relating pial surface area $A}_{t$, exposed surface area $A}_{e$ (the exposed surface area can be thought of as the surface area of a piece of cling film wrapped around the brain; mathematically, for the remaining paper it is the convex hull of the brain surface), and average cortical thickness $T$:
This scaling law, relating powers of cortical thickness and surface area metrics, was shown to be valid across mammalian species (Mota and HerculanoHouzel, 2015) and within the human species (Wang et al., 2016), as well as to the structures and substructures of individual brains (Wang et al., 2019; Leiberg et al., 2021). Notably, across all these cases the dimensionless offset $k$ is shown to be near invariant. However, this universality, presumed to be arising from universal physical principles and evolutionarily conserved biomechanics, says little about what that cortical shape actually is, beyond a constraint binding three morphometric parameters. In this paper, we take Equation 1 as an empirical starting point to create a new and hierarchical way of expressing cortical shape. Specifically, we introduce a coarsegraining procedure that renders the cortex at different spatial scales, or resolutions. We show that coarsegrained primate cortices at each spatial scale can be understood as approximations of the same universal selfsimilar archetypal form, of which the observed scaling law Equation 1 can then be shown to be a direct consequence.
Besides revealing a symmetry in nature hidden under much apparent complexity, our results indicate a conservation of morphological relationships across evolution. We will show that these results further provide us with a new and powerful tool to express and analyse cortical morphology. As an example, we will calculate the effects of human ageing across spatial scales and show that the effects are highly scaledependent.
Mathematical background
The universal scaling law Equation 1 can be rewritten in a suggestive way
where the $A}_{0}=\frac{{T}^{2}}{{k}^{4}$ is a fundamental area element that defines the threshold between gyrencephaly (folded cortex) and lissencephaly (smooth cortex) when $A}_{t}={A}_{e}={A}_{0$. For a constant $k$ the value of $A}_{0$ is a multiple of $T}^{2$, indicating that cortical thickness determines the size of the smallest possible gyri and sulci.
This rewriting highlights a new perspective, or interpretation of the scaling law: it now suggests a relationship between intrinsic and extrinsic measures of cortical size (given the folded laminar structure of the cortex, areas are the more natural way of measuring its ‘size’), $A}_{t$ and $A}_{e$, respectively, measured in units of $A}_{0$. This is reminiscent of fractal scaling (Mandelbrot, 1983), where a complex shape reveals ever smaller levels of selfsimilar detail as it is probed in ever smaller scales (or equivalently, higher resolutions), represented here by $A}_{0$. The scaling, or power exponent between the measured intrinsic and extrinsic sizes is the socalled fractal dimension.
Although actual fractals are mathematical abstractions, they can often be defined as the limit of iterative processes. Many structures in nature, and in particular biology (Elston and Zietsch, 2005; Codling et al., 2008; Ionescu et al., 2009; Losa, 2011; Klonowski, 2016; Di Ieva, 2016; Reznikov et al., 2018), are good approximations of a fractal. Equation 2 is suggestive, but not proof, that cortices are among these forms, with a fractal dimension of $1.25\times 2=2.5$ (the factor 2 being the topological dimension of areas). Indeed, fractal scaling for various aspects of cortical morphology has often been postulated (Free et al., 1996; Kiselev et al., 2003), with a number of recent papers making use of MRI data (Marzi et al., 2021; Jao et al., 2021; Meregalli et al., 2022; Díaz Beltrán et al., 2024). Most recently published estimates of fractal dimension for the whole cortex are indeed close to 2.5 (King et al., 2010; Madan and Kensinger, 2016; Madan and Kensinger, 2017; Marzi et al., 2020).
Here, for the first time, we propose to directly construct morphologically plausible realisations of cortices at any specified spatial scale, or resolution. This is achieved through a coarsegraining method that removes morphological details smaller than a specified scale while preserving surface integrity. For example, at a set scale of 3 mm, sulcal walls that are less than 3 mm apart would be removed, and the neighbouring gyri would be fused. This method is a new systematic way of obtaining shape properties from the cortex in terms of a sequence of morphometric measurements as spatial scale varies. By examining how areas scale across coarsegrained versions of actual primate cortices, we will be able to directly verify cortical selfsimilarity.
Method
Coarsegraining method
As a starting point for a coarsegraining method, we suggest to turn to a wellestablished method that measures a fractal dimension of objects: the socalled boxcounting algorithm (). Briefly, this algorithm fills the object of interest (the cortex in our case) with boxes, or voxels of increasingly larger sizes and counts the number of boxes in the object as a function of box size. As the box size increases, the number of boxes decreases; and in a loglog plot, the slope of this relationship indicates the fractal dimension of the object. In our case, this method would not only provide us with the fractal dimension of the cortex, but, with increasing box size, the filled cortex would also contain less and less detail of the folded cortex. Intuitively, with increasing box size, the smaller details below the resolution of a single box would disappear first, and increasingly larger details will follow – precisely what we require from a coarsegraining method. We, therefore, propose to expand the traditional boxcounting method beyond its use to measure the fractal dimension, but to analyse the reconstructed cortices as different realisations of the original cortex at the specified spatial scale.
Concretely, our proposed method requires the bounding pial and white matter surfaces of the cortical ribbon as input. We obtained these surfaces based on reconstructions from magnetic resonance imaging data in 11 different primate species. Algorithmically, we then segment the space between the original pial and white matter surfaces into a 3D grid of boxes of the desired scale $\lambda$, where each box is a cube of dimensions $\lambda \times \lambda \times \lambda$. We also term the 3D grid of cubes ‘voxelisation,’ as it effectively captures the cerebral cortex as voxels in 3D space (Figure 1B bottom row). At any given scale, or voxel size, this process effectively erases morphological features (folds) that are smaller than the cube size. Visually, increasing the voxel size appear as if the cortex is ‘melting’ and ‘thickening’ (Figure 1, and videos: https://bit.ly/3CDoqZQ).
A more technical and detailed description and discussion of the algorithm is provided in Appendix 1. Note this method has also no direct dependency on the original MR image resolution, as the inputs are smooth grey and white matter surface meshes reconstructed from the images using strong (bio)physical assumptions and, therefore, containing more finegrained spatial information than the raw images (see also Appendix 2).
Rescaling coarsegrained outputs for analysis
Morphological properties, such as cortical thicknesses measured in our ‘melted’ brains are to be understood as a thickness relative to the size of the brain. Therefore, to analyse the scaling behaviour of the different coarsegrained realisations of the same brain, we apply an isometric rescaling process that leaves all dimensionless shape properties unaffected (more details in Appendix 3.1). Conceptually, this process fixes the voxel size, and instead resizes the surfaces relative to the voxel size, which ensures that we can compare the coarsegrained realisations to the original cortices, and test if the former, like the latter, also scales according to Equation 1. Resizing, or more precisely, shrinking the cortical surface is mathematically equivalent to increasing the box size in our coarsegraining method. Both achieved an erasure of folding details below a certain threshold. After rescaling, as an example, the cortical thickness also shrinks with increasing levels of coarsegraining, and never exceeds the thickness measured at native scale.
Independent morphological measures of shape
To better characterise the coarsegrained cortices in terms of their similarity in offset, we use a previously introduced (Wang et al., 2021) a set of independent measures, $K$, $I$, and $S$, that summarise the morphometry of the cortex in a natural and statistically robust way. In this framework, isometrically scaled copies of the same morphometry all map onto a line along the $I=\mathrm{log}{A}_{t}+\mathrm{log}{A}_{e}+\mathrm{log}{T}^{2}$ direction, which is perpendicular to a $K\times S$ plane that fully summarises their shape. $K=\mathrm{log}{A}_{t}\frac{5}{4}\mathrm{log}{A}_{e}+\frac{1}{4}\mathrm{log}{T}^{2}$ is the direction defined by the offset $k$ of the scaling law Equation 1, while direction $S=\frac{3}{2}\mathrm{log}{A}_{t}+\frac{3}{4}\mathrm{log}{A}_{e}\frac{9}{4}\mathrm{log}{T}^{2}$ captures the remaining information about shape, and can be regarded as a simple measure of morphological complexity. Wang et al., 2021 provides a detailed derivation and demonstration of the superior sensitivity and specificity of these new morphometric measures. The advantage of using this framework here is that we can assess the offset $K$ (and shape term $S$) without interference by isometric size effects, including any rescaling procedures.
Data and processing
With the exception of the marmoset data, all other cortical surface reconstructions were based on healthy individual brains.
Human data
To study healthy human adults, we used the Human Connectome Project (HCP) MRI data, available at https://db.humanconnectome.org/; Van Essen et al., 2012, obtained using a 3T Siemens Skyra scanner with 0.7 mm isotropic voxel size. We used the HCP minimally preprocessed FreeSurfer data output, which provided the pial and white matter surface meshes we required. We selected five random subjects in the age category 22–25 y.o. and show one example subject (103414) in the main text, and the remaining subjects in Appendix 3.
To study the alterations associated with human ageing, we used T1 and T2 weighted MRI brain scans from The Cambridge Centre for Ageing and Neuroscience (CamCAN) dataset (available at http://www.mrccbu.cam.ac.uk/datasets/camcan/ Shafto et al., 2014; Taylor et al., 2017). CamCAN used a 3T Siemens TIM Trio System with 1 mm isotropic voxel size (for more details see Shafto et al., 2014; Taylor et al., 2017). From the CamCAN dataset we retained 644 subjects that successfully completed preprocessing (with Freesurfer reconall) without errors. From these subjects, we selected all subjects between the ages of 17–25 inclusive (forming the 20 y.o. cohort, n=27); we also selected all subjects between the ages of 77–85 inclusive (forming the 80 y.o. cohort, n=86).
To confirm the ageing results, we also obtained an independent dataset from the Nathan Kline Institute (NKI)/Rockland sample (Nooner et al., 2012; http://fcon_1000.projects.nitrc.org/indi/pro/nki.html) using the same procedure as described for the CamCAN dataset.
The MR images of both CamCAN and NKI datasets were first preprocessed by the FreeSurfer 6.0 pipeline reconall, which extracts the greywhite matter boundary as well as the pial surface. These boundaries were then quality checked by visual inspection for particularly the young and old cohorts and manually corrected where needed.
For all three datasets, we obtained the pial and white matter surfaces for further analysis. In the current work, the analysis is always hemispherebased, as in our previous work (Mota and HerculanoHouzel, 2015; Wang et al., 2016). We did not perform a more regionalised analysis, which is also possible (Wang et al., 2019; Leiberg et al., 2021). Future work using the principle demonstrated here can be directly extended to derive regionalised measures across scales.
Nonhuman primate data
Macaque
Rhesus Macaque MRI scans were carried out at the Newcastle University Comparative Biology Centre. Macaques were trained to be scanned while awake and sat in a primate chair. Both T1weighted MPRAGE and T2weighted RARE sequences were acquired, using a vertical MRI scanner (Biospec 4.7 Tesla, Bruker Biospin, Ettlingen, Germany).
Scans were processed using a custom macaque MRI pipeline, incorporating ANTs, SPM, FreeSurfer, and FSL. Briefly, this involved the creation of precursor mask in SPM, denoising (ANTs DenoiseImage) and debiasing (ANTs N4BiasFieldCorrection, and Human Connectome BiasFieldCorrection script), creation of a final mask in SPM and then processing in FreeSurfer (using a modified version of the standard FreeSurfer processing pipeline). Reconstructed pial and white matter surfaces were visually qualitycontrolled in conjunction with the MR images.
Marmoset
The marmoset MRI structural scan was collected as part of the development of the NIH marmoset brain atlas (Liu et al., 2018; Liu et al., 2020). Data was collected ex vivo from a 4.5yearold male marmoset using a T2star weighted 3D FLASH sequence using a horizontal MRI scanner (Biospec 7 Tesla, Bruker Biospin, Ettlingen, Germany).
A total of 10 scans were collected and averaged into one final image. In combination with scans from other modalities, cortical boundaries were manually delineated on each coronal slice. Boundaries were then refined through comparisons with other atlases. Volumetric data was then converted to surfaces using a custom pipeline involving an intermediate generation of highresolution mesh data (Madan, 2015), decimation (Madan, 2016), and remeshing. Reconstructed pial and white matter surfaces were visually quality controlled in conjunction with the MR images.
Other nonhuman primates (NHPs)
The remaining NHP MRIs and subsequent brain surface extraction are detailed in Ardesch et al., 2022; Bryant et al., 2021, and provided to the authors in a processed format. Briefly, a range of specialised scanners were used to acquire optimal images for each species. FreeSurfer 6.0 with some modifications was used for surface reconstruction, complemented by FSL, ANTS, and Matlab. All surfaces were visually inspected for accuracy and consistency across datasets.
Comparative neuroanatomy data
The comparative neuroanatomy dataset for different mammalian species is the same as previously published (Mota and HerculanoHouzel, 2015). Note that for this dataset, we only had numerical values for the total and exposed surface area, as well a average cortical thickness estimates. We did not perform any analysis across scales in this dataset (hence surfaces were not required), but only used it as a reference dataset.
Ethics statement
All analyses were performed on anonymised data that were acquired previously as part of other studies/consortia with ethical approval from Newcastle University (reference: 22/SC/0016).
Statistical analyses
Briefly, linear regression is used in either a mixedeffect model to capture effects across individuals and species, or in simple fixedeffect settings to estimate regression slopes to obtain fractal dimension.
In the final part of the Results, we analyse the effect between a group of 20yearolds and 80yearolds. Effect size is calculated as Cohen’s D between the two groups.
Throughout the paper, statistical significance is not a crucial argument, and we report pvalues only for reference and completeness.
More details can be seen in the analysis code, from which the reader can directly reproduce all the main result figures.
Results
All primate brains follow the same scaling law across spatial scales
We have analysed cortices of 11 different primate species: Graybellied Night Monkey (Aotus lemurinus), Tufted Capuchin Monkey (Cebus apella), Blackandwhite Colobus (Colobus guereza), Senegal Bushbaby (Galago senegalensis), Woolly Monkey (Lagothrix lagotricha), Graycheeked Mangabey (Lophocebus albigena), Rhesus Macaque (Macaca mulatta), Common Marmoset (Callithrix jacchus), Chimpanzee (Pan troglodytes), Whitefaced Saki (Pithecia pithecia), and various cohorts of human subjects. We applied our coarsegraining procedure to their pial and white matter surfaces, and empirically determined (i) that all species followed a power law (linear regression ${R}^{2}>0.999$ for all species); (ii) the slope of said power law is $\alpha =1.255$ on a group level (CI: [1.254 1.256]) using linear mixed effect modelling (see Figure 2 for visualisation, and Appendix 3 for a detailed breakdown by species); and importantly, (iii) all species also show a similar offset $\mathrm{log}k\approx 0.65263$, with a standard deviation of intercept across species estimated at 0.02 from linear mixed effect modelling in $\mathrm{log}k$.
Taken separately, the scaling for each species is proof that their cortices are selfsimilar with the same scaling: they each approximate a fractal with fractal dimension ${d}_{f}=2.5$. Considering all species together, different species also overlap substantially (similar offset), and only differ from each other in the range of scales over which the approximation is valid (see Appendix 3.2 and Appendix 4). Thus, as Figure 2 illustrates, the data supports a universal scaling law across primate species and spatial scales:
with $\displaystyle k=\mathrm{0.2277.}$
Primate brains at different spatial scales are morphometrically similar to each other and other mammalian species
To better characterise the coarsegrained cortices in terms of their similarity in offset, we use a set of independent morphometric measures, $K$, $I$, and $S$, that summarise the morphometry of the cortex in a natural and statistically robust way. We can, therefore, assess the offset $K$ (and shape term $S$) without interference by the isometric size or rescaling.
We can measure $K$ and $S$ for any object, but a fuller expression is captured by the trajectory of said object as a function of coarsegraining in the $K\times S$ plane. This is a very convenient and informative way of summarising an object: selfsimilar objects correspond to straight trajectories as the $K\times S$ plane is in loglog space. In particular, objects without any folds or protrusion (i.e. convex, such as the box with finite thickness in Figure 3A) corresponds to the line $K=\frac{1}{9}S$, as $A}_{e}={A}_{t$ for all levels of coarsegraining. Horizontal trajectories (constant $K$) represent fractal objects with fractal dimension ${d}_{f}=2.5$ (Figure 3B). And finally, in the $K\times S$ plane, a group of objects can said to be ‘universal’ when their trajectories overlap, so that they can all be regarded as coarsegrained versions of one another (Figure 3C).
Primate cortices (Figure 3D) display a nearly invariant $K$ in all cases. But, over all levels of coarsegraining, $K$ also remains nearinvariant in all trajectories as $S$ decreases, resulting in a set of horizontal lines that largely overlap with each other and other mammalian species. The variance in $K$ across scales and all 11 species is $<0.01$, which is at least an order of magnitude lower than the variance in $S$. Primate brains, therefore, have all three characteristics of selfsimilarity, fractality (with ${d}_{f}=2.5$), and universality (invariant $K$ for all scales and species) at the same time.
Thus, coarsegrained primate cortices are morphometrically similar to, and in terms of the universal law, ‘as valid as’ actual existing mammalian cortices. Note, of course, that the coarsegrained brain surfaces are an output of our algorithm alone and not to be directly/naively likened to actual brain surfaces, e.g., in terms of the location or shape of the folds. Our comparisons here between coarsegrained brains and actual brains is purely on the level of morphometrics across the whole cortex. In contrast, we tested various nonbrain objects, and while e.g., the walnut, and bell pepper form (partially) straight lines, they vary in both K and S (see Appendix 5). These objects may have a fractal regime, but their fractal dimension is not 2.5, nor are they similar to primate or mammalian brains in terms of $K$ or $S$. Furthermore, Appendix 6 underscores the algorithmic and statistical robustness of these results using multiple realisations of the coarsegraining procedure on the same object. In this framework, our main result can thus be expressed simply: for all the cortices we analysed, and for none of the noncortices, coarsegraining will leave $K$ largely unaffected, while morphological complexity $S$ will decrease.
Morphometric measures as functions of scale reveal scalespecific effects of ageing
In the final part of our work, we show how our algorithm and the associated new understanding of brain morphology may become useful in applications. As an example, we will focus on how the ageing process affects human cortical morphology across scales. In Figure 4A, we compare the total surface area ${A}_{t}(\lambda )$ as a function of scale $\lambda$ for a young (20 yearold) vs. an old (80 yearold) group of human brains. The difference in $A}_{t$ between the groups takes a U shape, and the strongest effect is seen at approximately 2 millimeters (greatest effect size of –8.635 seen at scale 2.188 mm), where older subjects have higher $A}_{t$. For scales over ∼5 mm, and under ∼0.5 mm the differences become relatively small, suggesting the ageing process has less effect on the largest and smallest cortical morphological features. Finally, we reproduced these results in an independent dataset in Appendix 7.2.
In this particular example, the scaledependency of morphological measures can be visually and intuitively understood by looking at the reconstructed surfaces at each scale: Figure 4B shows some coronal slices of the cortical surface. On a scale of 0.27 mm, the gyri in the younger subjects are densely packed, but the older subjects show the expected widening between gyral walls and decrease in gyral surface area at the crown (see e.g. Jin et al., 2018; Madan, 2019; Kochunov et al., 2008 for recent investigations and references therein). At 1.86 mm, the younger cortices have already partially ‘melted,’ erasing most small sulci between the densely packed gyri. In the older humans, however, the gyri are less dense, the sulci is more open, and thus, at 1.86 mm, most gyri and sulci have not been erased yet. At scale 7.94 mm, both young and old brains have ‘melted’ down to similarly nearlissencephalic cortices. For a more detailed multiscale investigation over the entire human lifespan, please refer to our new preprint (Leiberg et al., 2024).
More broadly, one can regard the melting process as a way of determining how cortical area is allocated across different scales. As the cortex ‘melts,’ the contributions to the total area from features smaller than the cutoff scale are eliminated. In the example in Figure 4, for instance, we can say about half (${10}^{54.7}={10}^{0.3}\approx 2$) of the total area in the 80 y.o. cortices are present in features smaller than 4 mm.
Discussion
We have devised a new way of expressing the morphology of the mammalian cerebral cortex, as the flow in the values of morphometric measures over a range of spatial scales. This was achieved by coarsegraining cortical surfaces, erasing morphological features smaller than the specified scale, while preserving surface integrity. After applying this method to the cortices of 11 primate species, we have shown that all these diverse cortices are approximations of the same archetypal selfsimilar (fractal) shape. Most of their morphological diversity can be ascribed to the speciesspecific ranges of spatial scales over which each approximation is valid, with the smallest scale being an invariant multiple of cortical thickness. This was a proofbyconstruction of fractality, a step beyond the usual boxcounting approach, which also yields scaledependent morphometrics. As a proofofprinciple we showed that healthy human ageing has highly scaledependent effects in a range of morphometrics.
Advantages and advances
Compared to previous literature, we can summarise our main contribution and advance as follows: (i) We are showing for the first time that representative primate species follow the exact same fractal scaling – as opposed to previous work showing that they have a similar fractal dimension (Hofman, 1985; Hofman, 1991), i.e., slope, but not necessarily the same offset, as previous methods had no consistent way of comparing offsets. (ii) Previous work could also not show direct agreement in morphometrics between the coarsegrained brains of primate species and other nonprimate mammalian species. (iii) Demonstrating in proofofprinciple that multiscale morphometrics, in practice, can have much larger effect sizes for classification applications. This moves beyond our previous work where we only showed the scaling law across (Mota and HerculanoHouzel, 2014) and within species (Wang et al., 2016), but all on one (native) scale with comparable effect sizes for classification applications (Wang et al., 2016).
In simple terms: we know that objects can have the same fractal dimension, but differ greatly in a range of other shape properties. However, we demonstrate here, that representative primate brains and mammalian brains indeed share a range of other key shape properties, on top of agreeing in fractal dimension. This suggests a universal blueprint for mammalian brain shape and a common set of mechanisms governing cortical folding. As a practical additional outcome of our study, we could show that our novel method of deriving multiscale metrics can differentiate subtle morphological changes much better (four times the effect size) than the metrics we have been using so far at a single native scale.
Expressing cortical morphology as a function of scale is more detailed than a list of summary morphometric measures, and more informative than the mere listing of every sulcus and gyrus. We propose this new syntax as the basis for a more rigorous characterisation of brain morphology and morphological changes. A clear advantage is that some biological processes may only act on a specific spatial scales, leaving other scales untouched (ageing in our example). By disambiguation of the spatial scale, it allows for an extra dimension of understanding, and tracking of biological processes. In the future, this approach can also be extended to cortical development and to various degenerative (Wang et al., 2016) and congenital (Wang et al., 2021) neuropathic conditions, especially if combined with a regionalised versions of this method applied to specific cortical regions (Wang et al., 2019) or local patches (Leiberg et al., 2021).
Implications of universality
Empirically, the main result of this paper is the demonstration of a universal selfsimilar scaling (Equation 3) for primate, and presumably mammalian cortices. There are two aspects of this universality: first that for each and every cortex the value of $K$remains the same for all scales as one removes substructures smaller than a varying length scale (or equivalently, that the fractal dimension is almost exactly 2.5 in all cases). Second, that the value of $K$ for the cortices of different species is approximately the same, as previously observed (Equation 1) across species (Mota and HerculanoHouzel, 2015) and individuals (Wang et al., 2016). One could imagine a set of objects for which one but not the other aspect of the universality in $K$ holds true. However, the fact that both universalities hold true is significant. It suggests the existence of a single highly conserved mechanism for cortical folding, operating on all length scales selfsimilarly with only a few morphological degrees of freedom. It also hints at the possibility of deriving cortical scaling from some variational principle. Finally, this dual universality is also a more stringent test for existing and future models of cortical gyrification mechanisms at relevant scales, and one that moreover is applicable to individual cortices. For example, any models that explicitly simulate a cortical surface as an output could be directly coarsegrained with our method and the morphological trajectories can be compared with those of actual humans and primate cortices. The simulated cortices would only be ‘valid’ in terms of the dual universality, if it also produces the same morphological trajectories (Note, we do not suggest to directly compare coarsegrained brain surfaces with actual biological brain surfaces. As we noted earlier, the coarsegrained brain surfaces are an output of our algorithm alone and are not to be directly/naively likened to actual brain surfaces, e.g., in terms of the location or shape of the folds. Our comparisons here between coarsegrained brains and actual brains is purely on the level of morphometrics across the whole cortex).
The scaling itself does not imply or favour any particular proposed gyrification model (ours Mota and HerculanoHouzel, 2015 included), and all results in this paper are agnostic about this choice. Indeed, our previously proposed model (Mota and HerculanoHouzel, 2015) for cortical gyrification is very simple, assuming only a selfavoiding cortex of finite thickness experiencing pressures (e.g. exerted by white matter pulling, or by CSF pressure). The offset $K$, or ‘tension term,’ precisely relates to these pressures, leading us to speculate that subtle changes in $K$ correlate with changes in white matter property (Wang et al., 2016; Wang et al., 2021). In the same vein of speculation, the scaledependence of $K$ shown in this work might, therefore, be related to different types of white matter that span different length scales, such as superficial vs. deep white matter, or Ufibres vs. major tracts. However, there are also challenges to the axonal tension hypothesis (Xu et al., 2010). Indeed, white matter tension differentials in the developed brain may not explain the location of folds, but instead white matter tension may contribute to a wholebrain scale ‘pressure’ during development that drives the folding process overall. Aside from speculations about the biological interpretation, the simplicity of the highlighted scaling law parallels many complex phenomena in nature that displays simple and universal scaling that can be derived from first principles (Barenblatt, 1996; West et al., 1997; Gagler et al., 2022). In addition, recent results suggest simplicity and symmetry are generically favoured on statisticalensemble grounds by evolution (Johnston et al., 2022). Our model correctly predicts the scaling law (Equation 1), but a more complete explanation for cortical gyrification is probably far more complex (Quezada et al., 2020) than can be accounted by such a simple model.
One specific example of said complexity is the exact patterns, locations, depth, and features of gyri and sulci. We know such patterns to be, for example, variable but also somewhat heritable in humans, whilst in macaques such patterns are relatively preserved across the species. Our work does not explain any of these observations, nor are the coarsegrained versions of human brains is supposed to exactly resemble the location/pattern/features of gyri and sulci of other primates. The similarities we highlighted here are on the level of summary metrics, and our goal was to highlight the universality in such metrics points towards highly conserved quantities and mechanisms.
Biological plausibility and implications
The observation that with increasing voxel sizes, the coarsegrained cortices tend to be smoother and thicker is particularly interesting: the scaling law in Equation 3 can be understood as thicker cortices ($T$) form larger folds (or are smoother i.e. with less surface area $A}_{t$) when brain size is kept constant ($A}_{e$). This way of understanding has also been vividly illustrated by using the analogy of forming paper balls with papers of varying thickness in Mota and HerculanoHouzel, 2014: comparing two paper balls of the same size ($A}_{e$) will show that the one that uses thicker paper ($T$) will be smoother, have larger folds and a smaller total surface area ($A}_{t$), in comparison with the one using thinner paper. The scaling law can, therefore, be understood as a physically and biologically plausible statements and our algorithm yields results in line with the scaling law.
More broadly, the interaction between brain development and evolution may also benefit from a scalespecific understanding. This may be important in elucidating what in cortical morphology is selected for by evolution, what is determined by physics; what is specified by genes, and what is emergent. For example, one can estimate the number of structural features at each scale (${A}_{t}(\lambda )$ as multiples of ${A}_{0}(\lambda )$), and it will be interesting to correlate this number to other quantifiers of cortical structure, such as the number of neurons (Mota and HerculanoHouzel, 2014), the number of functional areas (Molnár et al., 2014), or the number of cortical columns (Kaas, 2012), possibly over different stages of development. Generally, the larger the cortical feature (i.e. from gyri to functional areas to lobes to hemispheres), the earlier during development it appears (Zilles et al., 2013; Garcia et al., 2018), and more broadly it is conserved over kinship (Pizzagalli et al., 2020) and phylogeny (Heuer et al., 2019; Valk et al., 2020). It thus seems likely that comparative neuroanatomical methods (Mars et al., 2014; Croxson et al., 2018) may be directly used to identify and contrast structures in coarsegrained cortices of more highly gyrified species with their analogues in lessgyrified species. One could then perhaps specify when evolution conserves and when it invents old and new cortical features.
From an application perspective, our final result illustrates clearly that the surface area difference between older and younger subjects at the ‘native’ scale (i.e. original free surfer surfaces) is negligible (effect size smaller than two standard deviations). However, in our analysis across scales, there is a clear optimal scale at ∼2 mm where the effect size is maximised between older and younger subjects (effect size is –8 standard deviations). For most classification applications in biology and medicine, the increased effect size and hence separability of groups in the scaledependent morphometrics represent a huge advance over the native scale.
Outlook
Our work here was limited to summary descriptors of entire cortical hemispheres, but future work will explore extensions of these methods to lobes and cortical areas, similarly to Wang et al., 2019; Leiberg et al., 2021. This will generate precise characterisations of the morphological differences between phylae and across developmental stages, and perhaps pinpoint the time and location of morphological changes leading to congenital and neurodegenerative conditions. Ultimately, we hope this new framework for expressing and analysing cortical morphology, besides revealing a hitherto hidden regularity of nature, can become a powerful tool to characterise and compare cortices of different species and individuals, across development and ageing, and across health and disease.
Appendix 1
Coarsegraining algorithm
The core coarsegraining algorithm underpinning our analyses takes pial and white matter surfaces from the cortical surface reconstruction as inputs and at any specified spatial scale $\lambda$ ‘voxelises’ the grey matter ribbon at the specified resolution. To achieve this, we set up a 3D voxel grid, where each grid voxel is of size $\lambda \times \lambda \times \lambda$.
In detail, we assign all voxels in the grid with at least four corners inside the original pial surface to the pial voxelisation. This process allows the exposed surface to remain approximately constant with increasing voxel sizes. A constant exposed surface is desirable, as we only want to gradually ‘melt’ and fuse the gyri, but not grow the bounding/exposed surface as well. We want the extrinsic area to remain approximately constant as we decrease the intrinsic area via coarsegraining; it is like generating iterates of a Koch curve in reverse, from more to less detailed, by increasing the length of the smallest line segment.
We then assign voxels with all eight corners inside the original white matter surface to the white matter voxelisation. This is to ensure integrity of the white matter, as otherwise white matter voxels in gyri may become detached from the core white matter, and thus artificially increase white matter surface area. Indeed the main results of the paper are not very sensitive to this decision using all eight corners, vs. e.g., only four corners, as we do not directly use white matter surface area for the scaling law measurements. However, we still maintained this choice in case future work wants to make use of the white matter voxelisations or derivative measures.
Finally, subtracting the white matter voxelisation from the pial voxelisation, we obtain the grey matter voxelisation (Figure 1B bottom row). The pial surface voxelisation is designed such that neighbouring gyri can fuse if their separation is smaller than the scale of interest, whilst not growing the cortex outwards. The white matter surface voxelisation allows the walls of each gyrus to thicken and fuse inwards. Visually, this process looks as if the cortex is ‘thickening’ inwards and smoothing on the outside.
From the voxelised cortex at each scale $\lambda$, we can then obtain coarsegrained pial and white matter surfaces, and extract estimates of global morphometric measures such as the average cortical thickness $T(\lambda )$, the total cortical surface area ${A}_{t}(\lambda )$, and the exposed surface area ${A}_{e}(\lambda )$. For each $\lambda$, we derived an outer isosurface equal to 0.5 for the pial voxelisation (voxel values are 1 within the voxelised pial surface and 0 otherwise). This isosurface is then defined as the coarsegrained pial surface at this scale. The surface area of this isosurface is used as ${A}_{t}(\lambda )$. The exposed surface area ${A}_{e}(\lambda )$ is subsequently derived from the convex hull of the pial isosurface. Finally, the cortical thickness $T(\lambda )$ is estimated as $\frac{{V}_{G}(\lambda )}{{A}_{t}(\lambda )}$, where ${V}_{G}(\lambda )$, which is the estimated grey matter volume, derived from the number of grey matter voxels, is multiplied by the voxel volume. These global morphometric measures are summary statistics of the brain at each particular scale, capturing information about both their intrinsic geometry (${A}_{t}(\lambda ),T(\lambda )$) and extrinsic geometry (${A}_{e}(\lambda )$). In this manner, for each cortex, we obtain not just one set of summary morphometric measures, but rather a set of measures as functions of $\lambda$. A detailed walkthrough of the coarsegraining algorithm and estimation of morphometric measures is provided on Github: https://github.com/cnnplab/CorticalFoldingAnalysisTools/blob/master/Scales/fastEstimateScale.m.
Note that although this process is inspired by the boxcounting algorithm, it is different from boxcounting and related convolutionbased algorithms: we do not simply apply successive convolutions with an increasing kernel size (or equivalent), which effectively would achieve a uniform and spatially isometric ‘smearing’ of the original cortical ribbon to a given scale, rather than a targeted erasure of surface details smaller than said scale. As a result, this method yields welldefined and wellbehaved white and gray matter surfaces. Thus, one can apply all the usual analytical tools to these realisations that are applicable to actual cortices. Note that our approach is more comparable in its principles to the calculation of the outer smoothed pial surface in FreeSurfer (Schaer et al., 2008), which utilises a dilation and erosion convolution to effectively erase details below a certain scale. Of course, our proposed procedure is not the only conceivable way to erase morphological details below a given scale; and we are actively working on related algorithms that are also computationally cheaper. Nevertheless, the current version requires no finetuning, is computationally feasible and conceptually simple, thus making it a natural choice for introducing the methodology and approach.
Given how broadly it has been verified, we expect the observed universality in cortical selfsimilar scaling to be robust to the details of the coarsegraining algorithm. It would be very informative to test this proposition in future. For example, an alternative method inspired by Yu et al., 2021 could be implemented, eschewing voxelisation and dealing only with the flow of nested surfaces with self and mutual avoidance explicitly implemented. Going in the opposite direction, more detailed models for the mechanisms of cortical folding (see e.g. Raznahan et al., 2011) can be regarded as a type of reverse melting, and could perhaps be implemented and tested in a similar fashion as finegraining procedures.
Appendix 2
Influence of original image resolution
In this section, we want to demonstrate the relatively weak effect of the original image resolution on our analysis outputs. To this end, we used five example HCP subjects, who were scanned at 0.7mm isotropic image resolution, and downsampled their images to 1mm isotropic images. We then proceeded with our analysis using two freesurfer outputs. (1) the HCP freesurfer pipeline output optimised for the 0.7mm resolution, and (2) a standard freesurfer pipeline output on the 1mm downsampled images. We proceeded with these two sets of surfaces in our analysis and show the resulting morphology measures in Appendix 2—figure 1.
We observe a relatively weak difference between these two sets of inputs/surfaces. Both sets largely follow the same trajectory across scales. Especially in the exposed area, the within and between subject differences are noticeably larger than between image resolutions. In total area and thickness, some small but systematic differences are seen in all subjects between the scales of 1 2mm. These are most likely differences in the cortical morphology reconstruction in the freesurfer surfaces using the two different resolution images as input. However, in none of these resulting morphology measures do we see an artifact specifically at 0.7 or 1mm. Any differences between image resolutions only result in subtle changes in the freesurfer meshes, which are smooth. Our analysis method, therefore, has no direct dependency on the image resolution, as our inputs are these smooth freesurfer meshes.
Appendix 3
Scaling properties by species
3.1 Obtaining the scaling law
From the coarsegrain procedure (Figure 1), we obtain surface meshes for the pial and white matter surface at each spatial scale. From those, we derive exposed area, total pial surface area, and average cortical thickness as described in Methods. However, the coarsegraining procedure, by itself, barely changes the exposed surface area, and only changes the total surface area minimally (Appendix 3—figure 1A). As the voxel size changes at each scale, it only starts to affect the exposed area at very large scales, where effectively the voxels are no longer a good description of the shape of the skull. Thus, if we scatter the raw data points from the coarsegraining procedure in the plane of the scaling law, there is barely any variance in the data (Appendix 3—figure 1B). Even after zooming in, we see an almost vertical line (Appendix 3—figure 1B). This is expected and is not evidence for or against our hypothesis.
Instead, we need to transform the data into a perspective that makes sense to be seen in the scaling law plane. Instead of thinking of the coarsegraining procedure as a process that rerenders the cortex at increasing voxel sizes, we can instead think of the procedure as rescaling the original cortical mesh into increasingly smaller sizes, and rerendering it with the same fixed voxel size. The analogy using Britain’s coastline is that instead of measuring the coastline with increasingly smaller rulers, we resize the map of the coastline to increasingly smaller sizes, but keep the size of the ruler the same. Both procedures are equivalent and produce the same fractal dimension.
To achieve this procedure, we rescaled our measurements of $A}_{e$, $T$, and $A}_{t$ by the voxel size ($\lambda$) and a fixed factor ($l}_{r$):
As we used isometric cubes as voxels, the voxel size refers to the length of a single side. E.g. $\lambda =1$ if we used a isometric $1\times 1\times 1m{m}^{3}$ voxel.
$l}_{r$ is a fixed factor for each cortical hemisphere, and does not change with $\lambda$. We use it to systematically shift all the data points within a range, such that the rescaled quantities are not larger than those from the original cortical meshes. One can easily verify that $l}_{r$ will not change the slope or offset of any scaling law, but simply represents a constant shift to all data points. In our data, some of the rescaled quantities would indeed be larger than those from the original cortical meshes, as the voxel size we choose is limited at the smaller end only by computational resources. In other words, we can use very small voxel sizes (relative to the mesh), which after resizing would yield very large values of $A}_{e}^{r$, $A}_{t}^{r$, and $T}^{r$. To avoid this, we chose $l}_{r$ simply as the ratio of the $I$ (isometric term) of the mesh at the smallest scale we used relative to the original mesh, divided by the $\lambda$ of the smallest scale:
where $\lambda}_{s$ is the smallest voxel size used for a particular cortex, $I}_{s$ is the corresponding $I$ for this cortex at the smallest voxel size.
Finally, $I}_{o$ is the $I$ term for the original cortical mesh. Indeed, the ratio of $\frac{{I}_{s}}{{I}_{o}}$ is always close to, but larger than one in our dataset. Thus, we can see in Figure 2A that most traces start very close to the original data point, indicating that our finest scale is reconstructing the original surfaces well.
Note that the rescaling is isometric, meaning that it only affects $I$, but does not change the data in the $K\times S$ plane (Figure 3). Rescaling by a factor proportional to $\lambda}^{2$ can, therefore, be understood as distributing the data points along the $I$ axis, while $l}_{r$ can be understood as fixing the position in the $I$ axis relative to the $I$ of the original mesh.
Figure 2A is produced with these rescaled quantities as described above. The only final step in producing Figure 2A, and in calculating the associated slopes, is the removal of artifactual data points where $\frac{{A}_{t}}{{A}_{e}}<0$, which can occur at very large voxel sizes relative to the cortical mesh. The algorithmic implementation in MATLAB can be found on Github: https://github.com/cnnplab/CorticalFoldingAnalysisTools/blob/master/Scales/fastEstimateScale.m, as part of our Cortical Folding Analysis MATLAB package https://github.com/cnnplab/CorticalFoldingAnalysisTools/, which also has been recently updated with a graphical user interface.
3.2 Speciesspecific details
In the following, we will show the detailed data for each species in terms of their scaling behaviour in Appendix 3—figure 2. Videos for all species, showing the pial surface at each scale can be found under https://bit.ly/3CDoqZQ for review purposes. Final versions of all underlying data, analysis code, and videos will be published on Zenodo and eLife upon acceptance of the paper.
Given the demonstrated overlap between species, and if their cortices are all approximations of the same form, how can one tell apart their cortices? The answer is that all approximations have a range of validity, which varies between cortices. For a gyrified cortex of area $A}_{t$, (When no dependence on $\lambda$ is indicated then we are referring to the values for the original cortex), the coarsegraining will remove details of everincreasing scale, lowering ${A}_{t}(\lambda )$ until attaining lissencephaly for ${A}_{e}({\lambda}_{lys})={A}_{t}({\lambda}_{lys})$. Indeed, it follows from rewriting Equation 3 in the form of Equation 2 that a given cortex’ shape will be comprised of selfsimilar structures with areas ranging from $A}_{0}={g}^{5}{A}_{t$ at their smallest to $A}_{t$ at their largest, where $g$ is the gyrification index $g=\frac{{A}_{t}}{{A}_{e}}$. Whenever this selfsimilar scaling is valid, $A}_{0$ acquires a further interpretation as the typical size of the smallest structures in a cortex: patches smaller than that must be approximately smooth. Consequently, $N}_{structures}\simeq \frac{{A}_{t}}{{A}_{0}}={g}^{5$ estimates the number of morphological features in each cortex (and adds a new interpretation for the gyrification index). For example, we estimate the human in our dataset to have about 105 morphological features in each cortical hemisphere, and the galago to have about four such morphological features. The corresponding ${N}_{structures}(\lambda )$ estimates how this number changes over coarsegraining. Appendix 4 later provides a more detailed discussion on this topic.
Appendix 4
Morphometric relations and the range of validity of the fractal approximation
There are many geometrical quantities that describe aspects of entire cortical hemispheres. But most information contained in such summary morphometric measures is captured, exactly or to good approximations, by the three used in Equation 1: the pial (or total) surface area $A}_{t$, exposed surface area $A}_{e$, and average cortical thickness $T$. Other quantities, like the gyrification index $g=\frac{{A}_{t}}{{A}_{e}}$, or the Grey Matter and Total volumes, can be expressed as products of power laws of these quantities; and the logarithm of products of power laws are linear combinations of the logarithms of the constituting variables, with their exponents as coefficients. Thus, the morphology of a given cortex can be fairly summarised as a point in a threedimensional morphometric space with components given by $\mathrm{log}{A}_{t}$, $\mathrm{log}{A}_{e}$, and $\mathrm{log}{T}^{2}$ (the last exponent guarantees all axes have the same dimension of $\mathrm{log}[Area]$).
This logarithmic structure guarantees that the quantities denoted by products of power laws of the original variables correspond to vectors in the morphometric space: the increase and decrease of their values occur along these vectors, and the planes perpendicular to each vector denote all configurations with the same value for the associated variable. For example, $G=\mathrm{log}g=\mathrm{log}{A}_{t}\mathrm{log}{A}_{e}$, the point in the morphometric space associated with the gyrification ratio value $\frac{{A}_{t}}{{A}_{e}}$, defines a displacement along a vector with coefficients $\{1,1,0\}$.
The quantities $K$, $S$ and $I$ used in Sec. 3.2 likewise corresponds to displacements along specific directions in morphometric space, which are all orthogonal to one another. These can be normalised, so that the sum of their squared coefficients equals unity. In this case, expressing cortical shape in the morphometric space in terms of the original variables or the new variables amounts simply to a rotation: a change of orthonormal base.
Figures 2 and 3 can be then regarded as particular 2D snapshots, taken from a certain direction, of a 3D set of points.
Of particular note, those vectors with a sum of coefficients equal to zero correspond to dimensionless variables. Any such vector is perpendicular to the $I$ direction, and lie on the $K\times S$ plane.
In this framework, we can picture the stepwise coarsegraining of a cortex, as described in Sec. 3.1, as a trajectory in this morphometric space. Empirically, Sec. 3.1, as codified by Equation 3, finds that the trajectories in morphometric space for all studied primates are linear, and largely overlap along a line of constant $K$. For cortical morphology, the first fact implies scaleinvariance, the second implies universality.
Aging, development, and the progression of neurodegenerative conditions likewise each correspond to different morphometric trajectories.
Now consider a morphometric trajectory, corresponding to the coarsegraining of a cortex up to the point of lissencephaly, i.e., to the spatial scale $\lambda}_{lys$ such that ${A}_{t}({\lambda}_{lys})={A}_{e}({\lambda}_{lys})$. In the idealised case (but very close to empirical, as seen in Sec. 3.1), all points along the trajectory should follow Equation 3, starting at the original values for the morphometric measures and ending at the intersection of the $K=constant$ coarsegraining trajectory and the $g=1$ line that defines the limit for lissencephaly (points with $g<1$ are geometrically possible, but are largely avoided by the trajectories of actual cortices, although not by those of all objects; see Appendix 5—figure 1). This trajectory spans a range of values for each morphometric measure: for dimensionless morphometric measures, such span is fully specified by Equation 3. For all others, one needs also to specify the voxel rescaling. The scheme described in Sec. Appendix 3.1 simply guarantees that, for a constant $k={A}_{t}{A}_{e}^{\frac{5}{4}}{T}^{\frac{1}{2}}$, the fundamental area element $A}_{0}=\frac{{T}^{2}}{{k}^{4}}=\frac{{A}_{e}^{5}}{{A}_{t}^{4}$ (and thus also the average thickness $T$) is kept constant for all realisations of the coarsegraining process. We can regard these as approximations of the same original cortex, all drawn with the same resolution, but in ever smaller isometric sizes.
This choice of spatial scale then allows our identification of $A}_{0$ as the typical size of the smallest morphological features in a given cortex. Equivalently, it is also the area of the largest possible lissencephalic cortex for a given average cortical thickness, just at the cusp of the onset of gyrencephaly.
This same choice also enable us to easily compute the expected span of the various morphometric measures over coarsegraining, assuming the universal scaling: each will range between its original value (expressed as a function of the original values for $A}_{t$, $A}_{e$, and $T$), and its value at spatial scale $\lambda}_{lys$, where $A}_{t}({\lambda}_{lys})={A}_{e}({\lambda}_{lys})={A}_{0$: we simply replace $A}_{t$ and $A}_{e$ by $A}_{0$, and $T$ by $\sqrt{{A}_{0}{k}^{4}}$. The ratios between the initial and final values of all morphometric measures will be given as powers of $g$, imbuing the gyrification index with a new significance.
We list below the expanse for all such morphometric measures ($S=\mathrm{log}s$, $K=\mathrm{log}k$ and $I=\mathrm{log}{v}_{I}$)
Morphometric measure  Original value =  Span ×  Value at lissencephalic limit 

$T}^{2$  $T}^{2$  1  $A}_{0}{k}^{4$ 
$A}_{0$  $\frac{A{e}^{5}}{A{t}^{4}}$  1  $A}_{0$ 
$A}_{e$  $A}_{e$  $g}^{4$  $A}_{0$ 
$A}_{t$  $A}_{t$  $g}^{5$  $A}_{0$ 
$V}_{total$  $\approx \frac{2}{9\sqrt{3\pi}}{A}_{e}^{\frac{3}{2}}$  $g}^{6$  $\approx \frac{2}{9\sqrt{3\pi}}{A}_{0}^{\frac{3}{2}}$ 
$V}_{GM$  ${A}_{t}T$  $g}^{5$  $A}_{0}^{\frac{3}{2}}{k}^{2$ 
$g$  $\frac{At}{Ae}$  $g$  1 
$N}_{features$  $\approx \frac{A{t}^{5}}{A{e}^{5}}$  $g}^{5$  1 
$k$  $A}_{t}{A}_{e}^{\frac{5}{4}}{T}^{\frac{1}{2}$  1  $k$ 
$s$  $A}_{t}^{\frac{3}{2}}{A}_{e}^{\frac{3}{4}}{T}^{\frac{9}{2}$  $g}^{\frac{21}{2}$  $k}^{9$ 
$v}_{I$  $A}_{t}{A}_{e}{T}^{2$  $g}^{9$  $A}_{0}^{3}{k}^{4$ 
$\hat{k}$  $A}_{t}^{\frac{4}{\sqrt{42}}}{A}_{e}^{\frac{5}{\sqrt{42}}}{T}^{\frac{2}{\sqrt{42}}$  1  $\hat{k}$ 
$\hat{s}$  $A}_{t}^{\frac{2}{\sqrt{14}}}{A}_{e}^{\frac{1}{\sqrt{14}}}{T}^{\frac{6}{\sqrt{14}}$  $g}^{\frac{18}{\sqrt{14}}$  $\hat{k}}^{3\sqrt{3}$ 
$\hat{v}}_{I$  $A}_{t}^{\frac{1}{\sqrt{3}}}{A}_{e}^{\frac{1}{\sqrt{3}}}{T}^{\frac{2}{\sqrt{3}}$  $g}^{3\sqrt{3}$  $A}_{0}{\hat{k}}^{\sqrt{14}$ 
Appendix 5
Validation with nonbrain objects
We chose to apply our coarsegraining procedure to a range of objects to validate our algorithm, and we included both negative and positive controls.
As a first positive control, we used a simple box with finite thickness (i.e. we simulated a inner ‘white matter’ surface and an outer ‘grey matter’ surface, both being a cube, one positioned inside the other). A schematic is shown in Appendix 5—figure 1. We know based on theoretical consideration that this box must be aligned exactly with the $K=\frac{1}{9}S$ line, as its $A}_{t$ must be the same as its $A}_{e$ at every stage of coarsegraining, which also implies a ‘fractal dimension’ of 1. This is exactly the case in our plot of $K$ against $S$, as all the grey data points from the coarsegraining algorithm align on $K=\frac{1}{9}S$ (thick black line).
As negative controls, we included three nonbrain objects, a bell pepper, two walnut halves, and a coarse outline of a bust/figurine known as Laurana. All three of these objects show a selfsimilar, or fractal, regime (corresponding to partially straight trajectories in $K\times S$) in our coarsegraining procedure, but their trajectories in $K\times S$ space are not flat (Appendix 5—figure 1), meaning their fractal dimensions are distinct from 2.5. Their trajectories also do not overlap with each other, indicating that they are fundamentally different shapes.
The human brain is included here as a reference and shows a clear flat trajectory (Appendix 5—figure 1) as presented in the main text. Also worth noting, after the human trajectory intersects the $K=\frac{1}{9}S$ line, it veers off to follow the line closely, indicating that the $A}_{e$ and $A}_{t$ remain the same with further steps of coarsegraining. Effectively, the human brain transitions to be lissencephalic convex structures once their fractal regime ends. For simplicity, we excluded these lissencephalic data points of extreme coarsegraining from the results in the main text.
Appendix 6
Validation with randomising grid for coarsegraining
To assess robustness of the coarsegraining algorithm, we ran 30 different realisations of the algorithm, but with a small random shift of the grid position relative to the surface meshes. The shift is chosen to be within the radius of $\lambda$. This allows subtle changes in the voxelisation at the boundary of grey matter and white matter. Interestingly, over five different human individuals, we observed subtle changes in the values $S$ over the 30 realisations, and to a lesser degree $K$ (Appendix 6—figure 1). Additionally, the variation between individuals, especially in $K$ at smaller spatial scales, is far greater than the variation introduced by the jittered realisations. These results suggest that the coarsegraining algorithm is extremely robust towards subtle changes at the grey or white matter boundary, especially far away from the lissencepahlic limit.
Appendix 7
Ageing process
7.1 All morphometric variables
In the main text, we used the ageing process as an example for scaledependent biological process. For completeness, here, we also show all the morphometric variables across scales in the 20 y.o. and 80 y.o. cohort in Appendix 7—figure 1, and the corresponding effect sizes in Appendix 7—figure 2. As expected, $A}_{e$ shows very little difference between the two cohorts. $I$, $S$, and $T$ show some effects at scales smaller than 4 mm, but the effect decreases monotonously for higher scales. $K$ demonstrates a more complex scaledependent effect: larger $K$ in younger subjects at small scales of 0.25 mm (in agreement with previous native scale analyses Wang et al., 2016), smaller $K$ in younger subjects at scales of approx. 2 mm, and a return to larger $K$ in younger subjects at large scales of more than approx. 5 mm. Note, however, that despite these large effect sizes, the actual change of values of $K$ is within the range of variation expected across species (seen in the main text). The range of variation of $K$ is approx 0.04 here, and at least an order of magnitude smaller than the range of variation in $S$ (approx 1.5).
7.2 Confirmation in independent dataset
To confirm that our observed ageing effect was not driven by sample or dataspecific properties, we also analysed an independent dataset (NKI) with the same methods. Here, we show the equivalent figures to Appendix 7—figure 1 and Appendix 7—figure 2 for the NKI data in Appendix 7—figure 3 and Appendix 7—figure 4. The same qualitative patterns can be seen in all plots, including the scalespecific effects in $K$ and $A}_{e$ at around 2 mm.
Data availability
The code for coarsegraining has been integrated into our MATLAB toolbox Cortical Folding Analysis Tools: https://github.com/cnnplab/CorticalFoldingAnalysisTools (copy archived at Wang et al., 2024), which now also includes a graphical user interface. Users will see the latest updates in this repository. The analysis code underpinning this paper is published on GitHub: https://github.com/cnnplab/2024_Folding_scales/ (copy archived at Wang, 2024). The postprocessing data (i.e. 'voxelisations' and derived metrics) are uploaded on Zenodo: https://doi.org/10.5281/zenodo.12820611. The data, together with the code, will allow readers to reproduce of our main results.

ZenodoData underlying "Neuroevolutionary evidence for a universal fractal primate brain shape".https://doi.org/10.5281/zenodo.12820610
References

BookScaling, selfsimilarity, and intermediate asymptotics: dimensional analysis and intermediate asymptotics. cambridge texts in applied mathematicsCambridge University Press.

Diffusion MRI data, sulcal anatomy, and tractography for eight species from the primate brain bankBrain Structure & Function 226:2497–2509.https://doi.org/10.1007/s0042902102268x

Random walk models in biologyJournal of the Royal Society, Interface 5:813–834.https://doi.org/10.1098/rsif.2008.0014

Structural variability across the primate brain: a crossspecies comparisonCerebral Cortex 28:3829–3841.https://doi.org/10.1093/cercor/bhx244

Fractal dimension analysis in neurological disorders: An overviewAdvances in Neurobiology 36:313–328.https://doi.org/10.1007/9783031476068_16

Scaling laws for branching vessels of human cerebral cortexMicrocirculation 16:331–344.https://doi.org/10.1080/10739680802662607

Evolution of neocortical folding: A phylogenetic comparative analysis of MRI from 34 primate speciesCortex; a Journal Devoted to the Study of the Nervous System and Behavior 118:275–291.https://doi.org/10.1016/j.cortex.2019.04.011

Size and shape of the cerebral cortex in mammalsBrain, Behavior and Evolution 27:28–40.https://doi.org/10.1159/000118718

The fractal geometry of convoluted brainsJournal Fur Hirnforschung 32:103–111.

A model of the lungs based on fractal geometrical and structural propertiesIFAC Proceedings Volumes 42:994–999.https://doi.org/10.3182/200907063FR2004.00165

Relationship between sulcal characteristics and brain agingFrontiers in Aging Neuroscience 10:339.https://doi.org/10.3389/fnagi.2018.00339

Is the brain cortex a fractal?NeuroImage 20:1765–1774.https://doi.org/10.1016/s10538119(03)00380x

The Fractal Geometry of the Brain413–429, Fractal analysis of electroencephalographic time series (EEG signals), The Fractal Geometry of the Brain, Springer, 10.1007/9781493939954_25.

Relationship among neuroimaging indices of cerebral health during normal agingHuman Brain Mapping 29:36–45.https://doi.org/10.1002/hbm.20369

BookLocal morphological measures confirm that folding within small partitions of the human cortex follows universal scaling lawIn: Cattin PC, Cotin S, Padoy N, Speidel S, Zheng Y, Essert C, editors. Medical Image Computing and Computer Assisted Intervention – MICCAI 2021. Cham: Springer International Publishing. pp. 691–700.https://doi.org/10.1007/9783030872342_65

BookFractals in biology and medicineJohn Wiley & Sons, Ltd.https://doi.org/10.1002/3527600906.mcb.201100002

Improved understanding of brain morphology through 3D printing: A brief guideResearch Ideas and Outcomes 2:e10398.https://doi.org/10.3897/rio.2.e10398

Testretest reliability of brain morphology estimatesBrain Informatics 4:107–121.https://doi.org/10.1007/s4070801600604

Robust estimation of sulcal morphologyBrain Informatics 6:5.https://doi.org/10.1186/s4070801900981

BookThe fractal geometry of natureNew York: W. H. Freeman and Comp.https://doi.org/10.1119/1.13295

Cortical complexity estimation using fractal dimension: a systematic review of the literature on clinical and nonclinical samplesThe European Journal of Neuroscience 55:1547–1583.https://doi.org/10.1111/ejn.15631

Evolution and development of the mammalian cerebral cortexBrain, Behavior and Evolution 83:126–139.https://doi.org/10.1159/000357753

How does your cortex grow?The Journal of Neuroscience 31:7174–7177.https://doi.org/10.1523/JNEUROSCI.005411.2011

A surfacebased approach to quantify local cortical gyrificationIEEE Transactions on Medical Imaging 27:161–170.https://doi.org/10.1109/TMI.2007.903576

Software2024_Folding_scales, version swh:1:rev:1d60b1fe5675767c61f7984f473771602bdfd008Software Heritage.

SoftwareCorticalFoldingAnalysisTools, version swh:1:rev:a51846299807038b7809158293dd60ca196f54cbSoftware Heritage.

Axons pull on the brain, but tension does not drive cortical foldingJournal of Biomechanical Engineering 132:4001683.https://doi.org/10.1115/1.4001683

Development of cortical folding during evolution and ontogenyTrends in Neurosciences 36:275–284.https://doi.org/10.1016/j.tins.2013.01.006
Article and author information
Author details
Funding
Engineering and Physical Sciences Research Council (EP/L015358/1)
 Yujiang Wang
 Karoline Leiberg
Engineering and Physical Sciences Research Council (EP/Y016009/1)
 Karoline Leiberg
UK Research and Innovation (MR/V026569/1)
 Yujiang Wang
UK Research and Innovation (MR/T04294X/1)
 Peter Neal Taylor
Instituto Serrapilheira (Serra170916981)
 Bruno Mota
Conselho Nacional de Desenvolvimento Científico e Tecnológico (PQ 2017 312837/20178)
 Bruno Mota
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
We thank members of the Computational Neurology, Neuroscience, and Psychiatry Lab (https://www.cnnplab.com) for discussions on the analysis and manuscript, and Dirk Jan Ardesch and Martijn van den Heuvel for helpful discussions and NHP brain surface data. PNT and YW are both supported by UKRI Future Leaders Fellowships (MR/T04294X/1, MR/V026569/1); YW and KL are further supported by the EPSRC (EP/Y016009/1, EP/L015358/1). B Mota is supported by Fundação Serrapilheira Institute (grant Serra1709–16981) and CNPq (PQ 2017 312837/20178).
Version history
 Preprint posted:
 Sent for peer review:
 Reviewed Preprint version 1:
 Reviewed Preprint version 2:
 Reviewed Preprint version 3:
 Version of Record published:
Cite all versions
You can cite all versions using the DOI https://doi.org/10.7554/eLife.92080. This DOI represents all versions, and will always resolve to the latest one.
Copyright
© 2023, Wang 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.
Metrics

 779
 views

 47
 downloads

 1
 citations
Views, downloads and citations are aggregated across all versions of this paper published by eLife.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading

 Evolutionary Biology
Euarthropods are an extremely diverse phylum in the modern, and have been since their origination in the early Palaeozoic. They grow through moulting the exoskeleton (ecdysis) facilitated by breaking along lines of weakness (sutures). Artiopodans, a group that includes trilobites and their nonbiomineralizing relatives, dominated arthropod diversity in benthic communities during the Palaeozoic. Most trilobites – a hyperdiverse group of tens of thousands of species  moult by breaking the exoskeleton along cephalic sutures, a strategy that has contributed to their high diversity during the Palaeozoic. However, the recent description of similar sutures in early diverging nontrilobite artiopodans means that it is unclear whether these sutures evolved deep within Artiopoda, or convergently appeared multiple times within the group. Here, we describe new wellpreserved material of Acanthomeridion, a putative early diverging artiopodan, including hitherto unknown details of its ventral anatomy and appendages revealed through CT scanning, highlighting additional possible homologous features between the ventral plates of this taxon and trilobite free cheeks. We used three coding strategies treating ventral plates as homologous to trilobitefree cheeks, to trilobite cephalic doublure, or independently derived. If ventral plates are considered homologous to free cheeks, Acanthomeridion is recovered sister to trilobites, however, dorsal ecdysial sutures are still recovered at many places within Artiopoda. If ventral plates are considered homologous to doublure or nonhomologous, then Acanthomeridion is not recovered as sister to trilobites, and thus the ventral plates represent a distinct feature to trilobite doublure/free cheeks.

 Evolutionary Biology
 Immunology and Inflammation
The incessant arms race between viruses and hosts has led to numerous evolutionary innovations that shape life’s evolution. During this process, the interactions between viral receptors and viruses have garnered significant interest since viral receptors are cell surface proteins exploited by viruses to initiate infection. Our study sheds light on the arms race between the MDA5 receptor and 5’pppRNA virus in a lower vertebrate fish, Miichthys miiuy. Firstly, the frequent and independent loss events of RIGI in vertebrates prompted us to search for alternative immune substitutes, with homologydependent genetic compensation response (HDGCR) being the main pathway. Our further analysis suggested that MDA5 of M. miiuy and Gallus gallus, the homolog of RIGI, can replace RIGI in recognizing 5’pppRNA virus, which may lead to redundancy of RIGI and loss from the species genome during evolution. Secondly, as an adversarial strategy, 5’pppRNA SCRV can utilize the m^{6}A methylation mechanism to degrade MDA5 and weaken its antiviral immune ability, thus promoting its own replication and immune evasion. In summary, our study provides a snapshot into the interaction and coevolution between vertebrate and virus, offering valuable perspectives on the ecological and evolutionary factors that contribute to the diversity of the immune system.