Hieu

Carlos

Yutong, Yue, Fangping, Oscar Mathilde Karim

Maxime Reda David

Jean-Rassaire Duong

Claire Mathieu

Nathan Guillaume

They

Mines Télécom Atlantique. I would like to express my special

Comparaison du développement du volume cérébral entre nos résultats de simulation, les résultats des simulations du modèle de Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] et les données de la littérature [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF] Axial force F is applied for a one-dimensional hypothetical element of tissue beam. In all cases, two equal and opposite forces are applied in the both ends of the beam, thus the resultant force is zero. (a) Under tension, beam is stretched and elongated. In the case of elastic deformation (middle panel), there exists a tensile stress (σ > 0). In the case of viscoelastic deformation or growing tissue (bottom panel), this stress tends to zeros at longtime limit due to the mechanical relaxation. (b) Under compression, beam is compressed and shortened. In the case of elastic deformation (middle panel), there exists a compress stress (σ < 0). In the case of growing tissue (bottom panel), this stress will be relaxed to zero. (c) Under sufficiently high compression, F > F c where F c is critical force, this beam may suddenly change shape and it shows up curvature. Such a phenomenon is called buckling or mechanical instability. Extracted from [START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF]. . . . . . . . . . . . 15 1.11 Wilhelm His' mechanical analogy: Left. Folding of a rubber tube; Right. Folding of a gut tube during morphogenesis. Adapted from [START_REF] His | On the Principles of Animal Morphology: Letter to Mr[END_REF]. . . . . . . . . . . . . . . 16 1.12 a. A schematic representation of the axonal tension-driven folding hypothesis [START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF] from [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF]. b. Actual distributions of axon tension based on dissection and histology data [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. c. Intracortical differential tangential growth model due to [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF]. d. This differential growth model predicts elastic buckling of the outer layer [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF].

Extracted from [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 1. [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF] Effects of cortical growth rate on wavelength, subcortical growth, and stress in a 2d model of cortical folding. Each column contains spatial maps of a different variable superimposed on the deformed geometry: radial growth G r , tangential growth G t ; radial stress σ r , tangential stress σ t . Each row corresponds to a different scaled ratio of cortical to subcortical growth rates. Extracted from [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF]. . . 21 1.14 Dynamic patterns of tangential growth in human brain development. (a) The regions of highest expansion migrate smoothly from the central sulcus and nearby regions 28-30 weeks (Top) to parietal 30-34 weeks (Middle) and finally to frontal and temporal regions 34-38 weeks (Bottom). (b) In healthy infants, areas of highest cortical expansion (red) are consistent with the trajectory from prenatal development. (c) Schematic illustrating the trajectory of the maximum growth region from primary motor, sensory, and visual cortices (labelled 'pre') to frontal, parietal and temporal lobes (labelled 'post'). Pre, prenatal; post, postnatal. Extracted from [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 1.15 Illustration of spatial discretization for 3 methods: (a) finite difference method (FDM), (b) finite element method (FEM), (c) finite volume method (FVM). . . . 25 1. [START_REF] Ben | Growth and instability in elastic tissues[END_REF] The flowchart of the interaction between MATLAB and FEBio. . . . . . . . . . . 27 1.17 Brain mesh at the initial time (22 weeks of gestation) and after deformation. Each tetrahedron has its own color because of different material properties: the shear modulus is defined according to the distance from the top surface, the tetrahedron close to the surface has a smaller shear modulus. . . . . . . . . . . . . . . . . . . 29 1.18 Diagram of the volume of an irregular tetrahedron after deformation. . . . . . . . [START_REF] Chang | Micropatterned stretching system for the investigation of mechanical tension on neural stem cells behavior[END_REF] 2.1 Brain is modelled as a soft elastic solid and a relative tangential expansion is imposed on the cortical layer as shown at left. The simulation starting from a smooth fetal brain shows gyrification as a result of uniform tangential expansion of the cortical layer. Extracted from [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]

2.13

Comparison of folding processes of the reference ellipsoid based on different global or regional expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . .

2.14

Left: the comparison of the curvatures for the linear and 1 st Gompertz growth regions. Right: the comparison of the 3D GI for the linear and 1 st Gompertz growth regions. "linear_overall" corresponds to the right semi-ellipsoid of 1 st row (in Figure 2.13) which is applied with global linear model, "gompertz_overall" corresponds to the left semi-ellipsoid of 2 nd row which is applied with global 1 st Gompertz model, "linear_half" corresponds to the right semi-ellipsoid of 3 rd row which is applied with regional linear model and "gompertz_half" corresponds to the left semi-ellipsoid of 3 rd row which is applied with regional 1 st Gompertz model. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 Step 1: The computation of the eigenfunctions of the Laplace-Beltrami operator of the cortical mesh provides the decomposition basis for the spectrum analysis.

Step 2: The curvature decomposition based on the eigenfunctions gives a series of coefficients that are the characteristics of curvature in the spectral domain. This step is reversible, which allows for spectral filtering of the curvature, and thus supports the spectral segmentation of cortical folding patterns.

Step 3: The production of the band power spectrum of curvature from the squared coefficients depends on the band design choices detailed in Figure 3.4 a). Extracted from [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF]. . . . . . . . . . . . . . . . . . . . . 3.4 Spectral banding. a) Spectral bands design: frequency intervals based on the hypothesis of doubling frequency. The spatial resolution achieved is illustrated by the relative size of the nodal domains of the last eigenfunction of the basis and the edges of the cortical mesh (median size brain of the database). Log-linear plot. min, median, max: brain of minimum, median and maximum size of the database. b) Spectral sizing: brain size, bandwidth and spatial resolution for the whole database. WL: wavelength, EV: eigenvalue/eigenvector. Extracted from [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF] 3 and cortical growth α t = 1.829t. 72 3.9 Comparison of average of local surface ratio computed with the sphere of a radius 25 mm on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t. . . . . . . . . . . . . 72 3.10 Folding power of B6 for different cortical thicknesses. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t. 73 3.11 (a) Flowchart and tools used to generate a tetrahedral mesh from MRI segmentation; (b) corresponding diagram of tetrahedral mesh generation, the first two set of images are from ITK-SNAP, the others are from ParaView (https: //www.paraview.org/). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 3.12 (a) Mesh for simulation (from Tallinen et al.); (b) Generated mesh in this section. [START_REF] Ioannis S Gousias | Magnetic resonance imaging of the newborn brain: manual segmentation of labelled atlases in term-born and preterm infants[END_REF] 4.1 Regions of highest cortical expansion change over time. Extracted from [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF]. . . . 78 4.2 Flow chart of relationship derivation from the cortical expansion defined by Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] to that defined in the brain growth model of Tallinen et al. [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. . . . 79 4.3 Pipeline of the application of a spatio-temporal cortical expansion map [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] Extracted from [START_REF] Zhao | A robust hole-filling algorithm for triangular mesh[END_REF] For the 1 st subject, three gestational periods are 28-30 wk, 30-33 wk and 33-38 wk; for the 2 nd subject, three gestational periods are 27-31 wk, 31-33 wk and 33-37 wk; for the 3 rd subject, three gestational periods are 27-29 wk, 29-33 wk and 33-36 wk; for the 4 th subject, three gestational periods are 27-30 wk, 30-34 wk and 34-36 wk. For example, LCE(28 -33wk) of the 1 st subject is calculated by LCE(28 -30wk) × LCE (30 -33wk)

Résumé étendu

La recherche dans cette thèse a été soutenue par le projet "Modélisation computationnelle biophysique pour l'étude de l'impact de l'AVC néonatal sur la croissance cérébrale", subvention DIC20161236453 de la Fondation pour la Recherche Médicale. Avec une prévalence de 1/2000 à 1/4000 naissances vivantes, l'accident vasculaire cérébral (AVC) ischémique périnatal est la forme la plus fréquente d'AVC infantile et constitue la principale cause de paralysie cérébrale unilatérale ou d'épilepsie chez les enfants nés à terme. [START_REF] Kirton | Life after perinatal stroke[END_REF]. Cependant, les causes de ces incapacités observées restent floues. Bien que les techniques basées sur l'image puissent fournir de nouvelles perspectives pour le diagnostic et la prédiction, elles ne peuvent pas fournir des informations sur les causes des incapacités sans informations biophysiques. Des revues récentes décrivent comment la modélisation informatique de la croissance cérébrale commence à relier les échelles (du niveau cellulaire vers le niveau cérébral et le fonctionnement au niveau cérébral) pour comprendre la croissance cérébrale et faire des prédictions efficaces [START_REF] Budday | Physical biology of human brain development[END_REF][START_REF] Goriely | Mechanics of the brain: perspectives, challenges, and opportunities[END_REF]. Les modèles basés sur la physique sont la pierre angulaire de ces études, en quantifiant par exemple les contraintes corticales et en prédisant les indices de gyrification. Comprendre le développement cortical (Figure 1) grâce à la modélisation biomécanique fait partie de nouveaux domaines de recherche prometteurs visant à comprendre la structure et la fonction du cerveau grâce à des approches basées sur la physique [START_REF] Goriely | Mechanics of the brain: perspectives, challenges, and opportunities[END_REF]. En particulier, un travail récent a montré que des modèles biomécaniques pourraient être utilisés en conjonction avec des observations IRM 3D in vivo pour simuler la croissance cérébrale précoce [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Les résultats de la simulation numérique ont démontré que l'expansion tangentielle relative du cortex cérébral contraint par la substance blanche a généré un stress de compression, résultant en des sillons et des gyri similaires à ceux du développement du cerveau foetal [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Le mécanisme d'expansion tangentielle relative et une simulation à partir d'un cerveau foetal de 22 semaines sont illustrés dans la Figure 2.

Malgré ces progrès récents, de nombreuses questions restent ouvertes concernant la morphogenèse des modèles de plissement, y compris les liens entre les paramètres physiques des modèles de simulation et les motifs de plissement observés dans les données IRM in vivo. Ainsi, cette thèse vise d'abord à étudier le modèle biomécanique de croissance du cerveau humain qui est basé sur l'hypothèse de croissance tangentielle différentielle [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], puis à étudier l'influence xv des paramètres physiques dans ce modèle, tels que la géométrie initiale, la croissance corticale et l'épaisseur corticale initiale, sur la morphologie de la surface.

De plus, nous constatons qu'il y a un manque de corrélations entre les résultats de simulation du cerveau foetal et les faits biologiques. En d'autres termes, la correspondance concernant l'augmentation du volume cérébral, le développement des motifs de plissement cortical et l'âge gestationnel n'est pas réaliste dans le modèle de croissance du cerveau humain. Pour résoudre ce problème, le deuxième objectif est d'améliorer le modèle de croissance du cerveau humain en assurant la précision du processus précoce d'expansion et de plissement du cerveau simulé par le modèle. L'expansion tangentielle corticale définie dans le modèle de croissance du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] est considérée comme un paramètre biophysique qui change linéairement avec le temps mais est spatialement invariant dans le cortex. Cependant, dans des études récentes [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Habas | Early folding patterns and asymmetries of the normal human brain detected from in utero mri[END_REF][START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF], les auteurs ont observé des différences de croissance régionales significatives à travers les surfaces corticales foetales et des changements de modèles de croissance pendant la période de gestation. Par conséquent, le troisième objectif est d'intégrer le mécanisme d'expansion corticale différentielle spatio-temporelle [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] dans le modèle biomécanique de croissance cérébrale [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], et d'appliquer une véritable carte d'expansion corticale au cerveau foetal pour tenter de connaître l'effet de l'expansion corticale spatio-temporelle sur les plis locaux.

Cette thèse est organisée de la façon suivante.

Le chapitre 1 présente le développement du cerveau humain qui comprend la structure du cerveau humain, la croissance, les modèles de croissance du cerveau, l'expansion corticale et le plissement cortical. Les concepts biomécaniques et la méthode des éléments finis (FEM) sont également brièvement présentés dans cette section pour introduire davantage les modèles biomécaniques et la modélisation du cerveau humain.

Le chapitre 2 présente le modèle biomécanique de plissement du cerveau humain basé sur l'hypothèse de la croissance tangentielle différentielle [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] en détail, ainsi que les études des effets de la croissance corticale, l'épaisseur corticale initiale et la géométrie initiale sur les motifs de plissement cortical. De plus, ce chapitre présente également diverses métriques pour quantifier la morphologie de la surface corticale, telles que la courbure moyenne, l'indice de xvi gyrification tridimensionnel basé sur la surface, la profondeur du sulcal. De plus, une nouvelle approche pour calculer l'orientation de plis et son anisotropie est proposée.

Le chapitre 3 propose un nouveau modèle de croissance de la longueur longitudinale du cerveau (BLL) et un modèle d'augmentation de l'épaisseur corticale qui sont intégrés dans le modèle de plissement du cerveau humain pour améliorer la précision de la croissance du cerveau pendant le processus de plissement cérébral. Sur la base du modèle combiné, l'impact de l'épaisseur corticale initiale du cerveau foetal sur les motifs de plissement est étudié à travers une analyse visuelle et des caractéristiques quantitatives. À la fin de ce chapitre, une enquête sur la génération d'un maillage tétraédrique de haute qualité à partir de segmentations IRM de cerveaux humains est démontrée.

Le chapitre 4 présente le processus de dérivation de la relation entre deux définitions d'expansion corticale différentes et un pipeline pour appliquer une véritable carte d'expansion corticale différentielle spatio-temporelle à un cerveau foetal de 22 semaines. Ensuite, ce chapitre montre une comparaison des résultats entre l'expansion corticale différentielle spatio-temporelle et l'expansion corticale linéaire spatialement invariante et une comparaison de la complexité de plissement entre différents lobes de l'expansion corticale différentielle spatio-temporelle. De plus, il existe une comparaison visuelle de l'orientation des plis primaires entre la surface simulée et une surface d'atlas cérébral foetal.

Le chapitre 5 conclut cette thèse et introduit des perspectives de travaux futurs, tels que l'exploration des méthodes de lissage géométrique et de correction automatique du maillage, l'investigation de l'effet d'autres changements géométriques (sauf le rapport d'allongement et la fissure longitudinale) et de variations spatiales d'épaisseur corticale sur les motifs de plissement. En plus, pour les différences entre les résultats de la simulation et les données réelles, l'idée sera d'utiliser une approche basée sur les données pour apprendre la dynamique de l'évolution corticale et obtenir le tenseur du gradient de déformation de croissance dans le modèle à partir des données de croissance de cerveaux humains, afin d'améliorer le modèle biophysique.

Le but ultime de cette thèse est de comprendre le développement précoce du cerveau humain et d'explorer les causes de plusieurs motifs de plissement cortical anormaux en utilisant la modélisation biomécanique. Les principales contributions de cette thèse sont :

• Une étude de l'influence de la croissance corticale, de la géométrie initiale et de l'épaisseur corticale initiale sur la morphologie de la surface basée sur un modèle biomécanique de croissance du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

• Une approche pour calculer les angles de plis sur la surface corticale et mesurer l'anisotropie de l'orientation de plis.

• Modélisation de la croissance du volume du cerveau humain pour permettre au modèle de plissement du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] de simuler le processus de croissance du volume du cerveau réaliste pendant le processus de plissement.

• Modélisation et intégration du mécanisme d'expansion corticale différentielle spatio-temporelle dans le modèle de plissement du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

Plus en détail, la première contribution de cette thèse, présentée au chapitre 2, consiste en une exploration du comportement d'un modèle biomécanique de croissance du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] en étudiant l'impact de la densité du maillage et les paramètres biophysiques (l'expansion corticale, la géométrie initiale et épaisseur corticale) sur la morphologie de la surface à l'aide d'ellipsoïdes. Nous montrons tout d'abord que lorsque la densité du maillage atteint un certain ordre de grandeur (10 6 tétraèdres/cm 3 ), les surfaces simulées peuvent atteindre une précision xvii de plissement suffisante, ce qui est utile pour les simulations ultérieures basées sur la croissance du cerveau humain à l'aide de ce modèle.

De plus, nous démontrons que le mode de croissance corticale affecte peu le degré de complexité de la morphologie de surface ; la variation de la géométrie initiale modifie l'orientation et la profondeur des plis, et en particulier, plus la forme est élancée, plus il y a de plis le long de son axe le plus long et plus les sillons deviennent profonds. De plus, plus l'épaisseur corticale initiale est fine, plus la fréquence spatiale des plis est élevée, mais plus les sillons deviennent peu profonds, ce qui est en accord avec les effets précédemment rapportés de l'épaisseur corticale [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Wang | On early brain folding patterns using biomechanical growth modeling[END_REF]. Ces résultats tendent à montrer que l'utilisation de tels modèles biomécaniques pourrait mettre en évidence les liens entre les troubles neurodéveloppementaux et les paramètres physiques.

La deuxième contribution de cette thèse également introduite au chapitre 2, est une nouvelle approche pour calculer les angles des plis à l'aide d'outils géométriques [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF][START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF] et mesurer l'anisotropie de l'orientation des plis en utilisant le Kullback-Leibler divergence. L'angle du pli est défini comme l'angle entre le gradient du vecteur de Fiedler [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF] et la direction principale de la courbure [START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF].

La troisième contribution est orientée vers l'amélioration de la précision de la croissance volumique du cerveau humain pendant le processus de déformation, qui est présentée dans le chapitre 3. Nous proposons un modèle de croissance de la longueur longitudinale du cerveau (BLL) dérivé de la mesure normative du cerveau foetal [START_REF] Kyriakopoulou | Normative biometry of the fetal brain using magnetic resonance imaging[END_REF] et l'introduisons dans le modèle biomécanique [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] pour former un modèle combiné, comme le montre la Figure 3. Ensuite, nous validons ce modèle BLL en comparant le volume cérébral simulé avec d'autres données de volume cérébral validées de la littérature [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF] dans la Figure 4. Sur la base du modèle combiné de croissance du cerveau humain présenté dans le chapitre 3, une étude de l'effet de l'épaisseur corticale initiale sur les motifs de plissement cortical cérébral humain est effectuée. Nous constatons que l'épaisseur corticale plus mince peut conduire à un pouvoir de plissement plus élevé correspondant aux plis cérébraux tertiaires. De plus, nous montrons que la morphologie avec l'épaisseur corticale initiale de 5,96 mm est similaire à la lissencéphalie et celle avec l'épaisseur corticale initiale de 0,74 mm ressemble à la polymicrogyrie. Cette observation peut être importante pour comprendre les causes de plusieurs troubles neurodéveloppementaux associés à des motifs de plissement cortical anormaux. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] et les données de la littérature [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF].

Étant donné que le modèle biomécanique de croissance du cerveau humain [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] a besoin d'un maillage tétraédrique de haute qualité comme point initial pour produire les motifs de plissement précis du cerveau et que nous n'avons généralement que des données de segmentation du cerveau humain, nous proposons donc également un pipeline pour générer un maillage tétraédrique à partir de la segmentation du cerveau humain dans le chapitre 3, comme le montre la Figure 5. Cependant, le pipeline doit être optimisé en termes de précision géométrique et d'efficacité de la correction de maillage. Les résultats dans cette thèse tendent à montrer que l'utilisation de tels modèles biomécaniques pourrait mettre en évidence les liens entre les maladies neurodéveloppementales et les paramètres biophysiques.

xx

List of Symbols

The following list describes several symbols that will be later used throughout this dissertation.

x (or X) Spatial coordinate t Time

θ angle u displacement L length S Surface V Volume M Mass λ Wavelength v L Longitudinal

General Introduction

Motivation

With a prevalence of 1/2000 to 1/4000 live births, perinatal ischemic stroke is the most frequent form of childhood stroke and constitutes the leading cause of unilateral cerebral palsy in term-born children [START_REF] Kirton | Life after perinatal stroke[END_REF]. Perinatal ischemic stroke is an umbrella term including several conditions that differ in pathophysiology, timing and thus in outcomes. Neonatal arterial ischemic stroke (NAIS) that refers to a perinatal ischemic stroke syndrome with neonatal signs related to an arterial infarct as revealed by brain imaging can lead to cerebral palsy or epilepsy.

Every case of NAIS is unique to the individual. Considering for instance NAIS leading to unilateral cerebral palsy, one person may have total paralysis and requires constant care, while another with partial paralysis might have slight movement tremors but requires little assistance. This is due to the type of injury and the timing of the injury to the developing brain. The prediction of long-term motor outcome requires new personalized approaches, i.e., patientspecific techniques to understand the causes of the observed disabilities. Although image-based techniques can provide new insights for diagnosis and prediction, it cannot provide insights on causes of observed disabilities without biophysical information.

Recent reviews describe how computational modeling of brain growth starts to bridge the scales (from the cellular level toward from and function at the brain level) to understand the brain growth and make efficient predictions [START_REF] Budday | Physical biology of human brain development[END_REF][START_REF] Goriely | Mechanics of the brain: perspectives, challenges, and opportunities[END_REF]. Understanding cortical development through biomechanical modeling is part of new promising research areas aiming at understanding the structure and function of the brain through physics-based approaches [START_REF] Goriely | Mechanics of the brain: perspectives, challenges, and opportunities[END_REF]. In particular, a recent work has shown that biomechanical models could be used in conjunction with 3D in vivo MRI observations to simulate the early brain growth [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. The numerical simulation results has demonstrated that the relative tangential expansion of the cerebral cortex constrained by the white matter generated compressive stress, resulting in cusped sulci and smooth gyri similar to those in developing fetal brains [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

Despite these recent advances, many questions remain open regarding the morphogenesis of folding patterns, including links between the physical parameters of simulation models and the folding patterns observed in in vivo MRI data. Thus this thesis first aims at studying the biomechanical human brain growth model which is based on differential tangential growth hypothesis [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], and then investigating the influence of physical parameters in this model, such as the initial geometry, the cortical growth and the initial cortical thickness, on surface morphology.

In addition, we find that there is a lack of correlation between simulation results of fetal brain and biological facts. In other words, the correspondence concerning the increase in brain volume, the development of cortical folding patterns and the gestational age is not realistic in the human brain growth model. To solve this problem, the second purpose is to improve the human brain growth model by ensuring the accuracy of the early brain development process simulated by the model.

The cortical tangential expansion defined in the human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is regarded as a biophysical parameter that changes linearly with time but is spatially invariant on the cortical surface. However, in recent studies [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Habas | Early folding patterns and asymmetries of the normal human brain detected from in utero mri[END_REF][START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF], the authors observed significant regional growth differences across fetal cortical surfaces and growth patterns changes during the gestation period. Therefore, the third aim is to integrate the spatio-temporal differential cortical expansion mechanism [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] into the biomechanical human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], and apply a cortical expansion map to a fetal brain in an attempt to know the effect of the spatio-temporal cortical expansion on local folding patterns.

Thesis overview

The ultimate goal of this thesis is to understand the early development of the human brain and explore the causes of several abnormal cortical folding patterns using biomechanical modeling. The main contributions of this thesis are:

• An investigate of the influence of the cortical growth, the initial geometry and the initial cortical thickness on surface morphology based on a biomechanical human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

• An approach for calculating the fold angles on cortical surface and measuring the anisotropy of the folding orientation.

• Modeling of the human brain volume growth to allow the human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] to simulate the realistic brain volume growth process during the folding process.

• Modeling and integration of the spatio-temporal differential cortical expansion mechanism [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] in the human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

Thesis organisation

This thesis is organised as follows: Chapter 1 introduces the human brain development which includes the human brain structure, growth, normative brain growth models, cortical expansion and folding, and quantitative description of cortical folding. The biomechanical concepts and finite element method (FEM) are also briefly presented in this section to further introduce the biomechanical models and modeling of the human brain.

Chapter 2 presents the biomechanical human brain folding model based on the hypothesis of differential tangential growth [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], as well as the investigations of the effects of the cortical growth, the initial cortical thickness and the initial geometry on folding patterns. In addition, this chapter also introduces various metrics for quantifying the cortical surface morphology, such as the mean curvature, the surface-based three-dimensional gyrification index, the sulcal depth. Moreover, a novel approach for calculating the folding orientation and its anisotropy is proposed.

Chapter 3 proposes a new brain longitudinal length (BLL) growth model and a cortical thickness increase model which are integrated into the human brain folding model to improve the brain growth accuracy during the cerebral folding process. Based on the combined model, the impact of the initial cortical thickness of the fetal brain onto cortical folding patterns is studied through visual analysis and quantitative features. At the end of this chapter, an investigation of the generation of a high-quality tetrahedral mesh from MRI segmentation of the human brain is demonstrated.

Chapter 4 presents the derivation process of the relationship between two different cortical expansion definitions and a pipeline for applying a real spatio-temporal differential cortical expansion map to a 22 weeks' fetal brain. Then this chapter shows a comparison of the simulation results between the spatio-temporal differential cortical expansion and spatially invariant linear cortical expansion, and a comparison of the folding complexity between different lobes of the spatio-temporal differential cortical expansion. Furthermore, there is a visual comparison of the orientation of the primary folds between the simulated surface and a fetal brain atlas surface.

Chapter 5 concludes this thesis and introduces perspectives for future works, such as the exploration of the methods of geometric smoothing and automatic mesh correction, the investigation of the effect of other geometric changes (except the elongation ratio and the longitudinal fissure) and areal variations in cortical thickness on folding patterns. One of the most significant differences between classical engineering materials and living creatures is that the latter can grow and remodel. Thus the growth is a critical development process for all organisms, such as the human brain. Even though the development of the human brain is intricate, and many things remain unknown, the anatomical structure and cortical morphogenesis have been well documented from the microscopic to macroscopic level [START_REF]Central nervous system[END_REF][START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF][START_REF] Lewis D Griffin | The intrinsic geometry of the cerebral cortex[END_REF][START_REF] Lemaitre | Normal age-related brain morphometric changes: nonuniformity across cortical thickness, surface area and gray matter volume?[END_REF][START_REF] Lewitus | Conical expansion of the outer subventricular zone and the role of neocortical folding in evolution and development[END_REF][START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF]. In addition, the leading biomechanical hypotheses and corresponding mechanisms are proposed and developed by many research teams. These mechanisms underlying the developmental process of cerebral cortical folding.

In this chapter, the key elements of the human brain development are first introduced, involving human brain structures (in Section 1.1.1), growth (in Section 1.1.2), normative brain growth models (in Section 1.1.3), cortical expansion and folding (in Section 1.1.4), and quantitative description of cortical folding (in Section 1.1.5). Moreover, important biomechanical concepts are briefly presented in Section 1.2.1 to better understand the existing biomechanical models of the human brain shown in Section 1.2.2. Finally, finite element method (FEM) and an attempt at biomechanical modeling using open source FEM codes are introduced in Section 1.3.

Human brain development 1.Human brain structure

The human brain consists of the brain stem, cerebellum and cerebrum, as shown in Figure 1.1. The cerebrum is divided into left and right halves by a deep groove called the longitudinal fissure [START_REF] Marian | The human brain coloring book[END_REF]. The cerebral cortex, which is also called cortical grey matter, is a thin folded sheet of neural tissue that forms the outermost layer of the cerebrum. It is a structure of six-layers which are organized from outermost pial surface to innermost white matter surface, and surrounds white matter [START_REF] Meyer | Genetic control of neuronal migrations in human cortical development[END_REF][START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF]. The layers from 1 st to 3 rd are transmitters and targets for corticocortical connexions; the 1 st layer is mainly made of fibers and contains only few interneurons; the 2 nd layer consists of many interneurons and few pyramidal neurons; the 3 rd layer is the main source of cortico-cortical efferents, made of interneurons and pyramidal neurons. The 4 th layer is the main target of afferent fibers in the thalamus and the target of intra-hemispheric cortico-cortical afferents; consists of interneurons and pyramidal neurons. The 5 th layer consists of large pyramidal neurons that send efferent axons to subcortical area. The 6 th layer contains few large pyramidal neurons and many small pyramidal neurons sending efferent fibers to the thalamus. The thickness of the human cerebral cortex is approximately 1-4.5 mm and the total average thickness is about 2.5 mm [START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF][START_REF] Lewis D Griffin | The intrinsic geometry of the cerebral cortex[END_REF]. Like other brain folded species, gyral regions are thicker than sulcal regions: the average thickness of gyral crown is 2.7 mm thick, while the average thickness of sulcal fundi is 2.2 mm [START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF]. The cortical surface area of the adult human brain is approximately 1600 cm 2 relative to brain volume of 1400 cm 3 [START_REF] Lemaitre | Normal age-related brain morphometric changes: nonuniformity across cortical thickness, surface area and gray matter volume?[END_REF]. The highly convoluted shape of the cortex allows the larger surface area to be fitted within the relatively small size of the cortical volume [START_REF] Lemaitre | Normal age-related brain morphometric changes: nonuniformity across cortical thickness, surface area and gray matter volume?[END_REF][START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF].

At the convoluted cortical surface, there are gyri (outward hills) and sulci (inward valleys) (Figure 1.2), which form during gestational weeks 16-40 [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF]. The positions of primary gyri and sulci are consistent among individuals, while the secondary gyri and sulci, which appear in the later stages of cortical development, are highly variable in location and appearance [START_REF] Lewis D Griffin | The intrinsic geometry of the cerebral cortex[END_REF]. The primary gyri and sulci could divide the two hemispheres into six major lobes: frontal lobe, parietal lobe, occipital lobe, temporal lobe, limbic lobe and insular lobe (Figure 1.3) [START_REF] Marian | The human brain coloring book[END_REF]. Each lobe has dedicated functions.

Growth

The growth problems are critical to all aspects of life, including cell division, morphogenesis, development, maintenance, cancer, and aging. All forms of life have experienced a certain degree of growth, and biological cells, tissues, organs and organisms exhibit a remarkable ability to grow throughout their entire lives [START_REF] Ambrosi | Growth and remodelling of living tissues: perspectives, challenges and opportunities[END_REF]. The classification of the growth is obtained by considering that the object is changed by its mass, or its material properties, or rearranging the relative positions of the material points [START_REF] Goriely | The mathematics and mechanics of biological growth[END_REF]. These three main processes are mass change, remodeling (property change) and morphogenesis (shape change).

The change in mass relates to a volume change or a change in density, i.e., it can occur either through a volume increase at constant density, as found during soft tissues development, or through a change in density at constant volume, as in the case of bone densification, or both found in a developing bone. Mathematically, growth theory must consider changes in mass, volume, and density, and must flexibly account for that mass penetrates the boundaries of the body, accumulates at the boundary, or occurs inside the body itself.

The term remodeling refers to the process of changing the characteristics of materials in a system without changing mass. It is well known that tissues can harden or soften during aging. In tissues, when old constituents are replaced by new ones through degradation and deposition, that have different orientations, diameters, or other characteristics (such as stiffness, etc.), changes in the internal microstructure can arise [START_REF] Ambrosi | Growth and remodelling of living tissues: perspectives, challenges and opportunities[END_REF]. A change in the microstructure can determine the overall behavior of the tissue [START_REF] Goriely | The mathematics and mechanics of biological growth[END_REF]. For example, a typical composition of soft tissue in many animals is a mixture of collagen fibers in an elastin matrix. Although elastin content remains largely unchanged for many years, there is a continuous conversion of collagen that depends on local biochemical and mechanical stimuli acting on the cells. The relative content of different types of collagen fibers and elastin determines the overall response of the tissue [START_REF] Jay | Cardiovascular solid mechanics: cells, tissues, and organs[END_REF]. This process can occur without changing the quality, but it is important to understand the tissue's response to mechanical loads. From a mathematical point of view, changes in material properties can either be modeled by considering a separate evolution of material parameters in a system, or by taking into account the evolution of separate tissue components at a lower scale.

Morphogenesis is a biological process that causes an organism to develop its shape. The process controls the organized spatial distribution of cells in the embryonic development of an organism, new tissues and organs are formed during this period. In this process, major reorganization and differentiation of cells occur after cell division, and it is important that there is a reorganization of material elements. This restructuring process takes place only when the adhesion between different constituents is weak enough so that they can separate and reattach. The observation of the restructuring process is important for modeling because the morphological tissues or organs undergoing morphogenesis exhibit rapid elastic stress relaxation and plasticlike flow [START_REF] Goriely | The mathematics and mechanics of biological growth[END_REF], which are usually described by modeling tissues as fluid, elastic or viscoelastic in mathematics [START_REF] Morin | Biomechanical modeling of brain soft tissues for medical applications[END_REF][START_REF] Nam | Strain-enhanced stress relaxation impacts nonlinear elasticity in collagen gels[END_REF][START_REF] Özkaya | Introduction to deformable body mechanics[END_REF][START_REF] Thibault P Prevost | Biomechanics of brain tissue[END_REF].

From a biologist's perspective, growth is mediated by gene activation and regulation. Simply put, different genes trigger different growth responses. For instance, the growth of the human brain cortex relies on the expansion of neural stem cells (NSCs), neural progenitors (NPs), and the subsequent generation of postmitotic neurons [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF]. The genetic understanding of growth is fundamental at the smallest scale, and it provides information at the local level on the changes in shape and volume of the growing components of a body. However, in many cases, volume change, remodeling and morphogenesis processes are largely dependent on mechanical factors and the related mechanobiological responses at the cellular level [START_REF] Ambrosi | Growth and remodelling of living tissues: perspectives, challenges and opportunities[END_REF]. Thus it's important to consider the physical, geometric and mechanical constraints in the development of soft tissues or organs and integrate them with genetic and biochemical signals to get a complete growth map.

Normative brain growth models

Considering that the computational modeling of human brain growth needs a proper model to describe the brain growth process, thus we also study the existing normative brain growth models. Several models that have been successfully applied to biological growth data can be used to represent normal brain growth, such as Power Law, Gompertz, Weibull and West [START_REF] Peterson | Normative human brain volume growth[END_REF].

Two-term power law

Power law is one in which two quantities change in proportional amounts, regardless of the initial size of these quantities, one quantity varies as a power of another. For the normative brain growth after delivery, considering the brain volume at birth which is not zero because of brain growth during gestation, the power law must be at least two-term. The function of the two-term power law for normative brain growth is:

α(t) = at b + c, (1.1)
where t is the time, a, b and c are the best-fit parameters.

Gompertz

The Gompertz distribution was initially described to fit survivals and played an important role in the early development of the insurance industry [START_REF] Bajzer | New dimensions in gompertzian growth[END_REF]. The cumulative distribution can be parameterized to have an S-shaped curve and fit animal growth. The Gompertz distribution can also simulate the slowing of tumor growth with necrosis [START_REF] Kristin R Swanson | Virtual and real brain tumors: using mathematical modeling to quantify glioma growth and invasion[END_REF] and track the growth of myelinated white matter in the brain of young children [START_REF] Sadeghi | Regional characterization of longitudinal dt-mri to study white matter maturation of the early developing brain[END_REF]. The model of Gompertz is:

α(t) = ae -e -b(t-c ) , (1.2)
where t is the time, a, b and c are the best-fit parameters.

Weibull

The Weibull function's cumulative distribution can also be parameterized to generate an Sshaped curve. Weibull originally designed this function to model chain failure [START_REF] Weihull | A statistical distribution function of wide applicability[END_REF], which occurs when any member of the chain links fail. Therefore, it is used in hazard and survival analysis, but it can also model growth, because any member's growth or division will lead to overall growth. The model of Weibull is:

α(t) = α m -ae -bt c , (1.3)
where t is the time, a, b and c are the best-fit parameters, α m is the maximum growth rate.

West

The West model of ontogenic growth was presented as a fundamental model for the growth of organisms and organs based on cell division and the energy requirements for cellular maintenance and division [START_REF] West | A general model for ontogenetic growth[END_REF]. It has been shown to fit the growth of a wide variety of animal species as well as the human brain [START_REF] Peterson | Normative human brain volume growth[END_REF]. The model of West ontogenic is:

( V V m ) 1/4 = 1 -(1 -[ V i V m ] 1 4
)e -at/4V 1/4 m , (1.4) where t is the time, a is the best-fit parameter, V is the volume, V m is the maximum volume and V i is the initial volume.

Peterson et al. identified the Weibull model as the most effective growth curve fit for both males and females, while the two-term power law model generated the worst scores among the four models [START_REF] Peterson | Normative human brain volume growth[END_REF].

Cortical expansion and folding

The important steps in the evolution of the human cerebral cortex are the expansion and folding of cortical surface [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF]. The increased growth and folding of the human cerebral cortex are important factors, which are related to the higher-order cognitive abilities [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF].

Cortical expansion

Cortical expansion takes place by the proliferation and expansion of neural progenitors and neural stem cells in the ventricular and subventricular regions during early brain development [START_REF] Chenn | Regulation of cerebral cortical size by control of cell cycle exit in neural precursors[END_REF][START_REF] Fernández | Cerebral cortex expansion and folding: what have we learned?[END_REF][START_REF] Lui | Development and evolution of the human neocortex[END_REF][START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF]. The significant increases in size of the developing brain occur during the third trimester of pregnancy. Many studies have shown that different regions of the cortex have different growth rates [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF]. Rajagopalan et al. [START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF] found that the local growth in the fetal cortical plate was either similar to or significantly greater than the overall cerebral rate through deformation tensor analysis. The deformation tensor model was used in the cortical plate because it is sensitive to focal changes in both cortical thickness and area.

Wright et al. [START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF] observed that regional differences in growth rate of fetal brains by parceling each fetal cortex using a nine regions anatomical atlas [START_REF] Ioannis S Gousias | Magnetic resonance imaging of the newborn brain: manual segmentation of labelled atlases in term-born and preterm infants[END_REF] (Figure 1.4), and they found that Gompertz model can fit the growth change in each lobar region. Specifically, the parietal and posterior temporal lobes showed the fastest growth, while the cingulate, frontal and medial temporal lobes developed more slowly. Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] provided specific spatio-temporal cortical growth maps for individual brains. They observed significant regional differences (spatial variations) in growth across the cortical surfaces of 30 preterm infants. Moreover, they found that these growth patterns change during the rapid cortical expansion period of 28-38 gestational ages (temporal variations) with the noninjured subjects following a highly consistent trajectory. This observation provides a detailed picture of the dynamic changes in cortical growth, linking knowledge about patterns of development at the microscopic (cellular) and macroscopic (folded) scales.

Cortical folding

Cortical folding is a hallmark of many but not all mammalian brains. The degree of folding increases as the size of the mammalian brain increases, but the scales are different between families. The cortical folding of the developing human brain takes place during gestational weeks 16-40 [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF], and especially during the third trimester of pregnancy, the brain undergoes dramatic cortical folding process [START_REF] Je G Chi | Gyral development of the human brain[END_REF][START_REF] Ma Rutherford | Magnetic resonance imaging of the brain in preterm infants 24 weeks' gestation to term[END_REF]. What causes the brain to fold has been the topic of much research. The main focus of interest is to reveal the underlying process by which the relatively consistent primary folds occur. Although the precise mechanism that causes coincident primary sulci but highly variable small sulci across individuals has not been fully understood and demonstrated, it has been shown that normal cortical growth and the resulting folding patterns are important for normal brain function [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF][START_REF] Sandra F Witelson | The exceptional brain of albert einstein[END_REF]. Defects in neurodevelopment (such as neuronal proliferation, migration and differentiation) cause disruption of cortical folding, which may occur at the early stage of the brain development and is related to a series of cognitive deficits in many genetic brain malformations and developmental disorders [START_REF] James Barkovich | A developmental and genetic classification for malformations of cortical development: update 2012[END_REF][START_REF] Fernández | Cerebral cortex expansion and folding: what have we learned?[END_REF].

In order to better understand the mechanisms of normal and abnormal cortical folding, it is important to observe early cortical folding patterns in the developing brain. The human cerebral cortex exhibits dramatic regional expansion and folding in fetal life, with the most prominent and dynamic genetic regulation [START_REF] Colantuoni | Temporal dynamics and genetic control of transcription in the human prefrontal cortex[END_REF][START_REF] Miller | Transcriptional landscape of the prenatal human brain[END_REF]. The primary cortical sulcal folding patterns are determined prenatally and are under strong spatio-temporal genetic control [START_REF] Hill | Similar patterns of cortical expansion during human development and evolution[END_REF][START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF]. The appearance of neocortical gyri and sulci can be divided into two stages: (1) the demarcation of primary gyri at GW (gestational weeks) 23-31; (2) the emergence of the secondary gyri and the growth of sulci in length and depth between the late stages of fetal development and early stages of postnatal life [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Je G Chi | Gyral development of the human brain[END_REF] (Figure 1.5). The formation of the primary cortical folds occurs at the early stage of the cerebral cortical growth, and their formations may be associated with the functional specialization of the cortex and the the map of cytoarchitectonic areas [START_REF] Cachia | A primal sketch of the cortex mean curvature: a morphogenesis based approach to study the variability of the folding patterns[END_REF]. The secondary and tertiary cortical folds appear to be more variable in appearance and in their relationship with functional areas [START_REF] Mohammed K Hasnain | Structure-function spatial covariance in the human visual cortex[END_REF][START_REF] Mohammed K Hasnain | Hemispheric asymmetry of sulcus-function correspondence: Quantization and developmental implications[END_REF]. Therefore, it is important to determine the hypothetical first cortical folds and verify their spatial distributions to further understand the anatomical and functional human brain development.

Figure 1.5: Chronology of gyrification in human. GW 31-32 divides phases 1 and 2. There is a gradual lack of consistency in cortical folding patterns among species in the final stages of gyrification, as minor developmental changes in gyri and sulci become more specialized to species and susceptible to local environmental and ultimate experiential variations. Extracted from [START_REF] Lewitus | Conical expansion of the outer subventricular zone and the role of neocortical folding in evolution and development[END_REF].

Regis et al. [START_REF] Régis | sulcal root" generic model: a hypothesis to overcome the variability of the human cortex folding patterns[END_REF] introduced sulcal roots to represent the first cortical folding locations. They proposed the shape and location of the sulcal roots of each person may be stable during the fetal period, and the variability during adulthood may be caused by the chaotic behavior of the folding process. It has been speculated that the primary cortical folds develop into the deepest local zones of sulci without spatial variance during development, which are called sulcal pits [START_REF] Im | Spatial distribution of deep sulcal landmarks and hemispherical asymmetry on the cortical surface[END_REF][START_REF] Lohmann | Deep sulcal landmarks provide an organizing framework for human cortical folding[END_REF]. The sulcal pits can be identified in a sulcal catchment basin through the structural information of the small gyri which are buried in depths of sulci.

Quantitative description of cortical folding

A variety of characteristics can be used to describe cortical folding. For example, global or local gyrification index, curvature, curvature-based spectral analysis, sulcal depth, sulcal length, sulcal area and sulcal fundus curves [START_REF] Philipp | Measures of folding applied to the development of the human fetal brain[END_REF][START_REF] Boucher | Depth potential function for folding pattern representation, registration and analysis[END_REF][START_REF] Fish | Influences of brain size, sex, and sex chromosome complement on the architecture of human cortical folding[END_REF][START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF][START_REF] Hill | Similar patterns of cortical expansion during human development and evolution[END_REF][START_REF] Im | Brain size and cortical structure in the adult human brain[END_REF][START_REF] Im | Fractal dimension in human cortical surface: multiple regression analysis with cortical thickness, sulcal depth, and folding area[END_REF][START_REF] Chiu-Yen | A geometric method for automatic extraction of sulcal fundi[END_REF][START_REF] Li | An automated pipeline for cortical sulcal fundi extraction[END_REF][START_REF] Lyu | Spectral-based automatic labeling and refining of human cortical sulcal curves using expert-provided examples[END_REF][START_REF] Mangin | A framework to study the cortical folding patterns[END_REF][START_REF] Pienaar | A methodology for analyzing curvature in the developing brain from preterm to adult[END_REF][START_REF] Schaer | A surface-based approach to quantify local cortical gyrification[END_REF][START_REF] Shi | Hamilton-jacobi skeleton on cortical surfaces[END_REF][START_REF] Shimony | Comparison of cortical folding measures for evaluation of developing human brain[END_REF][START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF][START_REF] Hyuk | Automated sulcal depth measurement on cortical surface reflecting geometrical properties of sulci[END_REF][START_REF] Zilles | The human pattern of gyrification in the cerebral cortex[END_REF]. The temporal changes of gyrification index were observed during fetal life [START_REF] Armstrong | The ontogeny of human gyrification[END_REF]. The spectral analysis of gyrification has been adopted in [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF] to parcellate the cortical surface into different spatial frequency bands according to curvatures. Cortical surface curvatures and sulcal depth have been measured at the overall level to quantify the degree of cortical folding from 22 to 39 weeks of gestational age [START_REF] Clouchoux | Delayed cortical development in fetuses with complex congenital heart disease[END_REF][START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF][START_REF] Dubois | Mapping the early cortical folding process in the preterm newborn brain[END_REF]. The changes in curvature and depth of cortical folding and local cortical expansion at the vertex level have also been observed using volume-or surfacebased registration techniques [START_REF] Clouchoux | Delayed cortical development in fetuses with complex congenital heart disease[END_REF][START_REF] Habas | Early folding patterns and asymmetries of the normal human brain detected from in utero mri[END_REF][START_REF] Schwartz | Modeling fetal cortical expansion using graph-regularized gompertz models[END_REF][START_REF] Scott | Volumetric and surface-based 3d mri analyses of fetal isolated mild ventriculomegaly[END_REF]. The methods based on local features (e.g., several curvature-based measures, sulcal depth and sulcal length) can capture local information, such as the depth or wideness of a sulcus or gyrus, while global measures (e.g., global gyrification index) can provide a overall description of the complexity of folding patterns.

Biomechanical framework of cortical folding

Without going back to a history of the understanding of the cortical morphogenesis, at the ontogenic level, we can nonetheless group together the different hypotheses and models according to some large families of explanation, by following for example the diagram proposed in In this section, we first present the concepts in biomechanics and then introduce the recent or well-known models of cortical folding development based on the biomechanical hypotheses of ( 2), ( 3) and (4).

Concepts in biomechanics

Biomechanics is the application of mechanical principles and methods for studying biological and medical problems. With the advances in non-linear continuum mechanics, modern continuum biomechanics emerged in the mid-1960s and rapidly grew with the development of computational methods and computer technologies, which enabled the solution of complex initial-boundary value problems and the implementation and interpretation of complex biomechanical experiments.

Deformation of biological tissues are governed by physical principles (the principle of force equilibrium), constitutive laws (material properties) and other non-mechanical processes [START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF]. In this section, we introduce some fundamental concepts in mechanics and the relation between mechanics and growth. It is necessary to keep these concepts in mind while solving biomechanical problems.

Force equilibrium

For a rigid body, Newton's second law of motion can be expressed as follows (when mass doesn't change in motion):

i f ext i = m a, ( 1.5)
where f ext i is external force. The total force (or resultant force) i f ext i determines the product of mass times acceleration. It is well known that if the total force equals to zero, this object remains at rest or continues to move at a constant velocity. If the total force doesn't equal to zero, this object will accelerate or decelerate with (a = i f ext i /m). For a deformable object, under external loading, some parts move more quickly than other parts and it ends up with changing the shape of this object. Thus for different parts of this object, their deformation could be different. It requires to apply Newton's equation for deformable bodies, i.e., the elastodynamic equation:

∇ • σ + F = ρ ¨ u, (1.6)
where σ is the Cauchy stress tensor, F is the body force per unit volume, ρ is the mass density and u is the displacement vector. The above equation is written on an infinitesimal cube inside the deformable body and is valid in the entire deformable body. Furthermore, in static case, there is no motion which means velocity and acceleration are zero. For a slow process where acceleration ( ¨ u → 0) and velocity are so small that can be neglected is called quasi-static case. The ratio of kinetic energy to internal energy can also be used to determine whether the system is quasi-static or not. When the ratio is small enough, the kinetic energy can be ignored and the system can be regarded as quasi-static.

The interest of solving dynamic problems (such as Equation 1.6) is that they tell the time history of structural deformation (such as brain growth) using time integration methods. One can get all the state quantities (position, velocity and acceleration) for each moment [START_REF] Lanczos | The variational principles of mechanics[END_REF].

Stress

Stress σ is a measure of force, normalized by the area σ = F A . The unit of stress is N m -2 in 3D or N m -1 in 2D. It is defined on an infinitesimal volume shown in Figure 1.7. Stress can be classified as normal stress (σ ii , with i = 1, 2, 3) and shear stress (σ ij , with i, j = 1, 2, 3 and i = j). σ ij means the force is applied on the surface i and along the direction j. It is noticed that stress tensor is symmetric σ ij = σ ji . This is so-called Cauchy stress tensor. In the next paragraph stress is Cauchy stress if we don't specify.

Strain and deformation gradient

Strain is a measure of deformation, which represents the displacement between particles in the body relative to a reference length. Strain is a non dimensional quantity. As shown in 1.8, the vector X defines the undeformed reference configuration and x defines the deformed configuration. The displacement u of any point inside the body can be defined as

u = x -X (1.7)
Then, the strain is the derivative of the displacement with respect to the reference configuration:

= ∂ u ∂ X = F -I (1.8)
where

F = ∂ x ∂ X
is the deformation gradient. When the deformation is small enough, the above strain is approximated as ij = 1 2 (u i,j + u j,i ) [START_REF] François | Mechanical Behaviour of Materials: Volume II: Fracture Mechanics and Damage[END_REF], where i, j are coordinates and u i,j means the partial derivative of i-direction displacement with respect to the direction j. However, the brain growth is far than a small deformation process. Therefore, large strain theory (finite strain theory) has been developed to deal with deformations in which strains and/or rotations are large enough beyond the infinitesimal strain theory. Depending on the deformation situation, several deformation tensors are used in the computational mechanics based on the deformation gradient F. Therefore, the deformation gradient tensor F is more commonly used.

Constitutive relations

A constitutive equation or constitutive relation is a relation between two physical quantities. In mechanics, constitutive relations are set up on the concept of stress and strain (or deformation gradient), it is also known as strain-stress relation. The simplest law is one-dimensional hook's law (linear elastic law), stress is proportional to strain: where E is Young's modulus. Here, the term linear means the proportional property, and the term elastic means the ability of a material to resume its normal shape after being deformed. In biomechanics, tissues are not simply linear elastic materials but non-linear elastic materials. The non-linear properties of tissues can be brought into the frame of hyperelastic materials for which the stress-strain relationship derives from a strain energy density function [START_REF] Gao | A finite strain nonlinear human mitral valve model with fluid-structure interaction[END_REF][START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF][START_REF] Ogden | Non-linear elastic deformations[END_REF][START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF].

σ = E (1.9)
For instance, for a compressible material, stress is given by

σ = 1 J ∂W ∂F F T , (1.10)
where W is strain-energy density and J = det(F).

In Figure 1.9, we can clearly see the difference of strain-stress relation between linear elastic and hyperelastic. Since the slope of the curve for the hyperelastic case is not always constant, meaning that applying same force can lead to different displacement values. The Neo-Hookean and Mooney-Rivlin solids are the first hyperelastic models derived by Ronald Rivlin and Melvin Mooney [START_REF] Mooney | A theory of large elastic deformation[END_REF][START_REF] Ogden | Non-linear elastic deformations[END_REF]. After that, other hyperelastic models have been given such as the Ogden model and the Arruda-Boyce model [START_REF] Ellen | A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials[END_REF][START_REF] William | Large deformation isotropic elasticity-on the correlation of theory and experiment for incompressible rubberlike solids[END_REF]. Apart from elasticity, brain tissue is also considered to be viscoelastic [START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF][START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF], which means that the existence of a velocity-depending force needs to be considered.

Stability

Force equilibrium can be stable or unstable. In the case of unstable equilibrium, a small disturbance will produce a large change. A typical unstable phenomenon is buckling as shown in Figure 1.10 (c), mechanical instability occurs on the surface of a uniformly compressed soft solid. This instability arises under sufficient compression (F > F c where F c is critical force) leading to the folding of the soft surface. In recent research, gyrification of brain is interpreted as a result of mechanical instability of a soft tissue that grows non-uniformly [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

a strain energy density function (the work per unit volume required to deform the material). Several models of brain folding assume brain tissue is elastic [START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF] or hyperelastic [START_REF] Chang | Micropatterned stretching system for the investigation of mechanical tension on neural stem cells behavior[END_REF][START_REF] Chenn | Regulation of cerebral cortical size by control of cell cycle exit in neural precursors[END_REF]. However, brain tissue is widely understood to be viscoelastic [START_REF] Je G Chi | Gyral development of the human brain[END_REF], such that sustained tension can lead to stress relaxation as the tissue lengthens passively. Stressinduced tissue growth can mimic viscoelastic relaxation. In figure 3a, elastic stretch (tension) leads to positive growth and a new lengthened configuration with no stress. Conversely, compression typically leads to negative growth and a shortened configuration with no stress. Some models of brain folding, such as the model in figure 5, incorporate viscoelastic behaviour explicitly to account for growth [START_REF] Cipolla | Visual motion of curves and surfaces[END_REF][START_REF] Clouchoux | Delayed cortical development in fetuses with complex congenital heart disease[END_REF][START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF].

Mechanical feedback: Cells both generate and respond to mechanical force [START_REF] Colantuoni | Temporal dynamics and genetic control of transcription in the human prefrontal cortex[END_REF]. For example contractile forces arise from interactions of myosin motor proteins and actin filaments. Actin polymerization itself can exert forces as filaments lengthen [START_REF] Courant | On the partial difference equations of mathematical physics[END_REF][START_REF] Dardenne | Variational tetrahedral mesh generation from discrete volume data[END_REF]. Proliferation, maturation of cells and accumulation of extracellular matrix material all tend to increase tissue dimensions (tissue growth), which will generate forces against constraining boundaries. On the other hand, applied forces can influence motor stepping rates, actin polymerization and gene expression, which in turn affect the processes of proliferation and growth.

Stability: Equilibria can be stable or unstable. For example, if a pure axial load is applied to a rod, it will shorten slightly but stay straight (figure 3b). If the load is large enough (above a critical level), the straight, compressed solution still exists, but lower-energy buckled states may also exist. Any slight perturbation will lead to a large deformation (an instability) and the rod will take a new, buckled shape (figure 3c). Importantly, the buckled shape is not different for every perturbation, but reflects a characteristic 'mode' of instability. However, the specific state that emerges, i.e. whether the rod buckles to the left or to the right, may depend on the perturbation. In addition, any imperfection in the rod (deviation from perfect straightness and symmetry) will affect the final configuration. In the context of the developing cerebral cortical sheet, deformations resulting from compressive stress in the cortical plate due to surface area expansion would be anticipated to adopt the more stable, low-energy configurations.

Models of cortical folding (a) Overview of hypothetical folding mechanisms

Four hypotheses have dominated research and speculation into the mechanisms of cortical folding: (1) the skull constrains growth of the brain and causes compressive stresses and buckling; (2) tension in axons connecting adjacent regions of the cortex draws those regions together to form gyri; (3) tangential expansion of outer cortical layers, greater than in inner layers, causes folding by a mechanical instability; (4) programmed patterns of growth cause more neurons to emerge, proliferate, reach the cortex and expand in some regions (gyri) more than others (sulci). As shown in figure 4, these models are not necessarily mutually exclusive. Axial force F is applied for a one-dimensional hypothetical element of tissue beam. In all cases, two equal and opposite forces are applied in the both ends of the beam, thus the resultant force is zero. (a) Under tension, beam is stretched and elongated. In the case of elastic deformation (middle panel), there exists a tensile stress (σ > 0). In the case of viscoelastic deformation or growing tissue (bottom panel), this stress tends to zeros at longtime limit due to the mechanical relaxation. (b) Under compression, beam is compressed and shortened. In the case of elastic deformation (middle panel), there exists a compress stress (σ < 0). In the case of growing tissue (bottom panel), this stress will be relaxed to zero. (c) Under sufficiently high compression, F > F c where F c is critical force, this beam may suddenly change shape and it shows up curvature. Such a phenomenon is called buckling or mechanical instability. Extracted from [START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF].

l * > 1 G(s) > 1 s > 0 s = 0 s < 0 s = 0 l * < 1 G(s) < 1 s o > 0 s m = ~0 s m = ~0F > F c l * m = ~1 s i < 0

Mechanics and growth

Scientific interest in using mechanics to understand the fundamental aspects of biological systems could be traced back at least to the beginnings of modern science itself. The idea that mechanics shapes organisms and should play a fundamental role in the description of growth was first proposed by Wilhelm His in 1888, as shown in Figure 1.11 [START_REF] His | On the Principles of Animal Morphology: Letter to Mr[END_REF][START_REF] Hopwood | Producing development: The anatomy of human embryos and the norms of wilhelm his[END_REF], then was echoed by D'Arcy Wentworth Thompson in 1942 [START_REF] Wentworth Thompson | On growth and form. On growth and form[END_REF] and Le Gros Clark in 1945 [START_REF] Le Gros | Essays on growth and form presented to d'arcy wentworth thompson[END_REF] while understanding brain morphology.

Stress influences the growth of soft tissues [START_REF] Goriely | The mathematics and mechanics of biological growth[END_REF]. The growth of soft tissues (changes in quality or material properties) often requires a cellular response. Mechanical transduction is a process by which mechanical signals applied in the bulk or at the boundaries of the body, are mediated down to the cellular level and then passed to the nucleus. This is a complex process that has not been fully understood. However, the influence of the mechanical environment on the regulation of the growth is undeniable, and the effects of stress on growth and physiological regulation of soft tissues are well documented [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF][START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Edward K Rodriguez | Stress-dependent finite growth in soft elastic tissues[END_REF][START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF].

In turn, stress is also influenced by growth [START_REF] Goriely | The mathematics and mechanics of biological growth[END_REF]. In many biological tissues, different parts of the tissue experience different growth rates due to a combination of cellular, chemical, and mechanical factors. The end result of this differential growth is that the tissue may be under stress even when unloaded, this stress field is called residual stress. The generation of residual stress through differential growth is a key feature of any mechanical theory of growing bodies [START_REF] Skalak | Compatibility and the genesis of residual stress by volumetric growth[END_REF]. For human brain, the differential growth of the grey matter constrained by the white matter generates residual stress in the grey matter layer and puts it into mechanical compression [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Essentially, as growth occurs, certain parts of the body need to be stretched or compressed to ensure integrity (no cavitation) and compatibility (no overlap).

Biomechanical folding hypotheses and models

Axonal tension-driven folding

The axonal tension hypothesis [START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF], which mechanical axonal tension along the cortico-cortical connections is the main driving force for the development of cortical folding, is illustrated in Figure 1.6 and 1.12 a. Strongly interconnected cortico-cortical connections, which generate tangential force components, could form gyri. This hypothesis is consistent with efficient wiringaxons that connect related areas will draw them together, thus decreasing the total length of neuron-to-neuron connections in the brain [START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF]. Disturbances in folding accompany alterations in axonal connectivity, as observed by [START_REF] David C Van Essen | Symmetry of cortical folding abnormalities in williams syndrome revealed by surface-based analyses[END_REF] who analyzed cortical folding in subjects with a genetic disorder known as Williams syndrome. [START_REF] Rakic | Specification of cerebral cortical areas[END_REF] noted similar effects in the brains of monkeys with disrupted visual systems. However, in both cases the causality and mechanism are not clear. The axonal tension hypothesis is particularly interesting because tension can be measured in individual axons [START_REF] Chada | Cytomechanics of neurite outgrowth from chick brain neurons[END_REF][START_REF] Steven R Heidemann | Tension as a regulator and integrator of axonal growth[END_REF][START_REF] Steven R Heidemann | Mechanical tension as a regulator of axonal development[END_REF][START_REF] Harish | Tension and compression in the cytoskeleton of pc 12 neurites[END_REF] as well as in tissue [START_REF] Xu | Residual stress in the adult mouse brain[END_REF][START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. However, the hypothesis that the walls of gyri are drawn together by axonal tension is not consistent with observed directions of stress in the ferret brain [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. Cuts of the folded brain along the radial axis of each gyrus remained closed, indicating that no tension pulls the related walls together. On the other hand, circumferential cuts perpendicular to axis of each gyrus were opened. In other words, axonal tension is along, not across gyri [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF] (Figure 1.12 b).

Differential tangential growth-driven folding

Stiff surface buckling on elastic foundations During folding, the cortical plate must expand faster than subcortical brain layers to attain its large surface area relative to brain Figure 1.12: a. A schematic representation of the axonal tension-driven folding hypothesis [START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF] from [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF]. b. Actual distributions of axon tension based on dissection and histology data [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. c. Intracortical differential tangential growth model due to [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF]. d. This differential growth model predicts elastic buckling of the outer layer [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF]. Extracted from [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF].

volume. The "differential tangential growth", meaning that outer layers grow tangentially at faster rates than inner layers of the brain, was proposed to explain cortical folding [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF] (Figure 1.12 c and d). The folding is based on two mechanical principles: (i) tangential expansion of an elastic layer, connected to an elastic foundation which does not expand, induces tangential compression in the expanding layer; (ii) a thin elastic layer under sufficiently large compressive stress will become unstable and buckle. If the thin layer is supported by an elastic foundation, the buckled shape will be sinusoidal, with the wavelength determined by the relative stiffnesses of the layer and foundation. Increasing the material stiffness, or elastic modulus, of the foundation leads to shorter wavelengths of the thin layer.

Richman et al. [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF] proposed the first quantitative model of folding based on differential growth. They analysed a model with two connected, growing layers of similar elastic material (equal elastic modulus in each layer, E1 = E2) connected to a soft elastic core (E core = E 1 /10). Growth was largest in the outer layer, so that compressive stress was induced in that layer. With these parameters, the authors obtained a theoretical prediction of wavelength that was consistent with the average wavelength (gyral peak-to-peak distance) of folding in the human brain.

Some studies have shown that the material properties in the outer and inner layers of the brain appear to be very similar. For instance, Xu et al. [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF] measured the regional mechanical properties of the ferret brain by indentation but found no significant differences between cortical and subcortical regions during the period of folding. Data from other studies of sheep [START_REF] Feng | Measurements of mechanical anisotropy in brain tissue and implications for transversely isotropic material models of white matter[END_REF] and pig brains [START_REF] Michael | Regional, directional, and age-dependent properties of the brain undergoing large deformation[END_REF] show also similar mechanical properties in grey and white matter.

Recent modeling works have expanded the set of results for buckling of stiff elastic layers on softer continua. For example, Hohlfeld and Mahadevan [START_REF] Hohlfeld | Unfolding the sulcus[END_REF] analyzed the surface instability in a soft, hyperelastic, incompressible material's strip with a thin, stiff skin when bent. This instability can result in furrows that resemble cortical folds. The buckling of the elastic shell on an elastic foundation has also been used to explain the formation of ridges and whorls during fingerprint formation [START_REF] Kücken | Fingerprint formation[END_REF]. Some buckling models have been proposed that do not attempt to clarify the physical processes of brain deformation, but rather rely on analogies to simpler behavior. A highly simplified, and the large-amplitude buckling model of cortical folding was introduced by Raghavan et al. [START_REF] Raghavan | A continuum mechanicsbased model for cortical growth[END_REF]. They represent the cortex as a thin elastic beam and predict the buckled shape by minimizing the energy density function, including the bending and tensile strain energy of the beam, the skull constraint, and other interaction terms that prevent self-contact and looping.

Elastic instability in growing continua

In the past two decades, biomechanical continuum models have been established and used to understand development and morphogenesis of biological tissues. Especially due to a theory of finite growth (a general continuum formulation for finite volumetric growth in soft elastic tissues) [START_REF] Edward K Rodriguez | Stress-dependent finite growth in soft elastic tissues[END_REF], which is based on a particular multiplicative decomposition of the deformation gradient, has been widely used to simulate the development and growth of cardiovascular tissue [2,[START_REF] Göktepe | A multiscale model for eccentric and concentric cardiac growth through sarcomerogenesis[END_REF][START_REF] Larry A Taber | Mechanics of cardiac looping[END_REF], skin [START_REF] Alexander M Zöllner | Growing skin: tissue expansion in pediatric forehead reconstruction[END_REF][START_REF] Alexander M Zöllner | On the biomechanics and mechanobiology of growing skin[END_REF] and brain [START_REF] Karel | Process for producing a chlorinated c3-6 alkane[END_REF][START_REF] Victor | Mechanics of head fold formation: investigating tissue-level forces during early development[END_REF][START_REF] Xu | Residual stress in the adult mouse brain[END_REF][START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. Their stress-dependent growth law is defined as a relation between a symmetric growth-rate tensor and a stress tensor, and this theory lays the foundation for several recent mathematical models of folding and creases, though these models rarely directly address the problem of cortical folding.

Using standard continuum mechanics terms, the position of the material element in the reference configuration is denoted as X, and the corresponding position of the same element in the deformed configuration is denoted as x. The deformation gradient tensor F = ∂x/∂X is expressed as the product of the growth tensor G and the elastic deformation gradient tensor F * :

F = F * • G.
(

Biological tissues are usually modeled as a hyperelastic material [START_REF] Dubois | Mapping the early cortical folding process in the preterm newborn brain[END_REF][START_REF] Shimony | Comparison of cortical folding measures for evaluation of developing human brain[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF][START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF][START_REF] David | A tension-based theory of morphogenesis and compact wiring in the central nervous system[END_REF], thus the Cauchy stress tensor σ is directly dependent on the elastic deformation according to the constitutive relationship:

σ = J * -1 ∂W ∂F * • F * T , (1.12)
where W is the strain energy density function of the material and J * = det(F * ) is the volume ratio of the elastic deformation. In general, an isotropic or almost incompressible isotropic constitutive behavior is assumed. The strain energy depends only on the elastic deformation, thus for example, in a nearly incompressible neo-Hookean material model, the material behavior is defined as

W = µ 2 (I * 1 J * -2/3 -3) + K 2 (J * -1) 2 . (1.13)
The strain energy depends on the shear modulus µ, on the bulk modulus K (K >> µ), on J * and on the first invariant (trace) of the right Cauchy-Green deformation tensor,

I * 1 = T r(C * ), where C * = F * T • F * .
The growth tensor G can be specified as a function of time or position, or it can be assumed that it depends on stress, thereby providing a mechanism for mechanical coupling by growth.

The growing elastic continuum model can show folding. Ben Amar and Goriely [START_REF] Ben | Growth and instability in elastic tissues[END_REF] used this approach to analyze the deformation of a growing incompressible hyperelastic spherical shell. They showed that growth leads to elastic instability with buckling modes of different wavelengths, depending on the growth parameters and the thickness of the shell. In this analysis, growth is specified (i.e., it does not depend on the stress in the material). Dervaux et al. [START_REF] Dervaux | Shape transition in artificial tumors: from smooth buckles to singular creases[END_REF] built a physical model of differential growth using circular slices of gel, in which the outermost layer expanded by swelling. This outer layer exhibited folding with wavelength determined by the thickness and modulus of the outer layer. Theoretical predictions agreed closely with the observed shapes. Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] extended experimental studies of differential growth into 3D, first using hemispheres of layered swelling gel (elastomer polydimethylsiloxane) in which the outer layer was expanded by swelling induced by immersion in hexane. By controlling the modulus and thickness of the swelling surface relative to the radius of the hemisphere, they obtained gyral and sulcal convolutions on the outer layer similar to those in fetal brains. Numerical experiments, based on a soft hyperelastic solid and in which the starting point approximated the smooth fetal human brain, resulted in 3D gyrification consistent with observations of fetal brains. These authors [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] also performed numerical simulations of growing hyperelastic materials and obtained shapes that closely matched their experimental observations.

All the studies in this thesis are based on the brain folding model proposed by Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], assuming that the differential tangential growth leads to the elastic instability in growing continua and produces the folds similar to those in developing fetal brains.

Stress-induced growth and viscoelastic instability

Models based strictly on uniform surface growth on an elastic foundation are not able to replicate two key features of mammalian cortical folding. Firstly, primary folds that are highly conserved between individuals have not been consistently and accurately reproduced. Secondly, models fail to recapitulate the observed growth of subcortical white matter. Developmental neuroanatomists emphasize that more superficial subcortical layers assume the folded shape of the cortex but deeper layers remain smooth [START_REF] Borrell | How cells fold the cerebral cortex[END_REF][START_REF] Smart | Gyrus formation in the cerebral cortex in the ferret. i. description of the external changes[END_REF]. White matter structure is essential for normal brain development, so useful models should be able to capture the structural development of the subcortical layers.

In vitro studies have shown that axons (the primary component of white matter) grow in response to tension and reach an equilibrium length that maintains a small, but finite, level of tension [START_REF] Bray | Axonal growth in response to experimentally applied mechanical tension[END_REF][START_REF] Chada | Cytomechanics of neurite outgrowth from chick brain neurons[END_REF][START_REF] Dennerll | The cytomechanics of axonal elongation and retraction[END_REF][START_REF] Steven R Heidemann | Tension as a regulator and integrator of axonal growth[END_REF]. Therefore, the white matter layer will probably respond to tensile stress by growing in the direction of tension and maintaining tension along the direction of axonal fibres. This prediction is consistent with the oriented tension observed in white matter [START_REF] Xu | Residual stress in the adult mouse brain[END_REF][START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. Since tension-induced axonal growth takes time, stress induced by relatively fast loading will relax with a characteristic time constant. This behaviour mimics the response of a viscoelastic material, where the relaxation time depends on the ratio of "viscous" resistance to elastic stiffness. If the inner subcortical tissue is modeled as a viscoelastic material, a distinct type of viscoelastic instability (buckling on a viscoelastic foundation) occurs in response to the growth of the cortical plate [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF]. Because a viscoelastic foundation appears stiffer when loading is applied faster, the occurrence of buckling and the wavelength of the resulting surface depends on the rate of cortical growth relative to the rate of relaxation in the subcortical region.

Toro et al. [START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF] incorporated stress-induced growth into a folding model by including discrete, radial elasto-plastic fibres (representing the subcortical tissue zones) connected to an expanding elastic ring (representing the cortical plate). Simulations reproduced many qualitative features of folding, especially the evolution of a wavelike folding pattern in the outer ring. They modulated the emergence and wavelength of folds by changing the growth rate, cortical thickness and stiffness of the outer layer. Folds were produced at specific locations through geometrical perturbations (small changes in the initial shape) or spatial variations in model parameters (for example, thickness or stiffness). This model clearly showed that simple physics can produce a wealth of folding behaviour, and illustrated the qualitative influences of physical parameters.

Moving beyond one-dimensional, discrete elements, continuum mechanical models can more accurately represent the behaviour of tissue in 2D or 3D, as they capture stresses and deformations in multiple directions. By representing the sub-cortical region as a continuous viscoelastic foundation, Bayly et al. [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF] obtained an analytical formulation for the effects of relative growth rate on the wavelength of cortical folds. The tangential dimension of the outer (cortical) layer was defined to grow at a rate of G T,0 (week -1 ). For example, if the tangential dimension dou-bles in one week then tangential growth G T = 2 (dimensionless) and G T,0 = 1week -1 . In the inner (subcortical) layer, growth was assumed to occur in both radial and tangential directions (G R and G T ), but only in response to radial and tangential stresses (σ R and σ T ). This is a form of stress relaxation, described by the rate constant, R. In this viscoelastic material framework, any initial stress (σ 0 ) will decrease to the equilibrium stress (σ e ) according to the form σ = σ e + (σ 0 -σ e )exp(-Rt) [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF]. Thus, the dimensionless ratio Γ G = G T,0 /R describes how fast the cortex grows (G T,0 ) relative to the subcortical, stress-dependent growth response (R). Since growth occurs slowly (on the time scale of days or weeks), inertial effects are negligible, and cortical folding is treated as a quasi-static process.

The formula for folding wavelength (λ), as a function of relative growth rate (Γ G ) cortical thickness (h), and ratio of elastic moduli (β = µ/µ f ) derived by Bayly et al. [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF] is:

λ = 2πh[ β(Γ + 1) 3Γ ] 1/3 , (1.14)
where Γ is a root of the characteristic polynomial: Γ 5 -(64/9)β 2 Γ 3 G (Γ 2 + 2Γ + 1) = 0. This formula predicts that wavelength increases with cortical thickness and modulus (µ) of cortex relative to modulus (µ f ) of subcortical regions. With similar stiffness in the cerebral cortex and subcortical regions, this model can account for the changes in wavelengths of folds and the stress fields observed in the developing brain. Furthermore, the growth rate of the cortex, which is relative to how fast the core grows in response to stress, largely determines the wavelength of cortical folds. 2D simulations confirmed that the larger the ratio of cortical to subcortical growth rates (Γ G ), the shorter the folding wavelength (Figure 1.13). The initial shape before tangential expansion and the spatial changes in tangential expansion itself may also affect the final shape. The effects of the cortical thickness, cortical growth rate and elastic moduli in the two layers are consistent with experiments and simulations from other groups (i.e., Dervaux et al. [START_REF] Dervaux | Shape transition in artificial tumors: from smooth buckles to singular creases[END_REF]; Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]; Budday et al. [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF]; Toro et al. [START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF]).

Budday et al. [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF] modeled the stress-dependent growth by a morphogenetically growing outer surface and a stretch-driven growing inner core. This method integrates two hypotheses for cortical folding: axonal tension and differential growth. Their model confirmed that the ratio between growth rates in outer and inner layers influenced folding wavelength, and predicted that the misbalance in cortical and subcortical growth (relative growth rate) and cortical thickness cause morphological abnormalities: a slower growing or thinner cortex generally enhances folding and a faster growing or thicker cortex reduces folding. These characteristics are in agreement with the classical pathologies of folding including polymicrogyria and lissencephaly. This model was further extended by Holland et al. [START_REF] Holland | Emerging brain morphologies from axonal elongation[END_REF] to cover subcortical anisotropy, because that more axon elongation occur in the direction of axonal orientations than in other directions. The simulation predicted that the orientation of the subcortical axons can change the location of gyri and sulci, proposing that tissue anisotropy is another important factor.

Cellular mechanisms of differential tangential growth

The mathematical modeling of differential tangential growth may be beyond the current understanding of cortical and subcortical growth at the cellular level, such as that the mechanism by which the cortical plate expands tangentially (in area) but not much in the radial direction (in thickness) is not clearly explained. Since the folding period is later than the birth, proliferation, and migration phases of most neurons, these events may lay the foundation for cortical expansion, rather than drive it directly.

Richman and co-authors [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF] initially believed that neuronal maturation, particularly dendritic arborization, may provide a main mechanism of tangential growth due to the consistence Figure 1.13: Effects of cortical growth rate on wavelength, subcortical growth, and stress in a 2d model of cortical folding. Each column contains spatial maps of a different variable superimposed on the deformed geometry: radial growth G r , tangential growth G t ; radial stress σ r , tangential stress σ t . Each row corresponds to a different scaled ratio of cortical to subcortical growth rates. Extracted from [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF].

with timing of cerebral cortical folding and the morphological differentiation of cells in the cortex [START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF]. The latest observations with diffusion-weighted MRI further support this possibility [START_REF] Nørhøj | Determination of axonal and dendritic orientation distributions within the developing cerebral cortex by diffusion tensor imaging[END_REF][START_REF] Christopher D Kroenke | Regional patterns of cerebral cortical differentiation determined by diffusion tensor mri[END_REF][START_REF] Wang | Folding, but not surface area expansion, is associated with cellular morphological maturation in the fetal cerebral cortex[END_REF]. The water diffusion on the cortical board is strongly anisotropic before folds appear in the cerebral cortex. Radial diffusion is less restricted than tangential surface diffusion, apparently because morphologically undifferentiated and tangentially oriented glial and neuronal cell processes limit diffusion perpendicular to their axes. As neurons mature and branch dendrites expand, these cell projections provide new ways for unrestricted diffusion in tangential directions, making diffusion more isotropic.

In ferrets, the timing of loss of water diffusion anisotropy [START_REF] Thomas | Cellular motility driven by assembly and disassembly of actin filaments[END_REF] coincides with morphological differentiation of pyramidal neurons, as assessed by Golgi staining approaches [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF][START_REF] Zervas | Ferret pyramidal cell dendritogenesis: changes in morphology and ganglioside expression during cortical development[END_REF]. Within-subject comparisons of water diffusion anisotropy and dendritic orientation complexity [START_REF] Nørhøj | Determination of axonal and dendritic orientation distributions within the developing cerebral cortex by diffusion tensor imaging[END_REF] further confirmed that the loss of anisotropy is associated with an increase in dendritic structure. Finally, longitudinal intrauterine diffusion tensor imaging measurements indicated that there is a correlation between the loss of diffusion anisotropy and folding in the fetal rhesus macaque brain [START_REF] Wang | Folding, but not surface area expansion, is associated with cellular morphological maturation in the fetal cerebral cortex[END_REF].

Another mechanism of cortical expansion is the late intercalation of previously born neurons in subcortical layers to the top of the cortical plate [START_REF] Striedter | Cortical folding: when, where, how, and why?[END_REF], which can be considered as a proximate force for tangential expansion that then leads to cortical folding. Several authors [START_REF] De | Coevolution of radial glial cells and the cerebral cortex[END_REF][START_REF] Kriegstein | Patterns of neural stem and progenitor cell division may underlie evolutionary cortical expansion[END_REF][START_REF] Reillo | A role for intermediate radial glia in the tangential expansion of the mammalian cerebral cortex[END_REF] noticed that in gyrencephalic animals, that is, the outer subventricular zone, a layer that is significantly different histologically is found. In this layer, basal radial glial cells (bRGCs) proliferate or produce intermediate progenitors, but bRGCs divide only once to directly generate two neurons in the lissencephalic mouse [START_REF] De | Coevolution of radial glial cells and the cerebral cortex[END_REF][START_REF] David V Hansen | Neurogenic radial glia in the outer subventricular zone of human neocortex[END_REF]. The increases in neuronal bRGCs and intermediate progenitors, can all produce neurons, increasing the final number of neuronal cells in the folded cerebral cortex. The bRGCs themselves seem to fan out [START_REF] Fernández | Cerebral cortex expansion and folding: what have we learned?[END_REF][START_REF] Reillo | A role for intermediate radial glia in the tangential expansion of the mammalian cerebral cortex[END_REF], perhaps providing a wider path for neuron migrations and allowing more to reach the cortex.

It is now important to determine the biomechanical implications of this biological difference lissencephalic and gyrencephalic species. Interestingly, anatomical experiments in mouse revealed tangential compressive stress in the cortex [START_REF] Xu | Residual stress in the adult mouse brain[END_REF][START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF], suggesting an increase in cortical expansion relative to the subcortex. Further researches are needed to understand how enhanced cortical neurogenesis mechanisms promote differential growth. Another consideration is that an increase in cortical thickness relative to the initial brain size increases the growth threshold of induced folding [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF][START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Mir | Role of mechanical factors in cortical folding development[END_REF].

Patterned growth and differential proliferation-driven folding

Heterogeneous patterned growth In addition to the growth differences between layers, the growth of the brain may be heterogeneous in each layer. Although the growth differences between layers are sufficient to cause folding, spatial patterns along a single layer have also been reported [START_REF] De | Genetic maps and patterns of cerebral cortex folding[END_REF][START_REF] De | Discrete domains of gene expression in germinal layers distinguish the development of gyrencephaly[END_REF][START_REF] Kriegstein | Patterns of neural stem and progenitor cell division may underlie evolutionary cortical expansion[END_REF]. If this spatial pattern is consistent among individuals of a given species, it may be based on consistently positioned primary folds that govern the geometry of the ferret brain and determine many features of human neuroanatomy.

In [START_REF] Smart | Gyrus formation in the cerebral cortex in the ferret. i. description of the external changes[END_REF][START_REF] Smart | Gyrus formation in the cerebral cortex of the ferret. ii. description of the internal histological changes[END_REF], the authors pointed out that the subcortical tissue expands heterogeneously in the radial direction and has the largest radial expansion under gyri, and gyrification is driven by non-uniform radial expansion of the subcortical layers. In contrast, in the differential tangential growth model of Bayly et al. [START_REF] Clouchoux | Delayed cortical development in fetuses with complex congenital heart disease[END_REF], the heterogeneous radial expansion of the subcortical tissue is caused by folding of the outer layer. How to distinguish whether the subcortical layers push the cortex outward to form gyri, or does the folding cortex pulls the inner layer outward through stress-induced growth? The presence of tension rather than compression along the radial axis of each gyrus provides direct physical evidence for the folding cortex that pulls the underlying gyral tissue [START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF].

The birth and proliferation of neurons and glial cells may be different under gyri and sulci, especially for primary folds. Kriegstein et al. [START_REF] Kriegstein | Patterns of neural stem and progenitor cell division may underlie evolutionary cortical expansion[END_REF] showed that proliferative layers was thicker under primary gyri than primary sulci. Later, the Kriegstein and Borrell groups [START_REF] David V Hansen | Neurogenic radial glia in the outer subventricular zone of human neocortex[END_REF][START_REF] Reillo | A role for intermediate radial glia in the tangential expansion of the mammalian cerebral cortex[END_REF] showed that proliferation in the outer subventricular zone (oSVZ) enriches neuronal density in the cortical plate. As mentioned above, Borrell and co-authors [START_REF] Reillo | A role for intermediate radial glia in the tangential expansion of the mammalian cerebral cortex[END_REF] have shown that basal radial glial cells (bRGCs) spread out in a tangential pattern in developing ferrets. These tangentially diffusing bRGCs are not found below sulci in the developing ferret brain or anywhere in the developing mouse brain. Substantial evidences suggest that there are indeed heterogeneous patterns of bRGC birth and proliferation, and that they play an important role in determining the location of primary folds.

Reconciling radial and tangential growth hypotheses

In the growth model of Bayly et al. [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF], the position of the primary gyrus was specified by a local increase in radial growth, but secondary folds evolved from mechanical instability. The model shows that the radial growth hypothesis and the differential tangential growth hypothesis are not mutually exclusive.

Another interesting possibility is that mechanical stress itself may trigger patterns of proliferation, differentiation, migration, and maturation. Studies in chicken embryos have shown that the precursor neuroepithelial progenitor cells of basal radial glial cells in the oSVZ proliferate in response to mechanical feedback [START_REF] Desmond | Focal adhesion kinase as a mechanotransducer during rapid brain growth of the chick embryo[END_REF][START_REF] Desmond | Internal luminal pressure during early chick embryonic brain growth: descriptive and empirical observations. The Anatomical Record Part A: Discoveries in Molecular[END_REF][START_REF] Kara | Contraction and stress-dependent growth shape the forebrain of the early chicken embryo[END_REF]. The mechanical properties of the substrate [START_REF] Adam J Engler | Matrix elasticity directs stem cell lineage specification[END_REF] and mechanical tension [START_REF] Chang | Micropatterned stretching system for the investigation of mechanical tension on neural stem cells behavior[END_REF] have been shown to influence the differentiation of stem cells into neurons. Franze and colleagues [START_REF] Franze | The mechanical control of nervous system development[END_REF][START_REF] David E Koser | Mechanosensing is critical for axon growth in the developing brain[END_REF] have shown that axons are mechanically sensitive, and some authors have observed that the appearance of neurites [START_REF] Chada | Cytomechanics of neurite outgrowth from chick brain neurons[END_REF][START_REF] Dennerll | The cytomechanics of axonal elongation and retraction[END_REF] and synapses [START_REF] Siechen | Mechanical tension contributes to clustering of neurotransmitter vesicles at presynaptic terminals[END_REF] are affected by tension.

For instability caused by uniform cortical growth, Bayly et al. [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF] showed that the mechanical stress becomes heterogeneous in the subcortex before any significant change in shape takes place. Local radial stress may be sufficient to induce radial growth below an expected gyrus. In this case, a primary fold may occur due to mechanical instability, but at a precise location, based on subtle factors, cause nonuniform stress in subcortical layers. As the geometry changes, so do the patterns of subcortical stress. Local stress below the cortex can cause radial growth and other secondary folds. In a similar manner, stress patterns, or folding patterns may be influenced by various physical factors including initial geometry, anisotropy of cellular process orientations [START_REF] Holland | Emerging brain morphologies from axonal elongation[END_REF], regional changes in mechanical properties [START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF], or regional variations in cortical growth.

On a larger scale, the MRI studies of ferret, macaque and human brains during folding reveal spatial-temporal changes in cortical area expansion, which is shown in Figure 1.14. Spatial changes appear to span multiple sulci and gyri, and patterns change over multiple days in ferret brains [START_REF] Andrew K Knutsen | Spatial and temporal variations of cortical growth during gyrogenesis in the developing ferret brain[END_REF] or several weeks in human brains [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF], consistent with the patterns of cortical maturation and functional development. By incorporating these patterns of tangential expansion, the 3D models could better generalize the primary folds of human brains [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], or some folding differences, which are observed in neurological diseases such as epilepsy, autism, schizophrenia and bipolar disorder.

Biomechanical modeling by Finite Element Method 1.3.1 Finite Element Method (FEM)

Introduction

Computers have changed our world and revolutionized the production mode. In the engineering field, traditional hand-made product design is replaced by computer-aided design (CAD). Not only because it avoid repeating work and improve accuracy, but also can it represent a final effect picture of the product before manufacturing a real one. Based on the CAD model, we can investigate other types of analyses. Finite element analysis is an important part of the CAD. It is known as a mature tool to solve engineering problems and mainly applied among structural mechanics, electromagnetism, acoustics, etc. Recently, FEM is more and more mentioned and employed in bio-mechanics science. Lots of frontier bio-mechanics research employed FEM as supplement, such as the human brain folding model proposed by Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Like other sciences, in biomechanics, the integration of experimental and analytical models is to gain knowledge and to understand mechanical response in vivo. Experience provide interpretable data in various situations. While these works are affected especially by achievements in the field of graphics processing and image processing technology, digital mechanics, genetics and molecular biology [START_REF] Naaji | The application of the finite element method in the biomechanics of the human upper limb and of some prosthetic components[END_REF]. On this point, advantages of FEM are as follows: FEM ensure the accuracy of the geometry representation; the possibility of parameter analysis; the possibility of remeshing and improvement of FE model meshes; the possibility of using FE models in other types of analysis.

FEM: a good candidate to solve PDE

Considering the mechanics response, we need to take into account the corresponding partial differential equations (PDE) such as Equation 1.6. How to solve PDE is the key point in these mathematical questions. I present here three basic numerical methods to solve PDE: finite difference method (FDM), finite element method (FEM) and finite volume method (FVM). They all involve subdividing the structure into a finite number of elements of simple geometry.

FDM is a simple method but it is far from enough facing to complex geometry. As shown in Figure 1.15 (a), discretization of equation must be rectangular stencil. Based on FDM, FEM and FVM are developed in the 20th century. FEM is very successful in the solid mechanics problems (most commonly treated in Lagrange reference frames) [START_REF] Belytschko | Nonlinear finite elements for continua and structures[END_REF][START_REF] Thomas | The finite element method: linear static and dynamic finite element analysis[END_REF], while FVM is believed more effective to treat the fluid/aerodynamic problems (most commonly treated in Eulerian reference frames) [4]. An element of FEM represents one part of the deformable solid, so it deforms and moves with the solid. An element of FVM represents a part of space. This element doesn't move, but fluid can flow in or out this element. The meshes of FEM and FVM are shown in Figures 1.15 (b) and (c) respectively. Density of elements will depend on the geometry complexity as well as the mechanical fields which implicit sometimes. For example, the mesh of FEM needs to consider the stress concentration due to the geometry. The FVM mesh near a airplane wing is shown in Figure 1.15 (c), mesh is denser and denser near the boundary layer in order to capturing the flow transition from laminar to turbulent.

Fundamental idea of the FEM is replacing a structure with a set of elements whose assemblage mesh can approximate its geometry [START_REF] Belytschko | Nonlinear finite elements for continua and structures[END_REF]. Inside each element, it is supposed that the solution is a linear combination of certain numbers of shape functions. The shape functions are well-defined for the standard elements (e.g., tetrahedron or hexahedron), and in practices, one can use them directly. The unknowns of the problem are the coefficients of these shape functions, called the degrees of freedom (DOF), which are, for example, the displacements (u) in our studies. As such, the elastodynamic equation (Equation 1.6) is discretized in space and can be solved efficiently.

For one dynamic problem (such as the brain growth), the solution of the equilibrium equation also depends on time. Therefore, we need to discretize the equilibrium equation in time and then solve it step by step [START_REF] Géradin | Mechanical vibrations: theory and application to structural dynamics[END_REF][START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Given the initial conditions at t 0 , one can get the results at the following steps t n+1 = t n + ∆t, n = 0, 1, 2, ... if we know the state of system at t n by writing a series of recursive functions. This procedure is called time integration. The choice of the time step ∆t is important to ensure the stability of calculation. ∆t of explicit integration schemes are subject to a restrictive condition Courant-Friedrichs-Lewy (CFL condition) [START_REF] Courant | On the partial difference equations of mathematical physics[END_REF]. CFL conditions require to take into account the spatial discretization of the system (mesh spacing) and the sound speed in the material, being related to the elastic modulus of the substance and the mass density. Because in a dynamic transient problem, force or displacement is transferred through the mechanical waves at sound speed (c) in the material. If ∆t is too large and c × ∆t could be longer than the distance between two adjacent points ∆L (i.e., c × ∆t > ∆L), it gives rise to the unreasonable solution and cause an computational instability of time integration.

Finite Element Method (FEM) codes

FEM codes are widely used to solve engineering problems. They can be classified into two categories [START_REF] Neupane | Comparison of some fem codes in static analysis[END_REF]:

• Commercial FEM codes • Open source FEM codes
Commercial FEM codes Commercial FEM codes are largely used in the companies for solving engineering problems and in many universities for educating. Most commercial codes have user friendly interface (GUI), CAD module, mesh generation module and exhaustive training and support, which make them popular among design office and industrial production. Nowadays, these commercial codes are not only limited solving engineering problems, but also used in the filed of physic, chemistry and biomedical engineering.

Commercial FEM codes guarantee a standard requirement of finite element analysis, they also provide a powerful solver and other new features, such as parallel computing. Because of copyright, one need to purchase a license before using them. The possibility of development of these commercial FEM codes is another problem which scientific researchers must take into account. We list here some commercial FEM codes:

• ANASYS • ABAQUS • COMSOL Multiphysics • Creo simulate • ADINA

Free or open source FEM codes

Free or open source FEM codes are also used by different academic institution, R&D companies. These source codes are written by groups of academic or individuals for specific purposes. These FEM codes are provided for free under GNU License and source codes are available to download. They are not suitable to solve wide range of engineering problems but could be used in some specific questions. They may haven't GUI, CAD module or mesh generation functions. Complicated geometries can't be modeled and discretized into meshes using these codes, thus complicated geometries are previously obtained by other CAD softwares and mesh generators. Documents of open source FEM codes can be accessed from their websites and one can also contact their developers directly in case of questions. We list here some free or open source FEM codes:

• CalculiX • Elmer • Gmesh • FEBio • Code_Aster • Cast3m
Among so many FEM codes, we did a comparison of performance (Table 1.1) via the site (http://www.feacompare.com/) and used FEBio (open source FEM code) to do the first experiment on the human brain growth. The main reason to use open source FEM code is that the code can be downloaded from internet for free and the sources is possible for development. FEBio is specifically designed for biomechanical applications, which solves the nonlinear finite element equations using a quasi-Newton method. It offers biologically relevant constitutive models and modeling scenarios. Compared with other open source codes, although Code_Aster (EDF) and Castem (CEA) are widely used in France, their application situation is traditional industry and the materials highlighted by their library are metal, alloy, beton, etc. While in FEBio, hyperelastic materials, such as Neo-Hookean materials, are emphasised which is often used to model the mechanical response of biological tissues. The FEBio software suite (Preview, FEBio, Postview) is a set of software tools for nonlinear finite element analysis in biomechanics and biophysics [START_REF] Maas | Febio: finite elements for biomechanics[END_REF]. One of the highlights of Preview is that it's possible to use and manipulate meshes of different formats generated by multiple softwares (such as .feb, .n, .inp, .k, .txt, .unv, .nas, .dxf, .stl, .hmascii, .surf, .msh, .byu, .mesh, .ele, .iges, .vtk, .raw, .mphtxt and .ply).

Code

Price CAD Large displacements Contact Hyperelastic ANASYS

> $4,000 Y Y Y Y ABAQUS $20,000 Y Y Y Y COMSOL $4000 Y Y Y Y CalculiX Free N Y Y Y Elmer Free N Y Y N FEBio Free N Y Y Y Code_Aster Free N Y Y Y Table 1.1: Comparison between different FEM codes

Biomechanical modeling

The first experiment aims to use FEBio to implement the biomechanical model of brain growth proposed in [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Since some parameters of the model need to be additionally defined (e.g., the relative tangential growth tensor), we also used the "GIBBON" (Geometry and Image-Based Bioengineering add-On), which is an open-source MATLAB toolbox and can be interfaced with FEBio. The flowchart of the interaction between MATLAB and FEBio is illustrated in Figure 1.16. The general definition of the deformation gradient is:

F = I + ∂u ∂X ,
where I is the unit tensor, u is the displacement of nodes and X is the generalized coordinates.

Especially, if the element used is linear triangles (2D) or linear tetrahedra (3D), this deformation gradient could be represented alternatively by the matrix of vectors describing the triangle or tetrahedron:

F = I + ∂u ∂X = I + (A -Â) Â-1 = I + A Â-1 -Â Â-1 = A Â-1 ,
where  is undeformed matrix of vectors and A is deformed matrix of vectors. For example, for 3D tetrahedrons elements, A is expressed as:

A =    x 1 -x 4 x 2 -x 4 x 3 -x 4 y 1 -y 4 y 2 -y 4 y 3 -y 4 z 1 -z 4 z 2 -z 4 z 3 -z 4   
However, the deformation gradient F defined in the model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is related to the growth tensor G as:

F = A(G Â) -1 = A Â-1 G -1
Therefore, we defined G by using the GIBBON in MATLAB and then added the term of G -1 in the definition of F in FEBio. At each time step, G should be imported into FEBio, and FEBio is responsible for calculating the deformation gradient (F ), volumetric strain energy density (W ), stress (σ) and deformed coordinates (x(u)). By using this flowchart (Figure 1.16), the simulations based on a 22 weeks' fetal brain mesh which was provided by Tallinen et al. have been performed. However, there were no obvious folds on the deformed cortical surface (Figure 1.17). In addition, we were not very clear how to solve the following three points. The first one is that does FEBio take into account the viscous damping in its dynamic system? Theoretically, the nodal resultant force of the brain mesh in a dynamic system should contain the viscous damping force to calculate nodal velocity and deformed coordinates. The second is that could FEBio solve a problem of self-contact between the folds on the cortical surface? The last is that FEBio sometimes detects negative Jacobian at certain elements after several iterations of deformation, it may be caused by the quality of the mesh, i.e., the degree of regularity of tetrahedra. The volumes of some irregular tetrahedra become zeros after certain iterations of deformation (Figure 1.18), which is described by the following equations:

h = L * sin θ, (1.15) if θ << π 2 , lim θ→0 h = lim θ→0 L * sin θ = lim θ→0 L * θ = 0, (1.16
)

lim θ→0 V = lim θ→0 1 3 * h * S = 0. (1.17)
V is the volume of an irregular tetrahedron after the deformation. The elements with close to zero or negative volumes cause difficulty in solving, and the deformed coordinates of mesh nodes tend to infinity. Last but not least, the interaction between FEBio and MATLAB takes a lot of time.

Considering the problems and difficulties mentioned above, we finally used pure Python (open source, code readability, concise syntax, extensive support libraries) to implement the biomechanical brain folding model proposed in [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], which is available at https://github. com/rousseau/BrainGrowth. The model is divided into three main parts, including an elastic force part, a contact force part and a displacement part. The calculation of each biophysical parameter is written in a separate function. More detailed descriptions of this model are shown in the following sections.

Conclusion

In this chapter, we first present the structure of human brain and the growth of living creatures, which serve to understand the expansion and folding process of the human brain. Then we briefly introduce the existing quantitative methods of cortical folding. In addition, we present the concepts in biomechanics and various models of cortical folding based on different biomechanical hypotheses. The biomechanical continuum model based on the hypothesis of differential tangential growth-driven folding can simulate gyral and sulcal convolutions on the cortical layer similar to those in fetal brains [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Therefore, all studies in this thesis are based on the model proposed in [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

In the model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the residual stress in the cortical layer that forms folds is caused by the differential growth of the cortex constrained by the white matter. The relation between the stress and the relative growth is expressed through the neo-Hookean constitutive relation. Considering that the biomechanical modeling in [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is based on the finite element method (FEM), thus we also outline FEM at the end of this chapter.

The biomechanical brain folding model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is implemented by us in Python, which is available at https://github.com/rousseau/BrainGrowth. Based on these FEM codes, we study the impacts of biophysical parameters onto folding patterns in the following sections.

Introduction

Human brain growth is accompanied by the folding of the cerebral cortex, which takes place in a hierarchical mode during gestational weeks 16-40 [START_REF] Sun | Growth and folding of the mammalian cerebral cortex: from molecules to malformations[END_REF], with primary folds forming the earliest and highly conserved, then secondary folds elaborating on these folds, etc [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF][START_REF] Christopher | How forces fold the cerebral cortex[END_REF]. Recent studies have shown that not only the molecular and cellular processes, but also mechanical forces play an important role in the formation of cortical convolutions [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Foubet | Mechanical morphogenesis and the development of neocortical organisation[END_REF][START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. It has been revealed that mechanical models based on the hypothesis of differential tangential growth could produce realistic folding patterns when they are applied on human fetal brain data [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. 3D numerical simulations of brain growth demonstrate that the relative tangential expansion of the cerebral cortex constrained by the white matter generates compressive stress, resulting in cusped sulci and smooth gyri similar to those in developing fetal brains [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

The cortical folding patterns are influenced by various physical parameters, e.g., the initial cortical thickness [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] David P Richman | Mechanical model of brain convolutional development[END_REF][START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the initial geometry [START_REF] Bohi | Global perturbation of initial geometry in a biomechanical model of cortical morphogenesis[END_REF][START_REF] Todd | A geometric model for the cortical folding pattern of simple folded brains[END_REF][START_REF] Toro | On the possible shapes of the brain[END_REF], the initial curvature of the surface [START_REF] Christopher | How forces fold the cerebral cortex[END_REF] and the relative growth [START_REF] Bayly | A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain[END_REF][START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF][START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Toro | A morphogenetic model for the development of cortical convolutions[END_REF][START_REF] Xu | Axons pull on the brain, but tension does not drive cortical folding[END_REF]. In addition to these recent observations, many questions are still open regarding the morphogenesis of folding patterns, including links between the physical parameters of simulation models and the folding patterns observed in in vivo MRI data. A deeper understanding of these parameters can significantly contribute to comprehend pathologies associated with characteristic changes in cortical folding. For instance, the neonatal arterial ischemic stroke may be related to abnormal early cortical development, and polymicrogyria, pachygyria, lissencephaly malformations can be accompanied by autism [START_REF] Auzias | Atypical sulcal anatomy in young children with autism spectrum disorder[END_REF][START_REF] Antonio Y Hardan | Increased frontal cortical folding in autism: a preliminary mri study[END_REF][START_REF] Wu Nordahl | Cortical folding abnormalities in autism revealed by surface-based morphometry[END_REF], schizophrenia [START_REF] Heidi | Pre-frontal lobe gyrification index in schizophrenia, mental retardation and comorbid groups: an automated study[END_REF][START_REF] Cachia | Cortical folding abnormalities in schizophrenia patients with resistant auditory hallucinations[END_REF][START_REF] Roger | Reduced cortical folding in individuals at high risk for schizophrenia: a pilot study[END_REF] or epilepsy [START_REF] Pang | Malformations of cortical development[END_REF]. These motivate us to study the effect of biophysical parameters on cortical folding patterns.

In this chapter, we use the cortical folding model proposed by Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] and build the simulations based on the tetrahedral meshes generated on ellipsoids. Since mesh density is a crucial issue in biomechanical modeling using finite element method, which is closely associated with the accuracy of the finite element model and determines its complexity degree, we first study the effects of variations in mesh density on surface morphology. Then based on an appropriate mesh density, we investigate the influence of the cortical growth, the initial geometry and the initial cortical thickness on surface morphology. From the simulation results, we first visually remark the difference in the appearance of folding patterns. Then we quantify these folds through various quantitative metrics, such as the mean curvatures [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF], the surface-based threedimensional gyrification index [START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF] and the sulcal depth. Besides, we introduce a novel approach to measure the anisotropy of the folding orientation, through geometric tools [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF][START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF] and the Kullback-Leibler divergence.

Specifically, this chapter attempts to answer the following questions: 1) What is the effect of the mesh density on the accuracy of surface morphology? 2) What is the impact of the temporal cortical growth model onto the folding patterns? 3) Does regional differential growth have a mutual influence on folding patterns? 4) What is the influence of the initial cortical thickness on the folding patterns of the brain? 5) Is there a relationship between the folding complexity (as measured by the average of the absolute values of mean curvatures, the surface-based threedimensional gyrification index and the sulcal depth) and the shape of the brain (the initial geometry)? 6) Does the orientation of the folds depend on the shape of the brain?

Biomechanical model of brain folding

Tallinen et al. proposed a human cortical folding model that can mimic a realistic brain folding process [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Brain growth is modeled by a relative tangential expansion of the cortical layer and the white matter layer, the cortical layer is assumed to grow more rapidly than the white matter layer. Based on the nonlinear stress-strain property of the human brain [START_REF] Zbigniew | Why does the human brain need to be a nonlinear system?[END_REF] and the bulk modulus that is assumed to be five times the shear modulus [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the human brain was brought into the frame of a modestly compressible Neo-Hookean material solid. The model is guided by the use of 3D MRI of a smooth fetal brain as an initial point. The relative tangential expansion mechanism and a simulation starting from the smooth fetal brain are shown in Figure 2.1. This model is based on an explicit dynamic solver for quasi-static equilibrium of the system and allows the simulation of the large strains and highly nonlinear mechanics involved in gyrification, but the brain solid should be discretized into high-density tetrahedral finite elements [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF]. The time step is defined as dt = 0.05a ρ/K to avoid computational instabilities [START_REF] Belytschko | Nonlinear finite elements for continua and structures[END_REF][START_REF] Thomas | The finite element method: linear static and dynamic finite element analysis[END_REF], where a is the 2.2. Biomechanical model of brain folding mesh spacing which should be set manually based on the average spacing in the mesh, ρ is the mass density and K is the bulk modulus.

Figure 2.1: Brain is modelled as a soft elastic solid and a relative tangential expansion is imposed on the cortical layer as shown at left. The simulation starting from a smooth fetal brain shows gyrification as a result of uniform tangential expansion of the cortical layer. Extracted from [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

The model uses free boundary conditions and initial conditions derived from the 3D MRI of a 22 weeks' fetal brain, with the initial displacement and the initial velocity which are zero. Two main forces are considered in this model. One is the elastic force, which is derived from the volumetric strain energy density of neo-Hookean and a deformation gradient. Another is the contact force, which takes place when a separation between a node and a triangle face at the brain surface is less than a threshold in order to prevent nodes from penetrating element faces.

The elastic traction of each deformed face of a tetrahedron is given by

s i = -σn i , (2.1)
where n i is the normal with the length proportional to the deformed areas of each tetrahedral face (i = 1, 2, 3, 4), σ is the Cauchy stress, i.e., the force per unit area in the deformed configuration, which is derived from the volumetric strain energy density W by

σ = 1 J ∂W ∂F F T (2.2)
The strain energy density functions establish the relation between the strain energy density of a material to the deformation gradient, and they usually separate the energy generated by isochoric deformations (shape change) and volumetric deformations (size change). The relation between the deformation and energy for each component can be determined and related to shear modulus (µ) and bulk modulus (K) for isochoric and volumetric deformations, respectively. For a compressible Neo-Hookean material the strain energy density W is defined as

W = µ 2 [T r (F F T )J - 2
3 -3] + K 2 (J -1) 2 , (2.3)
where K is assumed to be 5µ, corresponding to a modestly compressible material. F is the deformation gradient, J is equal to det(F ). Note that T r (F F T )J -2

3 -3 is greater than zero for shape changing ("shear") deformations but not changes in volume. On the contrary, (J -1) 2 is greater than zero for size changing ("bulk") deformations but not changes in shape. The deformation gradient F is defined in this model by

F = A(G Â) -1 , (2.4)
where A is the deformed configuration of a tetrahedron. Â is the stress-free initial configuration of a tetrahedron. It differs from the traditional definition of F = A Â-1 by integrating the relative tangential growth tensor G. The relative tangential growth tensor G, which describes the tangential expansion perpendicular to the normal vector n of the tetrahedron, is calculated by

G = gI + (1 -g)n ⊗ n, (2.5)
where g is the relative tangential expansion ratio of the grey matter to the white matter, which associates with the distance of a tetrahedron from surface in material coordinates (shown in Figure 2.1) and is given by the relation

g = 1 + α t 1 + e 10( y H i -1) , ( 2.6)
where α t controls the magnitude of expansion and varies linearly with t, t parametrizes time of model and has a non-linear relation to gestational age (GA) as t = 6.926 × 10 -5 GA 3 -0.00665GA 2 +0.250GA-3.0189 [START_REF] Wang | On early brain folding patterns using biomechanical growth modeling[END_REF], t ∈ [0, 1] corresponds to the gestational age ∈ [22weeks, adult], y is the distance from the top surface, which is calculated for four vertices of each tetrahedron and will be averaged, H i is the initial cortical thickness.

The contact force at the brain surface is obtained via penalty based vertex-triangle contact processing [START_REF] Ericson | Real-time collision detection[END_REF]. If a node is close enough to a triangle face at the surface, we should consider that there is a contact force f c for the node and f c for the face [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]:

f c = k c a 2 ( d -h h ) d, (2.7)
where k c is the contact stiffness defined as 10K in the model, a is the initial mesh spacing, d is the distance between the node and the face, and h is the contact offset which is manually set as a 5 .

The dynamic scheme in this model is derived from Newton dynamic and described as

v t+dt = v t + f t -γv t V n m dt, ( 2.8)
where γv t V n is the damping force (γ is the viscous damping, v t is the velocity, V n is the mesh volume calculated based on the configuration of the growing tetrahedron (G Â)), f t is the sum of the elastic force and the contact force, dt is the time step. m is equal to V n ρ, where ρ is the mass density. The position is calculated by

x t+dt = x t + v t+dt dt, (2.9)
x t+dt is used to calculate the matrix A of a deformed tetrahedron.

When the relative tangential growth tensor G is initialized, the brain solid starts to grow and the deformation gradient is formed, the corresponding elastic force can be calculated. The resultant force (the sum of the elastic force, the contact force and the damping force) is applied as the nodal force on each node of the mesh to produce the deformation of the brain solid.

Biophysical and numerical parameters

The biophysical parameters defined in the model, such as the initial geometry, the cortical growth (α t ) and the initial cortical thickness (H i ), may affect the surface morphology of soft solids.

To study the relation between these biophysical parameters and the folding patterns, we first investigate the effect of the mesh density on folding accuracy, then based on an appropriate mesh density, we study the impact of the cortical growth (α t ), the initial cortical thickness (H i ) and the initial geometry onto soft solids morphology. For the definition of the initial cortical thickness (H i ), we should consider the size of the initial geometry. Our study is based on an ellipsoid (hereinafter referred to as the reference ellipsoid). The two equatorial radius and polar radius of the reference ellipsoid are approximately 10, 9 and 7 mm respectively, thus the longitudinal length (LL) of this ellipsoid is 20 mm. For a normal fetal brain at 22 weeks, the brain longitudinal length (BLL) is approximately 60 mm [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF][START_REF] Kyriakopoulou | Normative biometry of the fetal brain using magnetic resonance imaging[END_REF], and the typical cortical thickness is 2.5 mm [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. In order to respect the ratio of the initial cortical thickness to the longitudinal length, which is described by

H brain BLL = H ellipsoid LL ellipsoid , ( 2.10)
thus the initial cortical thickness for the reference ellipsoid is 0.83 mm. It should be noted that the model has a coordinates normalization part (the three-dimensional coordinates will be ∈ [-1, 1]), thus the initial volume of the solid does not affect simulation results, and the initial cortical thickness will be normalized to 0.042. For the cortical growth (α t ), the cortical layer usually has an areal growth by a factor of g 2 = 8 (in Equation 4.4) [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], thus the linear cortical growth is originally defined as α t = ( √ 8 -1)t = 1.829t in this model.

Mesh density

Mesh density is the number of elements per unit volume in a volumetric mesh. In finite element analysis, mesh density is a crucial issue, which is closely associated with the accuracy of the finite element model and determines its complexity degree. In order to investigate the effects of variations in mesh density on surface morphology, based on the same surface mesh of the reference ellipsoid (its lengths in three dimensions are 20, 18 and 14 mm), we generate tetrahedral meshes of approximately 10 2 to 10 6 tetrahedra by using Netgen (an automatic 3d tetrahedral mesh generator). When the number of tetrahedral elements is relative small (10 2 tetrahedra), we increase the number of tetrahedra at an eight-time rate, and then the tetrahedra become more, and we decrease its increase multiple. They have 535, 4280, 34240, 200944, 1181216, 2314240, 3897088 and 5248576 tetrahedra, respectively. The normalized volume of these volumetric meshes is 2.5 cm 3 , therefore, the mesh densities of them are approximately 214, 1712, 13696, 80378, 472486, 925696, 1558835 and 2099430 tetrahedra/cm 3 .

Cortical growth

The cortical growth (defined by α t in the model) has an effect on surface morphology [START_REF] Bayly | Mechanical forces in cerebral cortical folding: a review of measurements and models[END_REF][START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF]. Different growth models, which are used to describe different change tendencies of the cortical growth, may also have an impact onto the folding patterns. To better understand it, we first use a linear growth model, which was initially defined in the brain folding model as α t = ( √ 8-1)t [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Secondly, considering that the Gompertz distribution can well represent the growth of human brains [START_REF] Peterson | Normative human brain volume growth[END_REF], thus we choose the Gompertz growth model to compare with the linear growth model, which is defined as:

α t = ae -e -b(t-c) , (2.11)
where a is the asymptotic value, b sets the growth rate and c sets the displacement along the t-axis. Since the purpose is to explore the effect of change tendency of the cortical growth, the initial and final values of the Gompertz growth model should be in agreement with those of the linear growth model, thus we assume a is a constant parameter ( √ 8 -1). With (b, c) equal to (6.6, 0.43) and (7.5, 0.19), we define two growing modes which correspond to the 1 st and 2 nd Gompertz models shown in Figure 2.2. In addition, another aim is to know how the extreme growth mode of non-growth for a long time and rapid growth at the end of the simulation to reach the same final growth will affect the folding patterns, thus the logistic model is adopted:

α t = a 1 + e -b(t-c) , (2.12)
where a is curve's maximum value, b is the logistic growth rate and c is the x value of the sigmoid's midpoint. With a = √ 8 -1, b = 50 and c = 0.9, the model allows to grow rapidly around time 0.9. Eventually, these growth models are integrated into the brain folding model, respectively. Based on each combined model, the simulation is performed on the reference ellipsoid with the mesh density 10 6 tetrahedra/cm 3 and the initial cortical thickness 0.83 mm. Furthermore, the regional differential growth may also affect the folding patterns between regions. To study the influence of regional growth, a simulation of the reference ellipsoid is executed based on the linear expansion model for the right semi-ellipsoid and the 1 st Gompertz expansion model for the left semi-ellipsoid.

Initial geometry

The pattern and location of folds can be influenced by initial geometry [START_REF] Bohi | Global perturbation of initial geometry in a biomechanical model of cortical morphogenesis[END_REF][START_REF] Todd | A geometric model for the cortical folding pattern of simple folded brains[END_REF][START_REF] Toro | On the possible shapes of the brain[END_REF]. For example, in ellipsoid models, most folds run either parallel or orthogonal to the ellipsoid's long axis [START_REF] Todd | A geometric model for the cortical folding pattern of simple folded brains[END_REF][START_REF] Toro | On the possible shapes of the brain[END_REF]. In order to understand the impact of the initial geometry more clearly, an affine transformation (elongated transformation) is applied to the initial geometry to determine whether the complexity of folding patterns and how the direction of folds will change. Therefore, we propose that, while keeping the volume and the y-axis length of the geometry unchanged, the reference ellipsoid is scaled in x and z directions to obtain a sphere or the ellipsoids with different elongation ratios. The elongation ratio is defined by the ratio of the x-axis length to z-axis length. The x-axis and z-axis lengths and the corresponding elongation ratios for these geometries are illustrated in Table 2.1. Based on each geometry with almost the same mesh density (10 6 tetrahedra/cm 3 ), we also vary the initial cortical thickness in the brain folding model from 0.03 to 1.63 mm to simulate the folding processes.

Furthermore, the significant difference in geometry of an ellipsoid and a human brain is the longitudinal fissure between two hemispheres. To consider the effect of the fissure and the induced contact force between the semi-ellipsoids onto the surface morphology, we also compare Length x (mm) Length y (mm) Length z (mm) Elongation ratio (x/z) the folding process of the complete ellipsoid (reference ellipsoid) and the corresponding ellipsoid with a fissure. The ellipsoid with a fissure is obtained by performing a Boolean operation on the reference ellipsoid and a torus using COMSOL (https://www.comsol.com/), as is illustrated in Figure 2.3. In addition, considering that the sharp boundaries close to the longitudinal fissure may cause changes in folding patterns, thus we smooth the boundaries of the two connected semi-ellipsoids. To achieve it, we first create a 2d quarter circle and select the vertex below to chamfer at its edge. Then we rotate it 360 degrees around its upper edge, and add a small cylinder connected in the middle. Finally, we scale it appropriately in x, y, z directions and mirror the semi-ellipsoid. The lengths in the three directions are approximately 20, 18 and 14 mm, respectively, which are almost the same as those of the reference ellipsoid. The demonstration of creating the ellipsoid with a smooth fissure is shown in Figure 2.4. The simulations are performed on these ellipsoids of almost the same mesh density (10 6 tetrahedra/cm 3 ) with the initial cortical thickness 0.83 mm and the linear growth model α t = 1.829t. To know the difference in folding between the ellipsoid with a fissure and its semi-ellipsoid and further study the effect of the contact force on surface morphology, we compare the folding processes of the ellipsoid with a fissure and its left semi-ellipsoid. The mesh density is almost the same for them, i.e., the number of tetrahedra of the ellipsoid with a fissure is almost double for that of its semi-ellipsoid. The simulations are also performed with the initial cortical thickness 0.83 mm and the linear growth model α t = 1.829t.

Initial cortical thickness

The cortical folding patterns can also be influenced by the initial cortical thickness [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF]. To understand the effect of the initial cortical thickness on surface morphology, based on ellipsoids of different elongation ratios with the linear growth model α t = 1.829t, we vary the initial cortical thickness in the brain folding model from 0.03 to 1.63 mm (0.03, 0.43, 0.63, 0.83, 1.03, 1.23 and 1.63 mm) to simulate the folding processes of ellipsoids. The cortical thicknesses from 0.43 to 1.23 mm are defined according to normative human cerebral cortex measurements [START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF] and the scale factor of the longitudinal length which is introduced in Section 2.3. The other two cortical thicknesses (0.03 and 1.63 mm) are the hypotheses for abnormal cortical thicknesses.

Quantitative methods

Curvatures on triangle meshes

The normal curvature on a 3D surface in some direction is the reciprocal of the radius of the circle that best approximates a surface normal slice in that direction [START_REF] Cipolla | Visual motion of curves and surfaces[END_REF]. The normal curvature for a smooth surface can be represented by the Weingarten matrix, i.e., the second fundamental tensor II, which is defined in terms of the directional derivatives of the surface normal:

II = D u n D v n = ∂n ∂u • u ∂n ∂v • u ∂n ∂u • v ∂n ∂v • v , ( 2.13)
where (u,v) are the directions of an orthonormal coordinate system in the tangent frame (the sign convention used here produces positive curvatures for convex surfaces with outward-facing normals). Multiplying this tensor by any vector in the tangent plane can get the derivative of the normal in that direction:

IIs = D s n. (2.14)
In our study, we compute the curvature based on Rusinkiewicz estimation [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF], which may be thought of as an extension of common methods, such as the curvature presented in Knutsen et al. [START_REF] Andrew K Knutsen | Spatial and temporal variations of cortical growth during gyrogenesis in the developing ferret brain[END_REF] and used by subsequent authors like Garcia et al. [START_REF] Garcia | Mechanics of cortical folding: stress, growth and stability[END_REF], for the purpose of estimating per-vertex normals by averaging adjacent per-face normals. This algorithm uses the "Voronoi area" weighting which can produce more accurate normal estimates of curvature than other weighting methods for triangles of varying sizes and shapes. In this algorithm, the per-face (per-triangle) curvature tensor is first computed by its three well-defined directions (the edges) together with the differences in normals in those directions (computed from the per-vertex normals). Then the algorithm performs a coordinate system transformation, because the curvature tensor represented in the coordinate system of a face should be averaged with the contributions from adjacent triangles. To this aim, they suppose that each vertex has its own orthonormal coordinate system in a plane perpendicular to its normal, and derive a coordinate change formula for converting the curvature tensor to the vertex coordinate frame. For the question of how much of the face curvature should be accumulated for each vertex, they take the weight w f,p to be the area of the face f divided by the squares of the lengths of the two edges that touch vertex p. They found that this produces more accurate normal estimates than other weighting methods and is exact for vertices which lie on a sphere.

Mean curvature of a vertex is defined by the average of the two principal curvatures (the maximal and minimal curvatures) of the vertex, and the principal curvatures are the eigenvalues of the vertex normal curvature tensor computed by Rusinkiewicz estimation:

II = u v K 1 0 0 K 2 u v , ( 2.15)
where K 1 and K 2 are the eigenvalues and (u , v ) are the principal directions, which are the directions in which the normal curvature reaches its minimum and maximum. Since surface curvature is useful to describe spatial variations in folding, thus for the overall folding complexity comparison, we compute the average (across all vertices on the mesh surface) of the absolute values of mean curvatures (at a vertex on the mesh surface) for each simulated surface. In the remainder of this chapter, we simply use the term curvature for the sake of clarity.

Three-dimensional gyrification index

Curvature-based features do not provide a complete description of the folding patterns. In order to describe globally the folding complexity by considering the depth and wideness of the cortical folding, we also use the surface-based three-dimensional gyrification index (3D GI). It is a global measurement which is defined as the ratio of the cortical surface area to the area of its smooth "convex hull" (the minimum surface area needed to completely enclose the brain) [START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF]:

3D GI = area of cortical surf ace area of convex hull .

(2.16)

To get the convex hull, we scale the initial smooth surface in three dimensions so that the three-dimensional lengths of the convex hull are equal to those of the simulated cortical surface.

Sulcal depth

Sulcal depth can be used as a quantitative marker of cortical morphology [START_REF] Shimony | Comparison of cortical folding measures for evaluation of developing human brain[END_REF] and the variation in sulcal depth is also considered to be related to some neurodevelopmental diseases, such as sulcal depth reduction in most cortical lobes in schizophrenia compared to healthy controls [START_REF] Lyu | Sulcal depthbased cortical shape analysis in normal healthy control and schizophrenia groups[END_REF].

Several approaches have been proposed to compute the sulcal depth [START_REF] Boucher | Depth potential function for folding pattern representation, registration and analysis[END_REF][START_REF] Im | Fractal dimension in human cortical surface: multiple regression analysis with cortical thickness, sulcal depth, and folding area[END_REF][START_REF] Hyuk | Automated sulcal depth measurement on cortical surface reflecting geometrical properties of sulci[END_REF] but a well-defined computation of depth remains an open question. In this work, we make use of an intuitive approach to calculate the sulcal depth by using the distance between the deformed mesh surface and the corresponding convex hull. Specifically, for each surface vertex of the deformed mesh, we find the intersection point on the convex hull by using the vector determined by the corresponding vertex of the initial mesh and this vertex and the method of traversing all triangles on the convex hull. Then we compute the distance between each surface vertex of the deformed mesh and its corresponding intersection point on the convex hull.

Folds orientation

For the purpose of describing and comparing the direction of the folds on the simulated surfaces, we calculate the angle between the gradient of Fiedler vectors [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF] and the principal directions of curvatures [START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF], which helps to understand whether the folds are isotropic. The Fiedler vector is the first non-constant eigenfunction of Laplace-Beltrami operator, represented by φ 1 in Equation 2.18 [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF]. The Laplace-Beltrami operator is defined as

∆ M = div • ∇ M = 1 det(g) i,j ∂ x j ( det(g)g i,j ∂ x i M ), (2.17)
where (M, g) is a compact 2-Riemannian manifold without boundary, det(g) := |det(g i,j )| is the absolute value of the determinant of the metric tensor, g i,j are the components of the inverse of the metric tensor. The (x i , x j ) are the coordinates of a local coordinate system x(p) = (x 1 , x 2 ), i.e., a local diffeomorphism M → R 2 around a point p.

The eigenvalues of -∆ M are λ 0 = 0 ≤ λ 1 ≤ . . ., and φ 0 , φ 1 , . . . are corresponding eigenfunctions, which satisfy

-∆ M φ i = λ i φ i .
(2.18)

The eigenfunctions are orthogonal in the sense of the scalar product < u, v > M = M uvdµ, where the volume form dµ is given by det(g)dx 1 dx 2 . Figure 2.5 gives a visual illustration of φ 1 , φ 2 , φ 3 .

Figure 2.5: Eigenfunctions φ 1 , φ 2 and φ 3 (from left to right). Colormap: blue(negative), red/yellow(positive). Extracted from [START_REF] Lefèvre | Spherical parameterization for genus zero surfaces using laplace-beltrami eigenfunctions[END_REF].

The Fiedler vector allows to describe the longitudinal extension of surfaces [START_REF] Lai | Intra-patient supine-prone colon registration in ct colonography using shape spectrum[END_REF][START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF][START_REF] Levy | Laplace-beltrami eigenfunctions towards an algorithm that" understands" geometry[END_REF][START_REF] Shi | Inverseconsistent surface mapping with laplace-beltrami eigen-features[END_REF]. The Fiedler's extrema are the most distant points, and its contour lines are slices in the elongation axis. The gradient of the Fiedler vector that is orthogonal to the contour lines gives the direction of elongation. The principal directions of curvatures are the corresponding eigenvectors of the principal curvatures (the eigenvalues of the Weingarten matrix). Based on the local scalar product between the gradient of the Fiedler vector and the principal directions of curvatures, we can obtain the fold angles.

In order to compare quantitatively the uniformity of the angular distribution of folds, we use the Kullback-Leibler (KL) divergence. The KL divergence, also called relative entropy, is used to measure how one probability distribution is different from a second reference probability distribution. For two discrete probability distributions P and Q defined on the same probability space, the KL divergence from P to Q is defined to be

D KL (P ||Q) = i P (i) log P (i) Q(i) . ( 2

.19)

For angular uniformity calculations, P corresponds to the fold angular distribution on the folded surface, and Q represents the theoretically uniform distribution of fold angles.

Results

Results

Mesh density

With the initial cortical thickness 0.83 mm and the linear cortical growth model α t = 1.829t, the simulation results of the brain folding model are shown in Figure 2.6. It can be observed that deformations first appear where the sizes of elements of entire ellipsoids are the smallest, i.e., the boundary areas of ellipsoids. Moreover, we note that the higher the mesh density is, the greater the number of folds becomes and the smaller the width of gyri seems to be. However, when the mesh density reaches approximately the order of 10 6 (925696 tetrahedra/cm 3 ), the size and the number of folds will rarely change anymore. That is to say, when the mesh density reaches a certain order of magnitude, further increases in mesh density will increase computational cost but cannot significantly change the spatial frequency of folding patterns. However, the folding patterns are different for the last three largest mesh densities, which is because the tetrahedral meshes of different densities that we use are different. These different perturbations in mesh can produce different patterns because they are the mechanism breaking the symmetry in the system [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF]. In short, when the density of the mesh reaches a certain order of magnitude, the size of the folds tends to be stable, but the folding patterns have random spatial variations according to each mesh.

The comparison of the curvature for surfaces of different mesh densities is shown in Figure 2.7 (left figure shows the dimensional curvature, right figure shows the dimensionless curvature). We can observe that, no matter we use the dimensional curvature or dimensionless curvature, the tendency does not change for different mesh densities. The difference is that the dimensional curvature tends to no longer increase and even slightly decreases after time 0.5, because the dimensional mean curvature based on Rusinkiewicz approach [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF] is inversely proportional to the surface area. Thus when the folding degree no longer becomes more complex, the dimensional mean curvature decreases with the growing surface area, resulting in a decrease in the average curvature. The dimensionless mean curvature is computed by multiplying the dimensional mean curvature by the square root of the surface area. In the following studies of this chapter, we only compute the dimensional curvature.

It can be seen that after time 0.3, except for the mesh densities of 925696 and 1558835, the higher the mesh density is, the greater the curvature becomes. For the last three largest mesh densities, the change in curvature is enough small (< 0.04), the error is less than 5% (0.04/0.9 ≈ 4.4%), which can be ignored. Therefore, we can understand it as when the mesh density reaches the order of 10 6 tetrahedra/cm 3 , the solution converges and the complexity degree of folding patterns no longer changes, the curvature is around 0.9.

The comparison of the 3D GI computed on surfaces of different mesh densities is shown in Figure 2.8. The higher the mesh density is, the greater the 3D GI becomes, but for meshes with densities greater than the order of 10 6 tetrahedra/cm 3 , the difference in the 3D GI is very small (< 0.1) and the tendence is not evident. Considering both the curvature and the 3D GI, the mesh with the density of 10 6 tetrahedra/cm 3 can already achieve sufficient folding accuracy using the brain folding finite element model, thus we use this mesh density to study the impact of the physical parameters onto the folding patterns.

Cortical growth

Global cortical growth

Based on the simulation results of the reference ellipsoid with the same initial cortical thickness 0.83 mm but different growth models, we first compute the 3D GI over time, as shown in Figure 2.9a. The increasing tendencies of the 3D GI for the different growth models are consistent with the curves of the cortical growth shown in Figure 2.2. The final 3D GI is almost the same for the models, showing that different growth models with the same initial and final growth will almost not affect the complexity of the final surface morphology.

For the 3D GI of 1.0, 1.4, 2.3 and 2.8, the surface morphology of different growth models is shown in Figure 2.9b. We can observe that the primary folds of different models appear almost at the same positions. However, the divergence of folding patterns takes place when the 3D GI reaches 2.3 for the linear model and the 1 st Gompertz model (in red frames), and it occurs when the 3D GI is 2.8 for the 1 st and 2 nd Gompertz models (in green frames), indicating that the same type of expansion model can make the difference in folds occur relatively late. For the logistic expansion model, the folds are relatively shorter compared to those of other models since the cortical growth is zero for a long time and suddenly increases close to the end. The final results (at time 1.0) of different growth models are shown in Figure 2.9c. The amount of the final folds is almost the same for these growth models, but the patterns of the folds are visually different. The impact of growth model onto the curvature is reported in Figure 2.10. For the 2 nd Gompertz growth model, the curvature increases and reaches its maximum earlier than that of other models. For the logistic model, the curvature suddenly increases towards the end of the simulation. The final value of the curvature (at time of model 1.0) is almost the same for these growth models. This measurement further confirms that different growth models with the same initial and final growth will almost not affect the complexity of the final surface morphology.

Furthermore, in order to quantify the correlation of the folding patterns generated by different growth models under the same surface folding complexity, therefore, at the same 3D GI, we calculate the Pearson correlation coefficient of the mean curvatures of all vertices on the surfaces produced by every two different expansion models. In addition, we also compute the Pearson

Results

correlation coefficient of the distances from all vertices to the convex hull between every two different expansion models. The results are shown in Figures 2.11a and 2.12. For the 3D GI of 2.3, the mean curvature of each vertex on surfaces is displayed on the undeformed surface, as shown in Figure 2.11b. The correlation coefficients of mean curvatures and distances to convex hull between different growth models are all higher than 0.6 and 0.7, respectively. Especially for the linear and 1 st Gompertz models, and the 1 st and 2 nd Gompertz models, the correlation coefficients of mean curvatures and distances are as high as 0.8 and 0.87 after all folds are formed. It can also be seen that before the 3D GI reaches 2.45, the correlation coefficients of mean curvatures and distances to convex hulls for the 1 st and 2 nd Gompertz growth models are higher than those of the others, which is consistent with the similar folding patterns of the two Gompertz models shown in Figure 2.9.

Regional cortical growth

The simulation of regional differential growth mode is shown in Figure 2.13. At the same time, when we compare each region of the surfaces of the regional differential growth mode (3 rd row in Figure 2.13) with two other global growth modes (1 st and 2 nd rows), it can be observed that from time 0.55, the regional differential growth mode leads to a mutual influence in the pattern of the folds especially near the contact area of the two semi-ellipsoids, which may be due to the discontinuous growth. In addition, the curvature and the 3D GI are computed on left and right semi-ellipsoids of the global and regional linear and 1 st Gompertz growth models, the comparisons are shown in Figure 2.14. In terms of regional growth, using the same growth model for the entire ellipsoid and half of the ellipsoid, the curvatures of the corresponding regions are almost the same over time, which indicates that the different regional growths do not affect each other the complexity of the folding patterns. For the 3D GI, for the regions of the same growth model (e.g., linear_overall and linear_half), the 3D GI values are almost the same over time, indicating that for a soft 2.5. Results tissue, the different regional growths don't have a mutual effect on the complexity degree of the folds.

Figure 2.14: Left: the comparison of the curvatures for the linear and 1 st Gompertz growth regions. Right: the comparison of the 3D GI for the linear and 1 st Gompertz growth regions. "linear_overall" corresponds to the right semi-ellipsoid of 1 st row (in Figure 2.13) which is applied with global linear model, "gompertz_overall" corresponds to the left semi-ellipsoid of 2 nd row which is applied with global 1 st Gompertz model, "linear_half" corresponds to the right semi-ellipsoid of 3 rd row which is applied with regional linear model and "gompertz_half" corresponds to the left semi-ellipsoid of 3 rd row which is applied with regional 1 st Gompertz model. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

Initial geometry and cortical thickness

Ellipsoids with different elongation ratios and cortical thicknesses

The simulation results for the geometries of the elongation ratios 1.0, 1.50 and 2.25, with the initial cortical thickness 0.83 mm and the cortical growth α t = 1.829t, are shown in Figure 2. [START_REF] Belytschko | Nonlinear finite elements for continua and structures[END_REF]. In terms of folding patterns, we can observe that after time 0.55 when most of the folds have already formed, the size and spatial frequency of the folds are almost the same for the three geometries. However, their folding patterns are different.

The simulated surfaces of the reference ellipsoid with the initial cortical thickness varying from 0.03 to 1.63 mm are shown in Figure 2.16. The surface with the thinnest initial cortical thickness 0.03 mm folds relatively late than the others. In addition to the surface of the thinnest initial cortical thickness, at time 0.32, other surfaces already showed some clear primary folds, and the thicker the cortex is, the more obvious the folds are; starting from time 0.55, most of gyri and sulci are formed for all of these thicknesses. As the initial cortical thickness increases, the gyri become larger and the folds become fewer, thus the adjacent sulci are more isolated; after time 0.55, the folds are almost no longer complex, but the volumes are still growing especially for the thicker initial cortical thicknesses. The final surface morphology of the thinnest cortical thickness (0.03 mm) resembles the folding patterns of polymicrogyria [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Budday | The role of mechanics during brain development[END_REF], while the surface morphology of the thickest cortical thickness (1.63 mm) is similar to the phenomenon of pachygyria.

For all ellipsoids of different elongation ratios and initial cortical thicknesses, the curvature is computed on each surface and the comparison is shown in Figure 2.17. It can be observed that the curvature does almost not depend on the elongation ratio (i.e., the initial geometry), but it depends on the initial cortical thickness. For the initial cortical thicknesses between 0.43 and 1.63 mm, the thinner the initial cortical thickness is, the more quickly the curvature increases

Results

after time 0.32, and the greater the curvature becomes. The decrease in the curvature after time 0.55 is due to the larger surface area increase relative to the degree of folding. For the thinnest initial cortical thickness 0.03 mm, the curvature is smaller than that of the others before time 0.55; Starting from time 0.55, it increases faster, and eventually becomes even greater than the curvature corresponding to the thicker initial cortical thicknesses.

To further verify the degree of folding, we also compute the surface-based three-dimensional gyrification index (3D GI) [START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF] on each surface. The comparison is shown in Figure 2.18. The 3D GI does not depend on the elongation ratio, but it depends on the initial cortical thickness. For the initial cortical thicknesses between 0.43 and 1.63 mm, the increment of 3D GI becomes smaller as the initial cortical thickness increases. When the initial cortical thickness is overly thin (0.03 mm), the 3D GI is smaller at the beginning due to the late generation of folds and the shallow primary folds; from approximately time 0.55, the increase in 3D GI begins to become faster but the increment of 3D GI is less than that of the other cortical thicknesses in the end. These quantitative measurements demonstrate that when the initial cortical thickness within a reasonable range, the thinner the initial cortical thickness is, the more complex the folding patterns are, which confirms the previously reported effects of cortical thickness [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Wang | On early brain folding patterns using biomechanical growth modeling[END_REF].

The folding process can also be quantified by studying the sulcal depth. The histograms of sulcal depth at different time with respect to different initial geometries and cortical thicknesses are shown in Figures 2.19 and 2.20. It can be seen that the geometry has an effect on sulcal depth. The simulated sulci on the ellipsoids are deeper than in the case of the sphere, the greater the elongation ratio is, the deeper the sulci become, and this tendency becomes more obvious over time. Combining the above results, we can conclude that the elongation ratio of the initial geometry does almost not change the surface curvature and 3D GI, but it has an impact on sulcal depth. Moreover, the thinner the cortex is, the shallower the sulci become, and this tendency also becomes more obvious over time. This measurement result is in agreement with the analysis of sulcation morphology in polymicrogyria [START_REF] James | Mri analysis of sulcation morphology in polymicrogyria[END_REF]. The visual depth maps at time 1.0 for different initial cortical thicknesses are shown on the deformed reference ellipsoids and corresponding initial ellipsoid in Figure 2.21. As the initial cortical thickness increases, the sulci become deeper and wider.

Furthermore, to understand the effect of the initial geometry on the folds orientation, we calculate the angle between the gradient of Fiedler vectors [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF] and the principal direction of curvatures [START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF] on the surface of each geometry. The angles distributions for the geometries of the elongation ratio 1.0, 1.50 and 2.25 with the initial cortical thickness 0.83 mm at time 0.32, 0.55 and 0.79 are shown in Figure 2.22. We can observe that at time 0.32, the distribution of the fold angles is almost uniform on the sphere and slightly nonuniform on the ellipsoid of the elongation ratio 1.50, while following a privileged direction on the ellipsoid of the elongation ratio 2.25. Since the direction of curvatures is perpendicular to the extension of folds, and the gradient of Fiedler Vector is along the longitudinal extension of surfaces, thus the peak appearing at around 90 • indicates the most of folds are along the longitudinal axis of the ellipsoid of the elongation ratio 2.25, which is consistent with previous observations of the location of folds in ellipsoid models [START_REF] Todd | A geometric model for the cortical folding pattern of simple folded brains[END_REF][START_REF] Toro | On the possible shapes of the brain[END_REF]. As time passes, the number of folds increases and the angular distribution becomes more and more uniform, especially for the ellipsoid of the elongation ratio 2.25, but it is still not as uniform as that of the sphere and the ellipsoid of the elongation ratio 1.50.

To compare quantitatively the angular uniformity degree of the folding patterns for the surfaces of different elongation ratios and cortical thicknesses, the KL divergence is computed on each surface, and the results are shown in Figure 2.23. The angular distribution of folds strongly depends on the elongation ratio, and the angular distribution of folds becomes nonuniform as

Results

the elongation ratio increases, the thicker the cortex is, the more obvious this tendency becomes. As time goes on, the fold angles become more uniform especially for the geometries with larger elongation ratios.

Ellipsoid with a fissure

With the cortical growth α t = 1.829t and the initial cortical thickness 0.83 mm, the simulation results are shown in Figure 2.24. In terms of folding patterns, we can observe that the size and spatial frequency of folds are very similar for the complete ellipsoid (reference ellipsoid) and the corresponding ellipsoid with a fissure, but the shape of folds are different, indicating that the difference in geometry does not change the size and amount of folds but the induced contact force can change the shape of folds. At time 0.27, the mechanical instability induced by the relative growth of the complete ellipsoid first appear at the boundary; however, for the ellipsoid with a fissure, the two semi-ellipsoids begin to contact each other, so that the folds are formed both at the boundary and in the contact area.

For the quantitative analysis, in order to consider only the complexity of folds in areas other than the longitudinal fissure, the curvature is calculated in the pink area on each surface, as shown in Figure 2.25. The comparison of the curvature between the two geometries is shown in Figure 2.25. It can be observed that, after approximately time 0.4 when most of the folds begin to form on their surfaces, the curvature of the ellipsoid with a fissure is slightly larger than that of the complete ellipsoid, showing that the folding patterns on the ellipsoid that with a fissure are slightly more complex than those on the complete ellipsoid. The comparison of the 3D GI between the complete ellipsoid and the ellipsoid with a fissure is shown in Figure 2.26. It can be observed that after around time 0.4, the difference in 3D GI gradually increases, and the 3D GI of the ellipsoid with a fissure is larger than that of the complete ellipsoid, which is consistent with the result of the curvature shown in Figure 2.25.

To understand the effect of the fissure on the folds orientation, we compute the fold angles in the pink area (in Figure 2.25) on the simulated surfaces, and the comparison of it between the two geometries is shown in Figure 2. [START_REF] Chada | Cytomechanics of neurite outgrowth from chick brain neurons[END_REF]

H i = 0.042.
of the angular distribution of folds for these two geometries, we compute the Kullback-Leibler (KL) divergence for the two pink areas, the results are shown in Table 2.2. It can be observed that the primary folds distribution on the complete ellipsoid is more uniform than that on the ellipsoid with a fissure, the fissure results in the primary folds formed in a privileged direction.

As the number of folds increases, the angular distribution becomes more and more uniform for the ellipsoid with a fissure and it tends to be almost the same for these two geometries. This experiment shows that a fissure on the soft solid and the induced contact between two hemispheres can lead to more complex folding patterns and result in nonuniform primary folds distribution.

Semi-ellipsoid

The comparison of surfaces is shown in Figure 2.28. We can see that the contact between two semi-ellipsoids still causes some changes in shape of folds, especially in the area where the two semi-ellipsoids are in contact with each other, compared with the left semi-ellipsoid, the folds of the two connected semi-ellipsoids seem smaller.

Figure 2.28: Comparison of surface morphology for the ellipsoid composed of two connected semi-ellipsoids (above) and the corresponding left semi-ellipsoid (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

For the quantitative analysis, to consider only the complexity of folds in the same areas except the longitudinal fissure, the curvature is calculated in the pink area on each surface, as shown in Figure 2.29. The comparison of the curvature between the two geometries is shown in Figure 2.29. It can be observed that, after time 0.4 when most of the folds begin to form on their surfaces, the curvature of the ellipsoid with a fissure is larger than that of the semi-ellipsoid, showing that the folds on the ellipsoid with a fissure are slightly more complex than those on the semi-ellipsoid.

The comparison of the 3D GI between the ellipsoid with a fissure and the semi-ellipsoid is shown in Figure 2.30. It can be observed that after around time 0.35, the 3D GI of the ellipsoid with a fissure is larger than that of the semi-ellipsoid, which is consistent with the result of the curvature shown in Figure 2.29.

Therefore, this experiment further confirms that the contact force between two hemispheres does affect the folds and can lead to more complex folding patterns.

Ellipsoid with a smooth fissure

With the cortical growth α t = 1.829t and the initial cortical thickness 0.83 mm, the simulation results are shown in Figure 2.31. In terms of folding patterns, we can observe that at time 0.27, more sulci appear on the surface of the ellipsoid which has a smooth fissure, the spatial frequency of folds seems to be slightly higher than that of the complete ellipsoid. But over time, the spatial frequency of folds becomes similar for both of them. Figure 2.31: Comparison of surface morphology for the complete ellipsoid (above) and the ellipsoid with a smooth fissure (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

For the quantitative analysis, in order to consider only the areas where the initial geometries are the same for the complete ellipsoid and the ellipsoid with a smooth fissure, thus the curvature is calculated in the pink area on each surface, as shown in Figure 2.32. The comparison of the curvature between the two geometries is shown in Figure 2.32. It can be observed that, after time 0.4 when most of the folds begin to form, the curvature of the ellipsoid with a smooth fissure is larger than that of the complete ellipsoid, indicating that the folds on the ellipsoid that with a smooth fissure are also slightly more complex than those on the complete ellipsoid.

The comparison of the 3D GI between the complete ellipsoid and the ellipsoid with a smooth fissure is shown in Figure 2.33. It can be observed that after around time 0.4, the 3D GI of the ellipsoid with a smooth fissure is larger than that of the complete ellipsoid, which is consistent with the result of the curvature shown in Figure 2.32.

This experiment demonstrates that, even if the longitudinal fissure between the two semiellipsoids is sufficiently smooth, the contact between them could lead to slightly more complex folding patterns.

Through the above four experiments, it can be seen that variations in geometry will change the folds. Even if an affine transformation (elongated transformation) of the initial geometry does almost not change the complexity of the folds, it could change the distribution and direction of the folds. And when the ellipsoid is very slender, the primary folds have a privileged orientation. In addition, a fissure on the soft solid and the contact force induced by the contact at the fissure will result in anisotropic primary folds and a little more complex final folding patterns.

Conclusion

The biomechanical model based on the differential tangential growth hypothesis has been used for the realistic simulation of the early expansion and folding process of the human cerebral cortex [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Therefore, such biomechanical model can be used to investigate the relationship between biophysical parameters and severe cortical folding malformations which are thought to be associated with neurodevelopmental diseases. In this chapter, we investigate the impacts of temporal cortical growth, regional differential growth, initial cortical thickness and initial geometry onto the surface morphology in an attempt to answer the questions raised in the introduction.

Since the mesh density can influence the accuracy of surface morphology, we first investigate the effect of the mesh density on folding patterns by varying the mesh density from 10 2 to 10 6 tetrahedra/cm 3 . We prove that when the density of the mesh reaches the order of 10 6 tetrahedra/cm 3 , simulated surfaces can already achieve sufficient folding accuracy, and further increases in mesh density will increase computational cost but cannot significantly change the spatial frequency of folding patterns. This information is useful for researches using this model to study cortical folding patterns.

Regarding the impact of the temporal cortical growth model onto the folding patterns, four growth modes are defined using linear, Gompertz and logistic models to simulate the folding process. The simulation results and quantitative indices (3D GI and mean curvature) demonstrate that, when all folds are formed on the surfaces, the different growth modes with the same initial and final growth will not cause noticeable changes of the complexity degree of the folding patterns. Nevertheless, the growth mode can affect the pattern of the folds. This may be due to the growth rate (the slope of the growth curves in Figure 2.2) at certain moments being too high to fulfill the quasi-static constraint. In addition, in a quasi-static case, even if the absolute growth rate should not affect folding, the relative growth rate of the cortex to the sub-cortical regions may have an effect on folding patterns.

Furthermore, the simulation of the differential growth mode of the left and right semi-ellipsoid shows that the difference in growth has a slight mutual effect on the patterns of the folds near the contact area, but it does not affect each other the complexity degree of folding. In recent studies, Garcia et al. observed significant regional differences in growth across the cortical surface of 30 preterm infants, which are consistent with the emergence of new folds [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF], but the effect of these regional differences on folding patterns has not been quantified. Therefore, choosing a proper local growth model may be a crucial step, which can contribute to the study of the regional differential growth of cortex.

To understand the influence of the initial cortical thickness on the folding patterns of the brain, we use five normal and two abnormal cortical thicknesses varying from 0.03 to 1.63 mm to study the surface morphology. The results show that, when the initial cortical thickness varies from 0.43 to 1.63 mm, the thinner the initial cortical thickness is, the higher the spatial frequency of the folds appears to be, but the shallower the sulci become, which is consistent with the reported effects of the cortical thickness in previous works [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Wang | On early brain folding patterns using biomechanical growth modeling[END_REF]. When the initial cortical thickness is overly thin (0.03 mm), it causes the primary folds to form late and the sulci to become shallower in depth, which is in agreement with the analysis of sulcation morphology in polymicrogyria [START_REF] James | Mri analysis of sulcation morphology in polymicrogyria[END_REF]. The final surface morphology of the cortical thickness of 0.03 mm is similar to the phenomenon of polymicrogyria, while the surface morphology of the cortical thickness of 1.63 mm resembles the phenomenon of pachygyria. These observations may be crucial for exploring the causes of autism, schizophrenia and epilepsy diseases which may be related to cortical malformations [START_REF] Auzias | Atypical sulcal anatomy in young children with autism spectrum disorder[END_REF][START_REF] Heidi | Pre-frontal lobe gyrification index in schizophrenia, mental retardation and comorbid groups: an automated study[END_REF][START_REF] Cachia | Cortical folding abnormalities in schizophrenia patients with resistant auditory hallucinations[END_REF][START_REF] Antonio Y Hardan | Increased frontal cortical folding in autism: a preliminary mri study[END_REF][START_REF] Roger | Reduced cortical folding in individuals at high risk for schizophrenia: a pilot study[END_REF][START_REF] Wu Nordahl | Cortical folding abnormalities in autism revealed by surface-based morphometry[END_REF][START_REF] Pang | Malformations of cortical development[END_REF].

To study whether there is a relationship between the shape of the brain and the folding complexity, the simulations are first performed on the initial geometry of the elongation ratio from 1.0 to 2.25. we change the elongation ratio of the initial geometry from 1.0 to 2.25 while keeping the volume and the y-axis length of the geometry unchanged. The quantitative results illustrate that the variation in the elongation ratio of the geometry has an impact on the sulcal depth, but not on the surface curvature and 3D GI. The slenderer the shape is, the deeper the sulci become. In addition, we simulate the folding process of the ellipsoid with a longitudinal fissure and of the corresponding semi-ellipsoid, and find that the longitudinal fissure (whether smooth or not) between the two hemispheres and the induced contact could cause the folding patterns to be slightly more complex.

In order to investigate whether the orientation of the folds depends on the shape of the brain, we propose to calculate the fold angles between the gradient of Fiedler vectors and the principal directions of curvatures on the surfaces of geometries with different elongation ratios or with a longitudinal fissure, and then use the Kullback-Leibler divergence to measure the anisotropy of the folding orientation. We found that the elongation ratio of the geometry can predict the orientation of the folds (at least primary folds). The slenderer the initial geometry is, the greater the number of primary folds along its longitudinal direction becomes. Moreover, the longitudinal fissure between the two hemispheres and the induced contact could result in nonuniform primary folds distribution.

In addition to the elongation ratio and the longitudinal fissure, other geometric changes may also affect the folding patterns, such as the initial curvature. Therefore, it's also important to study the effect of other geometric changes on folding. The calculation method of folding orientation can be used for future works to measure the orientation of the folds on surfaces with other geometric changes.

Introduction

Although the biomechanical brain folding model can produce sulci and gyri similar to those in developing fetal brains [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the correspondence concerning the increase in the human brain volume and gestational age is not realistic in this model. In this chapter, to make the increase in the brain volume more in agreement with the biological measurement results (e.g., [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF]), we propose a brain longitudinal length (BLL) growth model, which is integrated into the brain folding model [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], to improve the human brain growth accuracy during the cerebral folding process. In addition, we also model the increase in the average cortical thickness according to the typical average cortical thickness at several ages [START_REF] Lyall | Dynamic development of regional cortical thickness and surface area in early childhood[END_REF][START_REF] Wang | Developmental topography of cortical thickness during infancy[END_REF].

Based on the improved model, we study the effect of the initial cortical thickness on cortical folding patterns. The numerical simulations are built on the tetrahedral mesh of a fetal brain at 22 weeks of gestation (provided by Tallinen et al.) with different initial cortical thicknesses. From 3.2. Modeling of brain longitudinal length and cortical thickness the simulation results, we first visually remark the difference in the folding patterns. Then, we quantify these folding patterns by the mean curvatures [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF], the surface-based three-dimensional gyrification index, the local degree of folding [START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF] and the curvature-based spectral analysis of gyrification [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF] to compare the complexity of the folding patterns.

In addition, since the initial point of the biomechanical brain folding model [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is a tetrahedral brain mesh, we also investigate how to reconstruct a human brain volume from brain segmentation and then generate a high-quality tetrahedral mesh on the reconstructed human brain volume at the end of this chapter. Specifically, this chapter attempts to answer the following questions: 1) In the model, is the relation between the increase in brain volume, the development of folding patterns and gestational age realistic? 2) What is the influence of the initial cortical thickness on the early folding patterns of the human brain? 3) How to reconstruct a human brain volume from MRI segmentation of brain tissue and generate a high-quality tetrahedral mesh on the human brain volume?

Modeling of brain longitudinal length and cortical thickness

Modeling the growth of brain longitudinal length (BLL)

In order to explain the brain volume growth during the cortical folding process, Tallinen et al. proposed a brain longitudinal length (BLL) growth model which is defined by

BLL (mm) = 59 1 -0.55t , ( 3.1)
where t parametrizes time of model and has a non-linear relation to gestational age (GA) as t = 6.926 × 10 -5 GA 3 -0.00665GA 2 + 0.250GA -3.0189, t ∈ [0, 1] corresponds to the gestational age ∈ [22weeks, adult]. According to the longitudinal length, the size of the initial brain will be scaled isotropically (in three dimensions). After measuring the brain volume of the simulation results, we found that this BLL growth model lacked brain volume growth during the folding process and could not mimic a volume growth process consistent with the biological measurement results (e.g., [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF]). Therefore, we propose a new BLL growth model. Referring to the normative fetal brain longitudinal length data in [START_REF] Kyriakopoulou | Normative biometry of the fetal brain using magnetic resonance imaging[END_REF], we find that a second-order polynomial is good enough to fit them and thus define the brain longitudinal length (BLL) as a function of gestational age (GA) from 22 to 40 weeks as

BLL (mm) = -0.067GA 2 + 7.16GA -66.05, (3.2)
with R 2 = 0.9987. The comparison of the previous and new BLL growth models is shown in Figure 3.1. The brain mesh is then scaled isotropically according to the longitudinal length of each gestational age.

Then we combine the cortical folding model (including the contact process, the elastic process and the dynamic process) and the new BLL growth model. The structure of the combined model is shown in Figure 3.2. Based on this combined model, we then study the effect of parameters on human brain folding patterns.

Modeling the increase in cortical thickness (H)

In the cortical folding model (e.g., [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF]), the relative tangential growth tensor is calculated for each tetrahedron by G = gI + (1 -g)n ⊗ n, which is perpendicular to the normal vector n of the tetrahedron. g is the expansion ratio of the grey matter relative to the white matter, which is , where α t controls the magnitude of expansion, considering that the cortical layer has an areal growth by a factor of g 2 = 8, thus α t is defined as α t = ( √ 8-1)t [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. y is the distance from the top surface, which is calculated for four vertices of each tetrahedron and will be averaged. H i is the initial cortical thickness. Considering that the cortical thickness increases slightly during the expansion process, and seems to reach a peak during puberty according to the IQ [START_REF] Shaw | Longitudinal mapping of cortical thickness and clinical outcome in children and adolescents with attention-deficit/hyperactivity disorder[END_REF], we thus also propose to model the cortical thickness H as a linear function of time of model, as follows: H = H i + βt, where β determines the magnitude of thickening. For a typical cortical thickness of 22 weeks' fetal brains, which is 2.5 mm (after normalization by BLL is 0.042), the corresponding deformed cortical thickness is approximately 3 mm (after normalization is 0.051) at about 2 years of age [START_REF] Wang | Developmental topography of cortical thickness during infancy[END_REF], which is 97% of its adult value [START_REF] Lyall | Dynamic development of regional cortical thickness and surface area in early childhood[END_REF], we thus define β = 0.051/0.97 -0.042 = 0.01 in the simulations.

Quantitative methods

In addition to the quantitative methods introduced in Chapter 2, such as the mean curvatures [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF] and the surface-based three-dimensional gyrification index, we also use two other quantitative metrics in this chapter, namely the local degree of folding [START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF] and the curvaturebased spectral analysis of gyrification [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF].

Local degree of folding

Since the three-dimensional gyrification index is a global measurement, it may be difficult to capture local information, such as the depth or wideness of a sulcus or gyrus. The quantification such as local gyrification index [START_REF] Schaer | A surface-based approach to quantify local cortical gyrification[END_REF][START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF] has the ability to provide the degree of local folding by considering the depth/wideness of sulci and gyri in regions-of-interest. The local degree of folding proposed by Toro et al. [START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF] is defined as the local surface ratio (LSR) between the pial (cortical) surface contained in a small sphere of a radius r around each vertex x and that of a disc of the same radius r:

LSR x =
surf ace in sphere(x, r) area of disc(r) .

(

Considering that the sphere has to be sufficiently large as to encompass a few folds, but small enough to make the approximation of the lissencephalic area reasonable, and the longitudinal length of the fetal brain before birth is approximately 115 mm, thus we define that the radius of the sphere is 25 mm. We measured LSR of each vertex on the simulated surfaces of the fetal brain with different initial cortical thicknesses, then the average LSR is computed over each whole surface.

Spectral analysis of gyrification

To quantify the cortical folding patterns of our simulations, we adopt the method for the spectral analysis of gyrification [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF]. This method performs a spectral decomposition of the mean curvature of the brain surface mesh based on the Laplace-Beltrami operator eigenfunctions. It allows to extend Fourier Analysis to more general domains such as graphs or surfaces and to produce power spectra. The steps summary of spectral analysis of gyrification is shown in Figure 3.3.

The Laplace-Beltrami operator is defined as ∆ M = div • ∇ M , where M is a Riemannian manifold. The eigenvalues of -∆ M are λ 0 = 0 ≤ λ 1 ≤ . . . and φ 0 , φ 1 , . . . are associated Step 2: The curvature decomposition based on the eigenfunctions gives a series of coefficients that are the characteristics of curvature in the spectral domain. This step is reversible, which allows for spectral filtering of the curvature, and thus supports the spectral segmentation of cortical folding patterns. Step 3: The production of the band power spectrum of curvature from the squared coefficients depends on the band design choices detailed in Figure 3.4 a). Extracted from [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF]. orthonormal basis of eigenfunctions, which satisfy

-∆ M φ i = λ i φ i . ( 3.4)
The approach of the decomposition of the mean curvature (C) is described by

C = +∞ i=0 C i φ i , C i = M Cφ i , ( 3.5)
where C i is the Fourier coefficient of the curvature in the eigenfunctions basis φ i . A band power spectrum is defined by

BS C (0) = C 2 0 , BS C (k) = i k 2 i=i k 1 C 2 i , (3.6)
where k is the number of the band, B0 power corresponds to the first constant eigenfunction, other bands correspond to sums of powers in different frequency bands defined via the wave-

Results

lengths. i k 1 and i k 2 correspond to the numbers of eigenfunctions associated to the band interval limits.

The approach of the local dominant band segmentation, which is presented in [START_REF] Germanaud | Larger is twistier: spectral analysis of gyrification (spangy) applied to adult brain size polymorphism[END_REF], allows us to parcellate the cortical surface into 7 spatial frequency bands, from B0 to B6. The bands B1, B2 and B3 are related to global brain shape. The last three bands (B4, B5 and B6) are associated with fold-related variations of curvature, which reveal the number of 1 st (primary), 2 nd and 3 rd order folds, respectively. The 3 rd order folds locate where the absolute values of curvature are the biggest on the brain surfaces. In consideration that we want to quantify the folding variations induced by the difference in the initial cortical thickness, we focus thus on the folding power of B4, B5 and B6.

Results

Validation of brain volume with biological time

Several numerical simulations of a fetal brain at 22 weeks of gestation are performed by using the combined biomechanical model. Compared with the simulation results of the previous model, the folding patterns have not been changed. We then compute the fetal brain volume of the simulations for the initial cortical thicknesses varying from 1.48 to 2.98 mm. A comparison between the volume data of our BLL growth model, the volume data using Tallinen's BLL growth model and other validated data of the literatures [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF] is shown in Figure 3.5. Comparing with the volume data of Tallinen (cyan points), our volume data (blue points), fitted with a second-order polynomial: BV (ml) = 0.59GA 2 -12.77GA + 53.46, is much more consistent with other validated brain volume data of the literatures. At the same gestational age, the different points represent the brain volume of the simulations with different initial cortical thicknesses. It is noted that the initial cortical thickness has a slight impact on the brain volume which can be neglected. We can conclude that our volume data (blue fitted curve) is consistent with these validated brain volume curves during gestational weeks, which indicates that our brain longitudinal length (BLL) growth model can mimic a realistic brain volume growth process.

Impact of cortical thickness onto morphology of brain

To understand the effect of the initial cortical thickness on the cortical folding patterns, we change the initial cortical thickness in the combined model from 0.74 to 5.96 mm. The cortical thicknesses from 1.48 to 2.98 mm are defined according to normative fetal brain measurements. The other two extreme cortical thicknesses (0.74 and 5.96 mm) are our hypotheses for fetal abnormal cortical thicknesses. The simulation results are shown in Figure 3.6. At the 29 weeks of gestation, the brain surfaces begin showing some clear line-like sulci. As the cortical thickness increases, the number of sulci decreases, the sulci became more isolated. Starting from 32 weeks of gestation, smooth gyri begin to appear on cortical surface. With the increase in the cortical thickness, the gyri became wider and fewer. The qualitative analysis is consistent with previous works such as [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF] relating the cortical thickness and the folding patterns. It can be found that a greatly increased cortical thickness associates with a severe decrease in number of gyri and sulci, which corresponds to the phenomenon of lissencephaly. The morphology with the initial cortical thickness of 5.96 mm in our simulations is similar to lissencephaly. On the contrary, polymicrogyria relates to an overly convoluted cortex, i.e., an increase in number and a decrease in size of gyri and sulci. The morphology with the thinnest initial cortical thickness of 0.74 mm resembles polymicrogyria.

To analyze quantitatively the folding patterns, we first compute the average of the absolute values of mean curvatures on the surfaces of the simulated brain, and the comparison is shown in Figure 3.7. For the brain surface with the thinnest cortical thickness, the curvature increases rapidly after 29 weeks of gestation and finally exceeds the values corresponding to other cortical thicknesses, indicating that its primary folds appear late but the 2 nd and 3 rd order folds are rapidly formed, and even more complex than the folding patterns of the others. For the surfaces with other cortical thicknesses, as the cortical thickness increases, the curvature decreases, corresponding to a smaller number of folds. Secondly, we calculate the 3D GI on each surface, the results of the initial cortical thicknesses between 1.48 and 2.98 mm are shown in Figure 3.8. It can be seen that there is a clear increase of the 3D GI with the increase in gestational age, regardless of the cortical thickness. The difference in the 3D GI is small at the same gestational week. However, the effect of the cortical thickness on the 3D GI is irregular. The GI based on surface areas may not be able to describe the complexity of the sulci and gyri when the initial cortical thickness changes little. The 3D GI computed on cerebral hemispheres for 12 healthy fetuses [START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF] is plotted for comparison. We remark a good agreement on the 3D GI at early gestational ages, but bias appears at later gestational ages may be due to differences in individual brain development. In future work, we will use more fetal brain MRI data at different gestational ages and compare the simulation results with the data to further validate the model. Furthermore, we also measure the local degree of folding proposed by Toro et al. [START_REF] Toro | Brain size and folding of the human cerebral cortex[END_REF] on the surfaces of the simulated brain (the initial cortical thickness varies from 0.74 to 5.96 mm). The average of the local surface ratio is shown in Figure 3.9. We observe that the results are almost consistent with the visual illustration (in Figure 3.6). For the brain surface with the thinnest cortical thickness, it starts with a smaller average of the local surface ratio, but then it folds faster to catch up with the folding degree of the surfaces with other cortical thicknesses. For the surfaces with other cortical thicknesses, as the cortical thickness increases, the number of folds decreases, and the average of the local surface ratio decreases. Therefore, the average of the local surface ratio may be a good measure of the complexity of the folding patterns.

Reconstruction and mesh generation of fetal brain

In addition, for different initial cortical thicknesses (H i ) from 0.74 to 5.96 mm, we analyze the curvature spectra and compute the folding power of B4, B5 and B6 on the brain surfaces. The total folding power is shown in Table 3.1. At the same gestational week, increasing the initial cortical thickness decreases the amount of growth in cortex which results in smaller total folding power, i.e., less total folds (1 st , 2 nd and 3 rd order folds). Moreover, we find that, for the thinner cortical thicknesses, the folding power of B6 accounts for an increasingly large proportion of the total folding power with the increasing gestational age. At each gestational week, the difference in B6 folding power is larger than that in B4 and B5 folding power between these cortical thicknesses, indicating that the small changes in cortical thickness have significant effect on the number of tertiary brain folds. These trends of B6 folding power for different cortical thicknesses are shown in Figure 3.10. The thinner the initial cortical thickness, the larger the increase of the B6 folding power demonstrates.

Reconstruction and mesh generation of fetal brain

In this study, we use segmentation of cortical grey matter, sub-cortical white matter and lateral ventricle of a fetal brain atlas between 23 and 37 gestational ages [START_REF] Serag | Construction of a consistent high-definition spatio-temporal atlas of the developing brain using adaptive kernel regression[END_REF]. To reconstruct a fetal brain volume at 23 weeks of gestation, we first combine the segmentation of cortical gray matter, sub-cortical white matter and ventricle at 23 weeks by using the addition operator of ITKTools (https://github.com/ITKTools/ITKTools), and then the brain is binarized with an appropriate threshold by using the threshold operator of ITKTools. After that, the filling tool in ITK-SNAP is used to manually fill the openings of the brain's outer surface.

ITK-SNAP (http://www.itksnap.org/pmwiki/pmwiki.php) is a software application used to segment structures in 3D medical images and it also provides several manual image tools. Finally, we use a Python algorithm Scipy.ndimage.morphology.binary _fill_holes which allows to do a binary dilation to fill the inner holes for closed solids.

The accuracy and convergence speed of FEM-based simulation depends on the quality of the mesh, our FEM-based three-dimensional brain growth model requires high-quality tetrahedral meshes to ensure the accuracy of folding patterns. The quality of a tetrahedral mesh mainly refers to the degree of regularity of its tetrahedra, which is generally measured in terms of minimal dihedral angle [START_REF] Dardenne | Variational tetrahedral mesh generation from discrete volume data[END_REF][START_REF] Du | Tetrahedral mesh generation and optimization based on centroidal voronoi tessellations[END_REF]. For highly distorted tetrahedra (slivers), the degree of dihedral angles can be as low as 0.01 or as high as 179.99 [START_REF] Ea | Mesh optimization: how to obtain good unstructured 3d finite element meshes with not-so-good mesh generators[END_REF][START_REF] Freitag | Tetrahedral mesh improvement using swapping and smoothing[END_REF][START_REF] Pablo D Zavattieri | Optimization strategies in unstructured mesh generation[END_REF], these silver elements affect the accuracy of numerical solutions in the application.

Since the volumetric mesh is typically generated from the surface mesh, thus we first study how to generate a high-quality surface mesh on the reconstructed brain surface. There are several criteria used to measure the quality of the surface elements, such as minimum angle, element volume, minimum edge and shape quality, etc. (https://www.gidhome.com/documents/ referencemanual/PREPROCESSING/MeshMenu/Meshquality). The shape quality criterion measures the likeness of the triangle to an equilateral triangle (3 equal sides and 3 equal angles (60 • )), thus the surface mesh should have as many equilateral triangles as possible. To achieve it, we first use Laplacian filter in 3DSlicer (https://www.slicer.org/) to smooth the brain globally since the the three-dimensional reconstructed brain is arranged in layers and its surface is rough. However, if we over-smooth it globally, it will cause geometric distortions, such as shallower longitudinal fissure between the two hemispheres (as shown in the 4 th image in Figure 3.11(b)), but insufficient smoothing could result in some remaining bumps on the surface of the brain, especially at the boundary. Therefore, a balance needs to be made between excessive smoothness and rough surface. 3DSlicer also allows us to generate the triangular surface mesh on the surface of the brain using Delaunay triangulation and save it in a .stl format file. If there are many sharp (the degree of angles is too low or too high) or irregular (3 unequal sides and 3 unequal angles) triangles located on the bumps, we need to repair the vertices of these triangles manually by using Blender software (https://www.blender.org/) which allows to manually and locally change the elements of the mesh.

Secondly, based on the repaired brain triangular mesh, we generate the tetrahedral mesh by using Netgen, as shown in last image in Figure 3.11(b). It can be observed that the mesh is denser where the curvature is large or there are singularities. Netgen (https://ngsolve.org/) is a tetrahedral advancing front mesh generation tool and provides a hull mesh adaptation algorithm in order to improve the input hull mesh prior to the volume meshing step [START_REF] Schöberl | Netgen an advancing front 2d/3d-mesh generator based on abstract rules[END_REF]. There are some parameters defined in Netgen, such as mesh granularity, steps, etc. The mesh granularity can be chose according to the need for calculation accuracy or manually defined by changing the parameters of the max/min mesh-size and the mesh-size grading. To obtain more uniform mesh elements, we set the max and min mesh-size as close as possible. In addition, there is a refinement tool in Netgen which can be used to further refine the mesh. For meshing steps, we select all of the generation and optimization steps, which include analyze geometry, mesh edges, mesh surface, optimize surface, mesh volume and optimize volume. Other parameters are set as default. The flowchart and corresponding diagram of the pipeline that we use to generate volumetric meshes are shown in Figure 3.11. An example of generated geometries and tetrahedral meshes is shown in Figure 3.12(b). Compared with the geometry and the corresponding tetrahedral mesh from Tallinen et al. (Figure 3.12(a)), the longitudinal fissure between two hemispheres that we reconstruct is overly smooth and not deep enough. The boundaries of the occipital lobe are also over-smoothed, resulting in distortion of some geometric features. In addition, the generated mesh elements are not as regular as those provided by Tallinen et al., especially in the bumps on the surface, which may need to be repaired and corrected more finely. Since our brain mesh is distorted due to over-smoothing, we use the fetal brain mesh provided by Tallinen et al. in this and next chapters.

Conclusion

In this chapter, to know whether the relation between the increase in human brain volume and gestational age in the biomechanical cortical folding model is realistic, we measure the brain volume of the simulation results and compare it with the biological measurement results. We find that the simulated brains lack volume growth during the folding process. Thus we introduce a brain longitudinal length growth model derived from normative fetal brain measurement into the brain folding model. This leads to the model being able to simulate a more realistic early brain volume growth, which has been validated by previous works of the literature. In addition, it should be noted that this correction does not change the folding patterns.

Furthermore, to understand the effect of the initial cortical thickness on folding patterns and complexity, we vary the initial cortical thickness in the combined model from 0.74 to 5.96 mm to simulate the folding process of the fetal brain. We evaluate the surface morphology of different cortical thicknesses through visual analysis and quantitative features. The results show that thinner cortical thickness leads to higher B6 folding power that corresponds to tertiary brain folds. In addition, we find that the morphology with the initial cortical thickness of 5.96 mm is similar to lissencephaly and that with the thinnest initial cortical thickness of 0.74 mm resembles polymicrogyria.

At the end of this chapter, we use various tools, such as ITKTools, ITK-SNAP, Scipy, 3DSlicer, Blender and Netgen, to study how to reconstruct a human brain volume from brain segmentation and then generate a tetrahedral mesh on the reconstructed human brain volume. However, the smoothing of the reconstructed human brain changes its geometry, especially the area of the longitudinal fissure and the surface boundaries. Therefore, we should further study whether there is an approach to ensure the accuracy of the initial geometry while smoothing. In addition, the manual correction of triangular mesh is complicated and greatly increases the time cost, to improve the efficiency and robustness of generating high-quality meshes, in future works, we should explore whether there is an appropriate automatic surface mesh correction method that can replace manual correction.

Introduction

In the human brain growth model proposed by Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the tangential cortical expansion is regarded as a parameter that is spatially invariant on the cortical surface and increases linearly with time. However, many studies on the human brain have shown that different regions of the cortex have different growth rates [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Habas | Early folding patterns and asymmetries of the normal human brain detected from in utero mri[END_REF][START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF]. Through the deformation tensor analysis, Rajagopalan et al. [START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF] found that the local growth in the cortical plate was either similar to or significantly greater than the overall cerebral rate in the fetus from in utero MRI. The deformation tensor model was used in the cortical plate because it is sensitive to focal changes in both cortical thickness and area [START_REF] Vidya Rajagopalan | Local tissue growth patterns underlying normal fetal human brain gyrification quantified in utero[END_REF]. Habas et al. [START_REF] Habas | Early folding patterns and asymmetries of the normal human brain detected from in utero mri[END_REF] identified regions of the fetal brain surface that undergo significant folding changes during 20-28 gestational weeks, and extracted spatio-temporal patterns of early cortical folding based on measurements of local surface curvature. Wright et al. [START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF] observed that regional differences in growth rate of fetal brains by parceling each fetal cortex into nine lobar regions, with the parietal and posterior temporal lobes showing the fastest growth, while the cingulate, frontal and medial temporal lobes develop more slowly.

In addition, Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] also observed significant regional growth differences across the cortical surface of 30 preterm infants scanned from approximately 28-38 wk postmenstrual age. They demonstrated that these patterns of cortical growth and maturation are consistent with those in past literature [START_REF] Ball | Richclub organization of the newborn human brain[END_REF][START_REF] Amy R Deipolyi | Comparing microstructural and macrostructural development of the cerebral cortex in premature newborns: diffusion tensor imaging versus cortical gyration[END_REF][START_REF] Dubois | Mapping the early cortical folding process in the preterm newborn brain[END_REF]. Specifically, Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] observed that the cortical expansion is significantly higher in the lateral parietal, temporal and occipital regions (red in Figure 4.1) and significantly lower in the medial and insular regions (blue in Figure 4.1) during the period from 30 to 38 wk. These growth patterns are consistent across right and left hemispheres. Moreover, they found that these growth patterns change during the rapid cortical expansion period (28-38 gestational ages), with the noninjured subjects following a highly consistent trajectory. Specifically, they observed that the relative expansion in the primary sensory, motor and visual cortices, as well as in the insula, decreases over time (green in In this chapter, we explore the growth data provided by Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] and attempt to answer the following questions: 1) What kind of local folding patterns will be produced by this spatio-temporal differential cortical expansion? 2) Compared to the global (spatially invariant) linear cortical expansion originally defined in the model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], how does this spatio-temporal differential cortical expansion affect folding? For these purposes, we model and integrate the spatio-temporal differential cortical growth mechanism into the combined brain growth model presented in Chapter 3, allowing to apply a spatio-temporal cortical expansion map [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] to a fetal brain to simulate the folding process. Specifically, we first derive the relation between the cortical expansion defined by Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] and Tallinen et al. [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Then we develop a pipeline that applies the spatio-temporal cortical expansion map to a fetal brain model at 22 weeks of gestation for simulation. From the simulation results, we can visually remark the difference in the appearance of local folding patterns over time. Then we measure the surface morphology through the mean curvatures [START_REF] Rusinkiewicz | Estimating curvatures and their derivatives on triangle meshes[END_REF], the surface-base three-dimensional gyrification index [START_REF] Clouchoux | Quantitative in vivo mri measurement of cortical development in the fetus[END_REF] and the folding orientation as introduced in Chapter 2.

Derivation of relationship between cortical expansions

The purpose of this section is to find the relationship between the cortical expansion defined by Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] and Tallinen et al. [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] in order to apply the spatio-temporal cortical growth map to the brain growth model of Tallinen [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] to that defined in the brain growth model of Tallinen et al. [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

The spatio-temporal cortical growth map (from https://balsa.wustl.edu/study/show/ K65Z) includes the relative cortical expansion (RCE) of each vertex for each individual at different gestational ages. The relative cortical expansion (RCE) is defined as the ratio of the local cortical expansion (LCE) to the global cortical expansion (GCE) by Garcia et al. [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF], which can be mathematically represented by

RCE i (t GA ) = LCE i (t GA ) GCE(t GA ) , ( 4.1)
where i is the index of the vertex on the brain surface, t GA represents the gestational age. GCE(t GA ) is defined as a global areal expansion at t GA compared to the reference gestational age t GA -∆t GA , which can be expressed as

GCE(t GA ) = A(t GA ) A(t GA -∆t GA ) , ( 4.2)
where A(t GA ) = Nvertex i=1 A i (t GA ) is the total cortical area at the gestational age t GA . The data of GCE is provided by Garcia et al., thus we can obtain the values of local cortical expansion. The LCE is also defined as

LCE i (t GA ) = A i (t GA ) A i (t GA -∆t GA ) , ( 4.3)
where A i (t GA ) is the local area around vertex i at the gestational age t GA and A i (t GA -∆t GA ) could be regarded as the local area around vertex i at the reference gestational age. The relative tangential expansion ratio of the grey matter to the white matter defined by Tallinen et al. [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] is not related to the real time (gestational age), we suggest using the three-dimensional gyrification index (GI) to make the link between the time of model (t m ) and the gestational age (t GA ). To achieve it, we compare the GI computed on the simulated surfaces (GI simu ) of the cortical growth of α(t m ) = 1.829t m at different t m with that presented in Armstrong et al. (GI real ) which was measured on human brains at different t GA [START_REF] Armstrong | The ontogeny of human gyrification[END_REF]. For a given t GA , it is assumed that GI real can be expressed by a function of f :

GI real = f (t GA ), (4.5)
GI simu can be computed at each t m through a function of g:

GI simu = g(t m ). (4.6)
Thus for a given t GA the corresponding t m should minimize the error between g(t m ) and f (t GA ):

t m = q(t GA ) = arg min |g(t m ) -f (t GA )|. ( 4.7)
Since the Gompertz model can be used to represent the cortical growth with t GA and the cortical growth defined in the model varies linearly with t m , the relation between t m and t GA is fitted with a Gompertz curve as t m = 0.987e -e -0.134(t GA -29.433 ) . (4.8)

We show here how to calculate the cortical expansion defined by Tallinen et al. [START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] by using the local cortical expansion computed from the growth data of Garcia et al.. Since the cortical expansion g 2 (0, t m ) is an areal growth at t m relative to the victitious stress-free state (initial time t m0 in the model), for the i th vertex of the brain surface, we can obtain

g i (0, t m ) = A i (t m )/A i (t m0 ).
Thus the magnitude of expansion at t m is

α i (t m ) ≈ g i (0, t m ) -1 = A i (t m )/A i (t m0 ) -1, which is derived from Equation 4.4.
With the relation between t m and t GA , we compute the corresponding t m , t m -∆t m , . . . for t GA , t GA -∆t GA , . . . of the growth data of Garcia et al., then we can obtain the local cortical expansion of these time of model (LCE i (t m ), LCE i (t m -∆t m ), . . .). Thus the α i at t m can be calculated by using LCE i (t m ), LCE i (t m -∆t m ), ... as:

α i (t m ) = A i (t m ) A i (t m0 ) -1 = A i (t m ) A i (t m -∆t m ) A i (t m -∆t m ) A i (t m -2∆t m ) • • • A i (t m0 + ∆t m ) A i (t m0 ) -1 = LCE i (t m )LCE i (t m -∆t m ) • • • LCE i (t m0 + ∆t m ) -1.
(4.9)

Application on a fetal brain

To apply the cortical expansion maps corresponding to the noninjured semi-brain surface meshes (from https://balsa.wustl.edu/study/show/K65Z) to the tetrahedral mesh of an entire brain at 22 weeks of gestation (provided by Tallinen et al.), we define a pipeline as shown in Figure 4.3. Each step will be detailed in the following subsections. In order to close the mesh boundaries, we adopt a robust hole-filling algorithm for triangular mesh [START_REF] Zhao | A robust hole-filling algorithm for triangular mesh[END_REF]. This algorithm creates the new triangles start from the boundary vertex with the smallest angle, and then creates the new vertices from the bisector of the angle formed by the three vertices of the edge, as shown in Figure 4.4. If the created vertex is close enough to every related boundary vertex (less than the given threshold), these vertices will be merged. The stop criterion of the algorithm is until the whole region is patched by all newly created triangles. Finally, all new vertices are locally smoothed, the smoothed coordinates of each new vertex correspond to the barycenter of its neighbors.

Mesh cutting and boundary closure

Mesh registration

To register the left and right semi-brain meshes relative to the semi-brain meshes provided by Garcia et al. from https://balsa.wustl.edu/study/show/K65Z, we use the spherical demons algorithm (Fast Diffeomorphic Landmark-Free Surface Registration) which was proposed by Yeo et al. [START_REF] Bt Thomas Yeo | Spherical demons: Fast surface registration[END_REF][START_REF] Bt Thomas Yeo | Spherical demons: fast diffeomorphic landmark-free surface registration[END_REF]. They showed that a large class of regularizors for the modified demons objective function can be efficiently approximated on the sphere [START_REF] Bt Thomas Yeo | Spherical demons: fast diffeomorphic landmark-free surface registration[END_REF] by using the theory of spherical vector spline interpolation [START_REF] Glaunès | Landmark matching via large deformation diffeomorphisms on the sphere[END_REF] and other differential geometric tools (mean curvature, sulcal depth). The modified demons objective function is given by

(s * , c * ) = arg min s,c F -M • c 2 + 1 σ 2 x dist(s, c) + 1 σ 2 T Reg(s), ( 4.10)
where F is the fixed mesh, M is the moving mesh, c is the desired registration and s is a hidden transformation that acts as a prior on c. The fixed mesh F and warped moving mesh M • c are treated as N × 3 vectors. Generally, dist(s, c) = s -c 2 , encouraging the resulting transformation c to be close to the hidden transformation s and Reg(s) = ∇s 2 , i.e., the regularization penalizes the gradient magnitude of the hidden transformation s. σ x and σ T provide a tradeoff among the different terms of the objective function.

They demonstrated that by using a limited class of diffeomorphisms combined with the demons algorithm, the resulting registration is diffeomorphic, and the speed gain is more than one order of magnitude compared with other landmark-free invertible spherical registration methods, such as [START_REF] Fischl | High-resolution intersubject averaging and a coordinate system for the cortical surface[END_REF][START_REF] Bt Thomas Yeo | Effects of registration regularization and atlas sharpness on segmentation accuracy[END_REF]. In addition, they validated this algorithm by demonstrating an accuracy comparable to that of the popular FreeSurfer algorithm [START_REF] Bt Thomas Yeo | Spherical demons: fast diffeomorphic landmark-free surface registration[END_REF].

We show the steps of applying the cortical expansion provided by Garcia et al. to our brain mesh for simulation in Figure 4.6. The relative cortical expansion of 28 weeks GA is projected on the atlas surface mesh, on our registered surface mesh and also on the surface of the volumetric mesh of 22 weeks GA which is used for simulation. The registered semi-brain surface mesh after resampling has the same number of vertices as that of the atlas surface mesh. For each vertex on the surface of the volumetric mesh which is used for simulation, we use Kd-Tree to find its closest vertex on the registered surface mesh, and then assign the lobar label of this vertex to it. Since the surface of the volumetric mesh has more vertices than the registered surface, the cortical expansion map after using Kd-Tree will be blurred.

Lobes parcellation and labels assignment

Referring to the parcellation methods of human brain lobes proposed by Wright et al. and Vasung et al. [START_REF] Vasung | Quantitative and qualitative analysis of transient fetal compartments during prenatal human brain development[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF], and expecting the central sulcus area to be divided into a separate lobe, we manually split the surface of a brain atlas mesh at 30 weeks of gestation (from https: //balsa.wustl.edu/study/show/K65Z) into 11 lobes using the surfpaint in Anatomist of BrainVISA (http://brainvisa.info/web/index.html), and then we give each lobe a label, as shown in Figure 4.7. Finally, for each tetrahedron of the mesh for simulation, we use Kd-Tree to find the nearest surface vertex to its centroid, and assign the lobar label of this surface vertex to it, which will be used to allocate the cortical expansion of each vertex of the atlas mesh.

Curve-fitting of average cortical expansion

Average cortical expansion From https://balsa.wustl.edu/study/show/K65Z, there are 4 noninjured subjects with 3 relative cortical expansion (RCE) maps corresponding to 3 periods from approximately 28-38 weeks postmenstrual age (can also be called gestational age before birth) with a restriction of scans < 6 weeks [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF]. For the 1 st subject, three gestational periods are 28-30, 30-33 and 33-38 weeks; for the 2 nd subject, three gestational periods are 27-31, 31-33 and 33-37 weeks; for the 3 rd subject, three gestational periods are 27-29, 29-33 and 33-36 weeks; for the 4 th subject, three gestational periods are 27-30, 30-34 and 34-36 weeks. The global cortical expansion (GCE) of each period for the 4 subjects are shown in Table 4.1.

For each subject and in each period, we multiply the RCE of each vertex by the corresponding GCE in this period to obtain the local cortical expansions (LCE) of this period by using the Equation 4.1. In order to obtain LCE of a t GA to the initial t GA0 of this subject, we just need to do a product. For example, the LCE(28 -33wk) of the 1 st subject is calculated by LCE(28 -30wk)×LCE(30 -33wk). We then compute the average LCE of each lobe of the left and right semi-brains by using

N vertex i=1

LCE i Nvertex

, the results for the right semi-brain are shown in Figure 4.9.

Considering that the relative tangential expansion ratio (g) defined by Tallinen et al. is a relative length growth based on an initial point of the fetal brain at 22 weeks of gestation, thus g(t GA = 22wk) = 1. However, for the LCEs (local areal growth) of the 4 subjects provided by Based on the normalized values of the root of the average LCEs, the values of α(t m ) are calculated by the Equation 4.9. Figure 4.10 shows the values of α(t m ) of 4 subjects and the original global cortical expansion model (α(t m ) = ( √ 8-1)t m and t m = 0.987e -e -0.134(t GA -29.433 ) ). We can observe that as the gestational age increases, the normalized real growth data of the left and right semi-brains of the 4 subjects are slightly larger than the original growth. This may be because the original cortical growth model was obtained from an atlas, thus the cortical surface should be slightly smoother, i.e., it is a little less folded.

Curve-fitting of cortical growth with time of model Based on the relation between the time of model (t m

) and the gestational age (t GA ), we transform the average cortical expansion (α(t m )) related to t GA to that related to t m , the results are shown in Figure 4.11. After weighing the accuracy of the fitting and the simplicity of the model, we choose a linear model to fit the average cortical expansion (α(t m )) with the time of model (t m ), which is defined as In addition, for each lobe, we also use the linear model to fit all the 12 constraints of the 4 noninjured individuals to obtain the curves which can represent the cortical growth of all individuals, the results are shown in Figure 4.12. These fitted growth curves are used for simulation. The corresponding cortical growth rate (the slope of the growth curve presented in 4.12) of each lobe of the left and right semi-brains is shown in Figure 4.13. It can be observed that the growth rate is higher in the temporal and parietal lobes, while it is lower in the connected region, the outer ring of the limbic and anterior frontal lobes.

α(t m ) = at m + b, ( 4.11)

Results

Visual analysis of folding patterns

Based on the linear curve-fitting of all lobar average cortical expansions of the 4 noninjured individuals shown in Figure 4.12, the simulation result of a 22 weeks' fetal brain with the initial cortical thickness of 2.5 mm is shown in Figure 4.14. Comparing to the simulation result of the global linear growth mechanism (α(t m ) = 1.829t m ) at different time of model, we can observe that on the surface of the spatio-temporal growth, the primary folds are almost at the same position as those on the surface of the global linear growth. The folds seem to be slightly more complex on the surface of the spatio-temporal growth than those on the surface of the global linear growth at each time, which is consistent with the growth curves presented in Figure 4.12.

According to the fitted growth curves shown in 4.12, the cortical expansions of different postmenstrual ages are displayed on the simulated surfaces of the spatio-temporal growth, as shown in Figure 4.15. The cortical expansion of the left and right semi-brains is asymmetrical. In each corresponding lobe, the cortical expansion of the left semi-brain is smaller than that of the right semi-brain. The cortical expansion is higher in the lobe of the central sulcus compared to other lobes at 28.6 weeks of gestation, and then the higher expansion occurs in the parietal, occipital, medial and posterior temporal lobes starting at 34 weeks of gestation, which is consistent with the observations presented in [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF][START_REF] Wright | Automatic quantification of normal cortical folding patterns from fetal brain mri[END_REF]. increases than other lobes from the start to the time of model 0.7, which is consistent with the cortical expansion presented in Figure 4.15. We then compute the 3D GI for each lobe of the spatio-temporal cortical growth, the results of principal lobes are shown in Figure 4.17. The lobes 4, 6 and 10 also have higher 3D GI increases than other lobes in the whole time period, which is in agreement with the growth curves presented in 4.12. It indicates that the higher cortical expansion can lead to more complex folding patterns, especially in the medial and posterior temporal, parietal and central sulcus lobes.

Comparison of spatio-temporal vs global growth

To answer whether the spatio-temporal differential cortical expansion affects the folding complexity compared to the global linear expansion, we compare the average of the absolute values of dimensionless mean curvatures and the 3D GI of each lobe between the two growth mechanisms. The results of principal lobes are shown in Figures 4.18 and 4.19. We can observe that almost all the lobes show the higher curvature and 3D GI on the surface of the spatio-temporal cortical growth compared to that of the global linear growth, which is consistent with the cortical expansion curves of the two growth mechanisms shown in Figure 4.12. These measurements indicate that this spatio-temporal cortical expansion can lead to more complex folding patterns compared to the original global linear cortical expansion. Furthermore, in order to quantify the correlation of the local folding patterns generated by the two growth mechanisms under the same surface folding complexity, we first search for the same 3D GI of each principal lobe for the two growth mechanisms by shifting the time of model, because the 3D GI is not the same for them at the same time of model, which is shown in Figure 4.20. The 3D GI depends on the cortical expansion, i.e., depends on the curves presented in Figure 4.12. Then at each same lobar 3D GI (corresponding to different time of the two growth mechanisms), we calculate the Pearson correlation coefficient of dimensionless mean curvatures for all vertices in each principal lobe of the surfaces produced by the two growth mechanisms. The results of the left and right semi-brains are shown in Figure 4.21. A visualization of the correlation coefficient of dimensionless mean curvatures between the two growth mechanisms at maximum 3D GI is shown in Figure 4.22. It can be seen that the correlation between the left and right semi-brains is asymmetric, with higher correlation in the lobes 5 and 10 for the left semi-brain and in the lobes 4 and 5 for the right semi-brain at maximum 3D GI.

Comparison with a fetal brain atlas

We also compare the orientation of the primary folds on the simulated surface of the spatiotemporal growth with that on the cortical surface of a fetal brain atlas [START_REF] Gholipour | A normative spatiotemporal mri atlas of the fetal brain for automatic segmentation and analysis of early brain growth[END_REF], as shown in Figure 4.23. To do this, we first generate a triangular mesh on the cortical surface of each hemisphere of the fetal brain atlas at 29 weeks GA and use a Laplacian filter to slightly smooth the mesh. Then we compute the fold angle for each vertex of the two surfaces (fetal brain atlas and simulation) at 29 weeks GA using our method proposed in Chapter 2. Finally, with respect to the cortical surface mesh of the fetal brain atlas, the simulated surface mesh and corresponding fold angles are registered and resampled by using the spherical demons algorithm [START_REF] Bt Thomas Yeo | Spherical demons: Fast surface registration[END_REF][START_REF] Bt Thomas Yeo | Spherical demons: fast diffeomorphic landmark-free surface registration[END_REF]. It can be observed that for each hemisphere, the distribution of folds orientation in the lobes of central sulcus (red frames) and the frontal lobes (black frames) of the two surfaces is similar, although there are still differences in the appearance of the folds.

Conclusion

In this chapter, in order to study the effect of the regional cortical expansion on folding patterns, we replace the spatially invariant linear cortical expansion with a spatio-temporal differential cortical expansion in the biomechanical human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. The derivation of the relationship between two different cortical expansion definitions allows us to use real spatiotemporal differential cortical expansion maps in the model. In addition, our proposed pipeline for mapping the spatio-temporal cortical expansion data corresponding to a brain atlas mesh to a fetal brain mesh for simulation could also be used to apply other cortical expansion maps to the human brain.

In order to know the effect of the spatio-temporal differential cortical expansion on local folding patterns, we compute the curvature, the three-dimensional gyrification index and compare them between different lobes. The results show that the complexity degrees of the folding patterns in the medial and posterior temporal, parietal and central sulcus lobes are higher than that of other lobes. In addition, we point out that in the frontal lobe and the lobe of central sulcus, the orientation of the primary folds on the simulated surface is similar to that on the cortical surface of a fetal brain atlas.

To understand whether this spatio-temporal differential cortical expansion will affect the folding compared to the spatially invariant linear cortical expansion, we also compare the curvature and three-dimensional gyrification index for each lobe of the two growth mechanisms. The results show that almost all the lobes of the spatio-temporal cortical expansion have higher curvatures and 3D GI than those of the global linear cortical expansion, indicating that this spatio-temporal cortical expansion can lead to more complex folding patterns compared to the global linear cortical expansion.

In order to evaluate the simulation results, in addition to the folds orientation, we should use various quantitative metrics to compare the simulated folding patterns with those on human brains in future works.

Conclusion

The objective of this thesis is to better understand the early development of the human brain and explore the causes of abnormal cortical folding patterns. Four main contributions are presented: 1) an investigation of the effect of the biophysical parameters (the cortical growth, the initial geometry and cortical thickness) on the surface morphology based on the biomechanical model proposed in [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], 2) an approach for computing the fold angles and the anisotropy of the folding orientation, 3) a modeling of the brain longitudinal length growth to improve the volume accuracy of the brain model, 4) a modeling of the spatio-temporal differential cortical expansion mechanism and approaches for applying of a cortical expansion map to a fetal brain. These contributions are motivated by the fact that the image-based techniques cannot provide insights on causes of neurodevelopmental disorders and disabilities without biophysical information, and the computational modeling is a powerful tool to help understand the causes and improve predictive power. In addition, since the biomechanical model based on the hypothesis of the differential tangential growth can produce the realistic patterns of convolutions on the brain surface over a developmental process [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the studies of brain development in this thesis are based on this model. The first contribution of this thesis, presented in Chapter 2, consists of an exploration of the behavior of the biomechanical model by investigating the impact of the mesh density and the biophysical parameters (the cortical expansion, the initial geometry and cortical thickness) onto the surface morphology using ellipsoids. We first show that when the mesh density reaches a certain order of magnitude, the simulated surfaces can achieve sufficient folding accuracy, which is useful for the further simulations based on the growth of human brains using this model.

In addition, we demonstrate that the cortical growth mode does almost not affect the complexity degree of surface morphology; the variation in the initial geometry changes the folds orientation and depth, and in particular, the slenderer the shape is, the more folds along its longest axis could be seen and the deeper the sulci become. Moreover, the thinner the initial cortical thickness is, the higher the spatial frequency of the folds is, but the shallower the sulci become, which is in agreement with the previously reported effects of cortical thickness [START_REF] Budday | A mechanical model predicts morphological abnormalities in the developing human brain[END_REF][START_REF] Wang | On early brain folding patterns using biomechanical growth modeling[END_REF]. These results tend to show that the use of such biomechanical models could highlight the links between neurodevelopmental disorders and physical parameters.

The second contribution of this thesis is also introduced in Chapter 2, is a novel approach to calculate the fold angles and measure the anisotropy of the folding orientation through geometric tools [START_REF] Lefevre | Fast surface-based measurements using first eigenfunction of the laplace-beltrami operator: Interest for sulcal description[END_REF][START_REF] Sylvain Petitjean | A survey of methods for recovering quadrics in triangle meshes[END_REF] and the Kullback-Leibler divergence, which could be used for future works to measure the folding orientation on any folded surface.

The third contribution is oriented toward improving the accuracy of the volume growth of the human brain during the deformation process, which is presented in Chapter 3. We propose a brain longitudinal length (BLL) growth model derived from normative fetal brain measurement [START_REF] Kyriakopoulou | Normative biometry of the fetal brain using magnetic resonance imaging[END_REF] and introduce it into the biomechanical model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. Then we validate this BLL model by comparing the simulated brain volume with other validated brain volume data of the literatures [START_REF] Armstrong | The ontogeny of human gyrification[END_REF][START_REF] Chang | The assessment of normal fetal brain volume by 3-d ultrasound[END_REF][START_REF] Qy Gong | Fetal and fetal brain volume estimation in the third trimester of human pregnancy using gradient echo mr imaging[END_REF][START_REF] Hsu | Quantitative analysis of normal fetal brain volume and flow by three-dimensional power doppler ultrasound[END_REF].

Based on the combined human brain growth model presented in Chapter 3, an investigation of the effect of the initial cortical thickness on human cerebral cortical folding patterns is performed. We find that the thinner cortical thickness can lead to higher folding power that corresponds to tertiary brain folds. In addition, we show that the morphology with the initial cortical thickness of 5.96 mm is similar to lissencephaly and that with the initial cortical thickness of 0.74 mm resembles polymicrogyria. This observation may be important for understanding the causes of several neurodevelopmental disorders associated with abnormal cortical folding patterns.

Considering that the biomechanical human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] needs a highquality tetrahedral mesh as the initial point to produce the accurate brain folding patterns and we usually only have human brain segmentation data, thus we also propose a pipeline for generating a tetrahedral mesh from segmentation of the human brain at the end of Chapter 3. However, the pipeline needs to be optimized in terms of the geometric accuracy and the efficiency of mesh correction.

The last contribution of this thesis, introduced in Chapter 4, is a modeling of the spatiotemporal differential cortical growth mechanism and an integration of it in the biomechanical human brain growth model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF]. We first show a relationship derivation of two different cortical expansion definitions, and then propose a pipeline for applying a spatio-temporal differential cortical expansion map [START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] to a fetal brain at 22 weeks of gestation. The simulation results show that the complexity degree of folding patterns in the posterior temporal, parietal and central sulcus' lobes is higher than that of other lobes. In addition, we point out that in the frontal lobe and the lobe of central sulcus, the orientation of the primary folds on the simulated surface is similar to that on the cortical surface of a fetal brain atlas [START_REF] Gholipour | A normative spatiotemporal mri atlas of the fetal brain for automatic segmentation and analysis of early brain growth[END_REF].

Perspective

Geometric smoothing and mesh correction

We define a pipeline for generating tetrahedral meshes from MRI segmentation of the human brain. However, it causes changes in the shape of the brain due to the global smoothing. In addition, it greatly increases the time cost because of the local manual correction. In order to ensure the accuracy of the initial geometry and improve the efficiency and robustness of generating high-quality meshes, it is necessary to explore whether there is a more suitable geometric smoothing method to preserve the characteristics of the longitudinal fissure between the two cerebral hemispheres and also the surface boundaries, and look for an automatic surface mesh correction method that can replace manual correction.

Initial geometry

For the study of the initial geometry, in addition to the elongation ratio and the longitudinal fissure, other geometric changes may also affect the folding patterns, such as the initial curvature. If the initial human brain surface is not smooth, what folding changes will occur? Therefore, another important future work would be to study the effect of other geometric changes on folding. The method that we propose to calculate the fold angles and measure the anisotropy of the folding orientation could be used for future works to measure the orientation of the folds on surfaces with other geometric changes.

Cortical thickness

In humans, the thickness of the cortex varies between brain regions. Gyral regions of the human brain are thicker than sulcal regions: the average thickness of gyral crown is 2.7 mm thick, while sulcal fundi are in average 2.2 mm thick [START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF]. Besides this local geometry-dependent change, cortical areas also have various thicknesses. For example, the thickness of the primary visual cortex is approximately 2.8 mm, while the primary auditory cortex is around 3 mm [START_REF] Fischl | Measuring the thickness of the human cerebral cortex from magnetic resonance images[END_REF]. Therefore, it's interesting to develop an atlas of cortical thickness that varies by region of the human brain in future works, and use it in the human brain growth model to study the relation between regional thicknesses and folding patterns.

Model improvement

In future works, we should further qualitatively and quantitatively compare the simulated folding patterns with the real folding patterns of the human brain. For the differences in folding patterns, the idea will be to use a data-driven approach to learn the dynamics of the cortical evolution and obtain the growth deformation gradient tensor in the model from growth deformation data of human brains, such as the approach presented in [START_REF] Wang | An inverse modelling study on the local volume changes during early morphoelastic growth of the fetal human brain[END_REF], to further improve the biophysical model. Abstract: Perinatal ischemic stroke constitutes the leading cause of unilateral cerebral palsy or epilepsy in term-born children. However, the causes of these observed disabilities remain unclear. In this context, computational modeling is a powerful tool to provide a better understanding of the early brain folding process. Recent studies based on biomechanical modeling have shown that mechanical forces play a crucial role in the formation of cortical convolutions. However, the effect of physical parameters in these models, and the correlation between simulation results and biological facts remain unclear. In this thesis, using a biomechanical brain folding model, we investigate the effect of the cortical growth, the initial geometry and the initial cortical thickness on cortical folding patterns. In addition, we improve the biomechanical model by adding a new brain longitudinal length growth model and a spatio-temporal differential cortical expansion mechanism to the model. Furthermore, in order to quantify the surface morphology of simulations, several descriptors of the folds are used such as curvaturebased measures, gyrification index and sulcal depth. Besides, we introduce a novel approach to measure the folding orientation through geometric tools.

5 1. 2 6 1. 3 6 1. 4 9 1. 5 10 1. 6 11 1. 7 13 1. 8 14 v 1 .

 52636495106117138141 . . . . . . . . . . . . . . . . . . . . . . . . . xix 5 (a) Organigramme et outils utilisés pour générer un maillage tétraédrique à partir de la segmentation d'IRM; (b) Diagramme correspondant à l'organigramme. . . . xix 6 Pipeline de l'application d'une carte d'expansion corticale spatio-temporelle [65] à un cerveau foetal de 22 semaines de gestation. . . . . . . . . . . . . . . . . . . . xx 1.1 Three major parts of the brain: cerebrum, brain stem and cerebellum. Extracted from https://www.myshepherdconnection.org. . . . . . . . . . . . . . . . . . . Gyri and sulci. Extracted from [22]. . . . . . . . . . . . . . . . . . . . . . . . . . Six major lobes. (a) Frontal lobe, temporal lobe, occipital lobe, parietal lobe, and insular lobe. Extracted from [1]. (b) Limbic lobe. Extracted from Wikipedia. . . Nine-region parcellation. Extracted from [189]. . . . . . . . . . . . . . . . . . . . Chronology of gyrification in human. GW 31-32 divides phases 1 and 2. There is a gradual lack of consistency in cortical folding patterns among species in the final stages of gyrification, as minor developmental changes in gyri and sulci become more specialized to species and susceptible to local environmental and ultimate experiential variations. Extracted from [116]. . . . . . . . . . . . . . . . . . . . . Different cortical folding models of the emergence of gyrification. Extracted from [167]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Components of stress in three dimensions. Extracted from Wikipedia. . . . . . . Mapping between the reference configuration and the deformed configuration. . . 14 1.9 (a) Strain-stress relation: difference between linear elastic and hyperelastic. (b) Stress-strain curves for various hyperelastic material models. Extracted from Wikipedia. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 Elastic, viscoelastic and unstable behaviour in response to mechanical loading.

  on a fetal brain at 22 weeks of gestation. . . . . . . . . . . . . . . . . . . . . . . . . . . 81 4.4 Rules for creating triangles: a)θ i ≤ 75 • ; b) 75 • < θ i ≤ 135 • ; c) θ i > 135 • .

  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 4.5 Cutting of the entire brain mesh into the semi-brain meshes and the closure of the semi-brain mesh boundaries. . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.6 Application of relative cortical expansion provided by Garcia et al. on the atlas, on the registered semi-brain surface meshes, and on the surface of the semi-brain volumetric mesh by using the spherical demons registration and the Kd-Tree method. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.7 Brain lobes parcellation is performed manually on the surfaces of the left and right semi-brains and each surface is split into 11 lobes (including 10 cortical lobes and 1 connected lobe): 1. connected lobe/region (gray), 2. medial and posterior frontal lobe (dark blue), 3. occipital lobe (navy blue), 4. medial and posterior temporal lobe (blue), 5. anterior temporal lobe (bright blue), 6. parietal lobe (blue-green), 7. insular lobe (green), 8 and 9. outer ring of the limbic lobe [gyrus fornicatus encompassing: cingulate gyrus with subcallosal area (bright green) and parahippocampal gyrus (gold)], 10. central sulcus lobe (red), 11. anterior frontal lobe (dark red). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 4.8 Normalization factors of the left and right hemispheres (LH and RH) of the 4 noninjured subjects. The relative cortical expansion model of Tallinen et al is g = 1 + 1.829t m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 4.9 Local cortical expansion (LCE) of each lobe on right hemisphere of each individual.

Figure 1 :

 1 Figure 1: Développement du cortex cérébral du foetus à l'adulte. Source : page d'accueil du laboratoire de Van Essen (http://brainvis.wustl.edu/wiki/index.php/Main_Page).

Figure 2 :

 2 Figure 2: Le cerveau est modélisé comme un solide élastique souple et une expansion tangentielle relative est imposée à la couche corticale comme illustré à gauche. La simulation à partir d'un cerveau foetal lisse montre une gyrification résultant d'une expansion tangentielle uniforme de la couche corticale. Extrait de [172].

Figure 3 :

 3 Figure 3: Composition du modèle combiné : le processus de contact, le processus élastique, le processus dynamique, le processus de croissance BLL et d'augmentation de volume. BLL est la longueur de l'axe le plus long du cerveau.

Figure 4 :

 4 Figure 4: Comparaison du développement du volume cérébral entre nos résultats de simulation, les résultats des simulations du modèle de Tallinen et al. [171,172] et les données de la littérature [5, 28, 72, 90].

Figure 5 :

 5 Figure 5: (a) Organigramme et outils utilisés pour générer un maillage tétraédrique à partir de la segmentation d'IRM ; (b) Diagramme correspondant à l'organigramme.

Figure 6 :

 6 Figure 6: Pipeline de l'application d'une carte d'expansion corticale spatio-temporelle [65] à un cerveau foetal de 22 semaines de gestation.

Figure 1 . 1 :

 11 Figure 1.1: Three major parts of the brain: cerebrum, brain stem and cerebellum. Extracted from https://www.myshepherdconnection.org.

Figure 1 . 2 :

 12 Figure 1.2: Gyri and sulci. Extracted from [22].

Figure 1 . 3 :

 13 Figure 1.3: Six major lobes. (a) Frontal lobe, temporal lobe, occipital lobe, parietal lobe, and insular lobe. Extracted from [1]. (b) Limbic lobe. Extracted from Wikipedia.

Figure 1 . 4 :

 14 Figure 1.4: Nine-region parcellation. Extracted from [189].

Figure 1 . 6 .

 16 There are four main hypotheses that dominate the speculation and study of cortical folding and are used to model the development of cortical folds: (1) external (Skull) constraint; (2) axonal tension; (3) differential tangential growth; (4) differential proliferation and patterned growth.

Figure 1 . 6 :

 16 Figure 1.6: Different cortical folding models of the emergence of gyrification. Extracted from [167].

Figure 1 . 7 :

 17 Figure 1.7: Components of stress in three dimensions. Extracted from Wikipedia.

Figure 1 . 8 :

 18 Figure 1.8: Mapping between the reference configuration and the deformed configuration.

Figure 1 . 9 :

 19 Figure 1.9: (a) Strain-stress relation: difference between linear elastic and hyperelastic. (b) Stress-strain curves for various hyperelastic material models. Extracted from Wikipedia.

Figure 3 .

 3 Figure 3. Elastic, viscoelastic and unstable behaviour in response to mechanical force. (a-c) Application of axial force (F) is shown for a hypothetical bar-shaped element of tissue. In all cases, the sum of forces on the bar is equal, denoted by black arrows in opposing directions. (a) Under pulling forces, the bar will stretch elastically (l* . 1), resulting in tensile stress (s . 0). In the case of either viscoelastic tissue or growing tissue, sustained tension may lead to permanent growth and relaxation of the original stress. (b) Conversely, under pushing forces, the bar will shorten elastically (l* , 1), resulting in compressive stress (s , 0) andpotentially-tissue shrinking to relieve stress. (c) Under sufficiently high compressive force, the same bar may buckle to a lower energy configuration. In the resulting fold, stress is tensile along the outer curvature region (red) and compressive along the inner curvature region (blue) due to bending. The neutral axis (green) represents the intermediate location where stresses due to bending are zero. (Online version in colour.)

Figure 4 . 4 Figure 1 . 10 :

 44110 Figure 4. Competing (or complementary) hypotheses for cortical folding. Dotted arrows (1 and 2) denote mechanisms not consistent with prior experiments. Additional experiments are needed to clarify the true mechanism(s) of cortical folding, which could include a combination of these hypotheses.

Figure 1 . 11 :

 111 Figure 1.11: Wilhelm His' mechanical analogy: Left. Folding of a rubber tube; Right. Folding of a gut tube during morphogenesis. Adapted from [86].

Figure 1 . 14 :

 114 Figure 1.14: Dynamic patterns of tangential growth in human brain development. (a) The regions of highest expansion migrate smoothly from the central sulcus and nearby regions 28-30 weeks (Top) to parietal 30-34 weeks (Middle) and finally to frontal and temporal regions 34-38 weeks (Bottom). (b) In healthy infants, areas of highest cortical expansion (red) are consistent with the trajectory from prenatal development. (c) Schematic illustrating the trajectory of the maximum growth region from primary motor, sensory, and visual cortices (labelled 'pre') to frontal, parietal and temporal lobes (labelled 'post'). Pre, prenatal; post, postnatal. Extracted from [65].

Figure 1 . 15 :

 115 Figure 1.15: Illustration of spatial discretization for 3 methods: (a) finite difference method (FDM), (b) finite element method (FEM), (c) finite volume method (FVM).

Figure 1 . 16 :

 116 Figure 1.16: The flowchart of the interaction between MATLAB and FEBio.

Figure 1 . 17 :

 117 Figure1.17: Brain mesh at the initial time (22 weeks of gestation) and after deformation. Each tetrahedron has its own color because of different material properties: the shear modulus is defined according to the distance from the top surface, the tetrahedron close to the surface has a smaller shear modulus.

Figure 1 . 18 :

 118 Figure 1.18: Diagram of the volume of an irregular tetrahedron after deformation.

Figure 2 . 2 :

 22 Figure 2.2: Curves of the cortical growth defined by different expansion models.

Figure 2 . 3 :

 23 Figure 2.3: Boolean operation is performed on an ellipsoid and a torus to form an ellipsoid with a fissure (which is composed of two connected semi-ellipsoids).

Figure 2 . 4 :

 24 Figure 2.4: Illustration of ellipsoid with smooth fissure.

Figure 2 . 6 :Figure 2 . 7 :Figure 2 . 8 :

 262728 Figure 2.6: Comparison of folding processes of ellipsoids with different mesh densities. Simulations are generated based on α t = 1.829t and H i = 0.042.

Figure 2 . 9 :

 29 Figure 2.9: a Evolution of 3D GI over time for surfaces of different expansion models. b Comparison of folding patterns on the reference ellipsoid at the same 3D GI based on different expansion models. In red frames: difference in folding patterns for linear and 1 st Gompertz models; in green frames: difference in folding patterns for 1 st and 2 nd Gompertz models. c Comparison of folding patterns on the reference ellipsoid at time 1.0 based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

Figure 2 . 10 :

 210 Figure 2.10: Comparison of the curvature for surfaces based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

Figure 2 .

 2 Figure 2.11: (a) Comparison of the Pearson correlation coefficients of curvatures for surfaces based on different expansion models. (b) Mean curvatures of vertices for surfaces of 3D GI of 2.3 (from left to right: linear, 1 st Gompertz, 2 nd Gompertz and logistic models). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

Figure 2 . 12 :

 212 Figure 2.12: Comparison of the Pearson correlation coefficients of distances to convex hulls for surfaces based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

22

 22

Figure 2 . 13 :

 213 Figure 2.13: Comparison of folding processes of the reference ellipsoid based on different global or regional expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042.

Figure 2 . 15 :

 215 Figure 2.15: Comparison of folding patterns on different geometries with the same initial cortical thickness. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 . 16 :

 216 Figure 2.16: Comparison of folding patterns on the reference ellipsoid for different initial cortical thicknesses. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 . 17 :

 217 Figure 2.17: Curvature over time of model for geometries with different initial cortical thicknesses and elongation ratios. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 .

 2 Figure 2.18: 3D GI over time of model for geometries with different initial cortical thicknesses and elongation ratios. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 . 19 :

 219 Figure 2.19: Probability density function (PDF) of sulcal depth over time for different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 . 20 :

 220 Figure 2.20: Probability density function (PDF) of sulcal depth over time for different initial cortical thicknesses (H i ). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 . 21 :

 221 Figure 2.21: Visualization of sulcal depth at time 1.0 for different initial cortical thicknesses.For each cortical thickness, the left column is the deformed reference ellipsoid and the right column is the corresponding initial ellipsoid. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 . 22 :

 222 Figure 2.22: Illustration of principal directions of curvatures, gradients of Fiedler vectors, and angles distributions for surfaces of sphere, ellipsoid of elongation ratio 1.50 and 2.25 at time 0.32, 0.55 and 0.79. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 . 23 :

 223 Figure 2.23: Kullback-Leibler divergence of uniformity of angular distribution for geometries with different elongation ratios and initial cortical thicknesses. Two scales are used here: the z-axis varies in each frame to allow a clear visualization of the behavior of KL divergence for a fixed initial cortical thickness (a)/time of model (b); the colorbar is the same across a row to show the evolution of KL divergence relative to uniform distribution. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t.

Figure 2 . 24 :

 224 Figure 2.24: Comparison of folding processes of the complete ellipsoid (above) and the corresponding ellipsoid composed of two connected semi-ellipsoids (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 .

 2 Figure 2.25: 1 st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 . 27 :

 227 Figure 2.27: Illustration of principal directions of curvatures, gradients of Fiedler vectors, and angles distributions for surfaces of complete ellipsoid and ellipsoid with a fissure at time 0.32, 0.55, 0.77 and 1.0. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042

Figure 2 .

 2 Figure 2.29: 1 st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 . 30 :

 230 Figure 2.30: Comparison of the 3D GI in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

Figure 2 .Figure 2 . 33 :

 2233 Figure 2.32: 1 st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042.

3 ChapterEarly

 3 Folding Patterns on a Human Brain Model with a Volume Correction 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 3.2 Modeling of brain longitudinal length and cortical thickness . . . . . . . . . . . . 65 3.2.1 Modeling the growth of brain longitudinal length (BLL) . . . . . . . . . . 65 3.2.2 Modeling the increase in cortical thickness (H) . . . . . . . . . . . . . . . 65 3.3 Quantitative methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 3.3.1 Local degree of folding . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 3.3.2 Spectral analysis of gyrification . . . . . . . . . . . . . . . . . . . . . . . . 67 3.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.4.1 Validation of brain volume with biological time . . . . . . . . . . . . . . . 69 3.4.2 Impact of cortical thickness onto morphology of brain . . . . . . . . . . . 70 3.5 Reconstruction and mesh generation of fetal brain . . . . . . . . . . . . . . . . . 73 3.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75

Figure 3 . 1 :

 31 Figure 3.1: Comparison of human brain longitudinal length growth models.

Figure 3 . 2 :

 32 Figure 3.2: Composition of the combined model: the contact process, the elastic process, the dynamic process, the BLL growth and volume increase process. BLL is the length of the longest axis of the brain.

Figure 3 . 3 :

 33 Figure 3.3: Process of spectral analysis of gyrification. Step 1: The computation of the eigenfunctions of the Laplace-Beltrami operator of the cortical mesh provides the decomposition basis for the spectrum analysis. Step 2: The curvature decomposition based on the eigenfunctions gives a series of coefficients that are the characteristics of curvature in the spectral domain. This step is reversible, which allows for spectral filtering of the curvature, and thus supports the spectral segmentation of cortical folding patterns. Step 3: The production of the band power spectrum of curvature from the squared coefficients depends on the band design choices detailed in Figure 3.4 a). Extracted from [68].

Figure 3 . 4 :

 34 Figure 3.4: Spectral banding. a) Spectral bands design: frequency intervals based on the hypothesis of doubling frequency. The spatial resolution achieved is illustrated by the relative size of the nodal domains of the last eigenfunction of the basis and the edges of the cortical mesh (median size brain of the database). Log-linear plot. min, median, max: brain of minimum, median and maximum size of the database. b) Spectral sizing: brain size, bandwidth and spatial resolution for the whole database. WL: wavelength, EV: eigenvalue/eigenvector. Extracted from [68].

Figure 3 . 5 :

 35 Figure 3.5: Comparison of the brain volume development between our simulation results, simulations results of Tallinen's model and literature data.

Figure 3 . 6 :

 36 Figure 3.6: Cortical morphology with increasing gestational weeks for different initial cortical thicknesses. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t.

Figure 3 . 7 :

 37 Figure 3.7: Average of the absolute values of mean curvatures computed on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t.

Figure 3 . 8 :

 38 Figure 3.8: 3D gyrification index computed on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t.

Figure 3 . 9 :

 39 Figure 3.9: Comparison of average of local surface ratio computed with the sphere of a radius 25 mm on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t.

Figure 3 . 10 :

 310 Figure 3.10: Folding power of B6 for different cortical thicknesses. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t.

Figure 3 .

 3 Figure 3.11: (a) Flowchart and tools used to generate a tetrahedral mesh from MRI segmentation; (b) corresponding diagram of tetrahedral mesh generation, the first two set of images are from ITK-SNAP, the others are from ParaView (https://www.paraview.org/).

Figure 3 .

 3 Figure 3.12: (a) Mesh for simulation (from Tallinen et al.); (b) Generated mesh in this section.

Figure 4 . 1 )

 41 , but it increases in the lateral temporal lobe over time (yellow in Figure4.1).

Figure 4 . 1 :

 41 Figure 4.1: Regions of highest cortical expansion change over time. Extracted from [65].

  et al.. We show the flow chart of the relationship derivation process in Figure 4.2.

Figure 4 . 2 :

 42 Figure 4.2: Flow chart of relationship derivation from the cortical expansion defined by Garcia et al.[START_REF] Kara E Garcia | Dynamic patterns of cortical expansion during folding of the preterm human brain[END_REF] to that defined in the brain growth model of Tallinen et al.[START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF].

Figure 4 . 3 :

 43 Figure 4.3: Pipeline of the application of a spatio-temporal cortical expansion map [65] on a fetal brain at 22 weeks of gestation.

Figure 4 .

 4 Figure 4.5 shows the cutting of the whole brain surface mesh into the semi-brain surface meshes and the closure of the semi-brain surface mesh boundaries. The cutting of the brain mesh into the semi-brain meshes is performed according to the surface coordinates of left and right semi-brain mesh, and the faces whose 3 vertices are on the boundary are removed.In order to close the mesh boundaries, we adopt a robust hole-filling algorithm for triangular mesh[START_REF] Zhao | A robust hole-filling algorithm for triangular mesh[END_REF]. This algorithm creates the new triangles start from the boundary vertex with the smallest angle, and then creates the new vertices from the bisector of the angle formed by the three vertices of the edge, as shown in Figure4.4. If the created vertex is close enough to every related boundary vertex (less than the given threshold), these vertices will be merged. The stop criterion of the algorithm is until the whole region is patched by all newly created triangles. Finally, all new vertices are locally smoothed, the smoothed coordinates of each new vertex correspond to the barycenter of its neighbors.

Figure 4 . 4 :

 44 Figure 4.4: Rules for creating triangles: a) θ i ≤ 75 • ; b) 75 • < θ i ≤ 135 • ; c) θ i > 135 • . Extracted from[START_REF] Zhao | A robust hole-filling algorithm for triangular mesh[END_REF].

Figure 4 . 5 :

 45 Figure 4.5: Cutting of the entire brain mesh into the semi-brain meshes and the closure of the semi-brain mesh boundaries.

Figure 4 . 6 :

 46 Figure 4.6: Application of relative cortical expansion provided by Garcia et al. on the atlas, on the registered semi-brain surface meshes, and on the surface of the semi-brain volumetric mesh by using the spherical demons registration and the Kd-Tree method.

Figure 4 . 7 :

 47 Figure 4.7: Brain lobes parcellation is performed manually on the surfaces of the left and right semi-brains and each surface is split into 11 lobes (including 10 cortical lobes and 1 connected lobe): 1. connected lobe/region (gray), 2. medial and posterior frontal lobe (dark blue), 3. occipital lobe (navy blue), 4. medial and posterior temporal lobe (blue), 5. anterior temporal lobe (bright blue), 6. parietal lobe (blue-green), 7. insular lobe (green), 8 and 9. outer ring of the limbic lobe [gyrus fornicatus encompassing: cingulate gyrus with subcallosal area (bright green) and parahippocampal gyrus (gold)], 10. central sulcus lobe (red), 11. anterior frontal lobe (dark red).

4 Figure 4 . 8 :

 448 Figure 4.8: Normalization factors of the left and right hemispheres (LH and RH) of the 4 noninjured subjects. The relative cortical expansion model of Tallinen et al is g = 1 + 1.829t m .

Figure 4 . 9 :

 49 Figure 4.9: Local cortical expansion (LCE) of each lobe on right hemisphere of each individual. For the 1 st subject, three gestational periods are 28-30 wk, 30-33 wk and 33-38 wk; for the 2 nd subject, three gestational periods are 27-31 wk, 31-33 wk and 33-37 wk; for the 3 rd subject, three gestational periods are 27-29 wk, 29-33 wk and 33-36 wk; for the 4 th subject, three gestational periods are 27-30 wk, 30-34 wk and 34-36 wk. For example, LCE(28 -33wk) of the 1 st subject is calculated by LCE(28 -30wk) × LCE(30 -33wk).

Figure 4 . 10 :

 410 Figure 4.10: Comparison of the normalized average cortical expansion of Garcia's data with the original expansion model (α(t m ) = 1.829t m and t m = 0.987e -e -0.134(t GA -29.433 ) ) defined in the brain growth model.

Figure 4 . 11 :

 411 Figure 4.11: Linear modeling of the normalized average cortical expansion (α(t m )) for each lobe of 4 individuals. 1 st column corresponds to the left hemisphere of 4 individual brains, 2 nd column corresponds to the right hemisphere of 4 individual brains.

Figure 4 . 14 :

 414 Figure 4.14: Comparison of the folding patterns during the human brain development and at the final simulation time for different growth mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 15 : 4 . 4 . 2 Quantitative analysis of folding patterns 4 . 4 . 2 . 1

 4154424421 Figure 4.15: Cortical tangential expansion of different postmenstrual ages displayed on simulated surfaces of the spatio-temporal growth. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 16 :Figure 4 . 17 :

 416417 Figure 4.16: Comparison of the average of the absolute values of dimensionless mean curvatures between different lobes of the left and right semi-brains of the spatio-temporal cortical expansion. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 19 :

 419 Figure 4.19: Comparison of the 3D GI of each principal lobe between the global linear and the spatio-temporal cortical expansion mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 20 :

 420 Figure 4.20: Relation of the 3D GI at the same time of model for each principal lobe between the global linear and the spatio-temporal cortical expansion mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 21 :

 421 Figure 4.21: Pearson correlation coefficient of the dimensionless mean curvatures for each principal lobe of the left and right semi-brains between the global linear and the spatio-temporal cortical growth mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 22 :

 422 Figure 4.22: Visualization of the Pearson correlation coefficients of the dimensionless mean curvatures between the global linear and the spatio-temporal cortical growth mechanisms at maximum 3D GI. Left: the left semi-brain; Right: the right semi-brain. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Figure 4 . 23 :

 423 Figure 4.23: Comparison of the primary folds orientation between the simulated surface of the spatio-temporal cortical growth and the cortical surface of a fetal brain atlas [69]. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm.

Titre:

  Modélisation et Caractérisation du Plissement Cortical Mots clés : Modèles biomécaniques, Plissement cérébral, Épaisseur corticale, Expansion corticale, Géométrie initiale, Orientation des plis. Résumé : L'AVC ischémique périnatal constitue la principale cause de paralysie cérébrale unilatérale ou d'épilepsie chez les enfants nés à terme. Cependant, les causes de ces incapacités observées restent floues. Dans ce contexte, la modélisation informatique est un outil puissant pour fournir une meilleure compréhension du processus de plissement cérébral précoce. Des études récentes basées sur la modélisation biomécanique ont montré que les forces mécaniques jouent un rôle crucial dans la formation des convolutions corticales. Cependant, l'effet des paramètres physiques dans ces modèles et la corrélation entre les résultats de la simulation et les faits biologiques restent flous. Dans cette thèse, en utilisant un modèle biomécanique de plissement cérébral, nous étudions l'effet de la croissance corticale, la géométrie initiale et l'épaisseur corticale initiale sur les patterns des plis corticaux. En outre, nous améliorons le modèle biomécanique en ajoutant un nouveau modèle de croissance de la longueur longitudinale du cerveau et un mécanisme d'expansion corticale différentielle spatiotemporelle au modèle. De plus, afin de quantifier la morphologie de surface des simulations, plusieurs descripteurs des plis sont utilisés tels que les mesures basées sur la courbure, l'indice de gyrification et la profondeur sulcale. Nous introduisons également une nouvelle approche pour mesurer l'orientation des plis à l'aide d'outils géométriques. Title: Modeling and Characterization of Cortical Folding Keywords: Biomechanical models, Brain folding, Cortical thickness, Cortical expansion, Initial geometry, Folds orientation.

1

  Développement du cortex cérébral du foetus à l'adulte. Source: page d'accueil du laboratoire de Van Essen (http://brainvis.wustl.edu/wiki/index.php/ Main_Page). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xv 2 Le cerveau est modélisé comme un solide élastique souple et une expansion tangentielle relative est imposée à la couche corticale comme illustré à gauche. La simulation à partir d'un cerveau foetal lisse montre une gyrification résultant d'une expansion tangentielle uniforme de la couche corticale. Extrait de [172]. . . . . . xvi 3 Composition du modèle combiné: le processus de contact, le processus élastique, le processus dynamique, le processus de croissance BLL et d'augmentation de volume. BLL est la longueur de l'axe le plus long du cerveau. . . . . . . . . . . . xviii 4

  Curves of the cortical growth defined by different expansion models. . . . . . . . 36 vi List of Figures 2.3 Boolean operation is performed on an ellipsoid and a torus to form an ellipsoid with a fissure (which is composed of two connected semi-ellipsoids). . . . . . . . . 2.4 Illustration of ellipsoid with smooth fissure. . . . . . . . . . . . . . . . . . . . . . 2.5 Eigenfunctions φ 1 , φ 2 and φ 3 (from left to right). Colormap: blue(negative), red/yellow(positive). Extracted from [112]. . . . . . . . . . . . . . . . . . . . . . . 2.6 Comparison of folding processes of ellipsoids with different mesh densities. Simulations are generated based on α t = 1.829t and H i = 0.042. . . . . . . . . . . . . 2.7 Comparison of the curvatures for the surfaces of different mesh densities (md). Left: dimensional curvature; Right: dimensionless curvature. Simulations are generated based on α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . 2.8 Comparison of 3D GI for the surfaces of different mesh densities (md). Simulations are generated based on α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . 2.9 a Evolution of 3D GI over time for surfaces of different expansion models. b Comparison of folding patterns on the reference ellipsoid at the same 3D GI based on different expansion models. In red frames: difference in folding patterns for linear and 1 st Gompertz models; in green frames: difference in folding patterns for 1 st and 2 nd Gompertz models. c Comparison of folding patterns on the reference ellipsoid at time 1.0 based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042. . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . . . 33 2.2 2.10 Comparison of the curvature for surfaces based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.11 (a) Comparison of the Pearson correlation coefficients of curvatures for surfaces based on different expansion models. (b) Mean curvatures of vertices for surfaces of 3D GI of 2.3 (from left to right: linear, 1 st Gompertz, 2 nd Gompertz and logistic models). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.12 Comparison of the Pearson correlation coefficients of distances to convex hulls for surfaces based on different expansion models. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and H i = 0.042. . . . . . . . . . . . . . . . . .

  and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Comparison of folding patterns on the reference ellipsoid for different initial cortical thicknesses. Simulations are generated based on mesh density 106 tetrahedra/cm 3 and α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 2.17 Curvature over time of model for geometries with different initial cortical thicknesses and elongation ratios. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . 50 2.18 3D GI over time of model for geometries with different initial cortical thicknesses and elongation ratios. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . 51 2.19 Probability density function (PDF) of sulcal depth over time for different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 2.20 Probability density function (PDF) of sulcal depth over time for different initial cortical thicknesses (H i ). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 2.21 Visualization of sulcal depth at time 1.0 for different initial cortical thicknesses. For each cortical thickness, the left column is the deformed reference ellipsoid and the right column is the corresponding initial ellipsoid. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 2.26 Comparison of the 3D GI in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 2.27 Illustration of principal directions of curvatures, gradients of Fiedler vectors, and Comparison of surface morphology for the ellipsoid composed of two connected semi-ellipsoids (above) and the corresponding left semi-ellipsoid (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.29 1 st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.30 Comparison of the 3D GI in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.31 Comparison of surface morphology for the complete ellipsoid (above) and the ellipsoid with a smooth fissure (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . 2.32 1 st column: the pink area is used to calculate the mean curvatures; 2 nd column: the comparison of the curvatures in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.33 Comparison of the 3D GI in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α

2.15 Comparison of folding patterns on different geometries with the same initial cortical thickness. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii 2.16 t = 1.829t. . . . . . . . . . . . . 53 2.22 Illustration of principal directions of curvatures, gradients of Fiedler vectors, and angles distributions for surfaces of sphere, ellipsoid of elongation ratio 1.50 and 2.25 at time 0.32, 0.55 and 0.79. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . 54 2.23 Kullback-Leibler divergence of uniformity of angular distribution for geometries with different elongation ratios and initial cortical thicknesses. Two scales are used here: the z-axis varies in each frame to allow a clear visualization of the behavior of KL divergence for a fixed initial cortical thickness (a)/time of model (b); the colorbar is the same across a row to show the evolution of KL divergence relative to uniform distribution. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 and α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . 55 2.24 Comparison of folding processes of the complete ellipsoid (above) and the corresponding ellipsoid composed of two connected semi-ellipsoids (below). Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 2.25 1 angles distributions for surfaces of complete ellipsoid and ellipsoid with a fissure at time 0.32, 0.55, 0.77 and 1.0. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and H i = 0.042 . . . . . . . . . . . . . . . . . . 57 viii 2.28 t = 1.829t and H i = 0.042. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1 Comparison of human brain longitudinal length growth models. . . . . . . . . . . 3.2 Composition of the combined model: the contact process, the elastic process, the dynamic process, the BLL growth and volume increase process. BLL is the length of the longest axis of the brain. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Process of spectral analysis of gyrification.

  . 3.5 Comparison of the brain volume development between our simulation results, simulations results of Tallinen's model and literature data. . . . . . . . . . . . . . Average of the absolute values of mean curvatures computed on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 3.8 3D gyrification index computed on each brain surface. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm

3.6 Cortical morphology with increasing gestational weeks for different initial cortical thicknesses. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and cortical growth α t = 1.829t. . . . . . . . . . . . . . . . . . . ix 3.7

  Comparison of the normalized average cortical expansion of Garcia's data with the original expansion model (α(t m ) = 1.829t m and t m = 0.987e -e -0.134(t GA -29.433 ) ) defined in the brain growth model. . . . . . . . . . . . . . . . . . . . . . . . . . . 4.11 Linear modeling of the normalized average cortical expansion (α(t m )) for each lobe of 4 individuals. 1 st column corresponds to the left hemisphere of 4 individual brains, 2 nd column corresponds to the right hemisphere of 4 individual brains. . . 4.12 Linear modeling of all the normalized average cortical expansion (α(t m )) of the 4 individuals for each lobe of left and right semi-brains. . . . . . . . . . . . . . . . 4.13 Cortical growth rate (the slope of the growth curve presented in 4.12) of each lobe of the left and right semi-brains. . . . . . . . . . . . . . . . . . . . . . . . . . 4.14 Comparison of the folding patterns during the human brain development and at the final simulation time for different growth mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.15 Cortical tangential expansion of different postmenstrual ages displayed on simulated surfaces of the spatio-temporal growth. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . 4.16 Comparison of the average of the absolute values of dimensionless mean curvatures between different lobes of the left and right semi-brains of the spatio-temporal cortical expansion. Simulations are generated based on mesh density 1.45 × 10 6 .21 Pearson correlation coefficient of the dimensionless mean curvatures for each principal lobe of the left and right semi-brains between the global linear and the spatio-temporal cortical growth mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . 4.22 Visualization of the Pearson correlation coefficients of the dimensionless mean curvatures between the global linear and the spatio-temporal cortical growth mechanisms at maximum 3D GI. Left: the left semi-brain; Right: the right semi-brain. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.23 Comparison of the primary folds orientation between the simulated surface of the Comparison between different FEM codes . . . . . . . . . . . . . . . . . . . . . . 2.1 Initial geometries with different elongation ratios . . . . . . . . . . . . . . . . . . 2.2 Kullback-Leibler divergence of angular uniformity . . . . . . . . . . . . . . . . . 3.1 Total folding power of B4, B5 and B6 . . . . . . . . . . . . . . . . . . . . . . . .

	List of Tables
	1.1
	xi xii

. . . . . . . . . . . . . . . . . . . . . . . 86 x 4.10 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . . . 4.17 Comparison of the 3D GI between different lobes of the left and right semi-brains of the spatio-temporal cortical expansion. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . 4.18 Comparison of the average of the absolute values of dimensionless mean curvatures of each principal lobe between the global linear and the spatio-temporal cortical expansion mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . 4.19 Comparison of the 3D GI of each principal lobe between the global linear and the spatio-temporal cortical expansion mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. 4.20 Relation of the 3D GI at the same time of model for each principal lobe between the global linear and the spatio-temporal cortical expansion mechanisms. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4spatio-temporal cortical growth and the cortical surface of a fetal brain atlas [69]. Simulations are generated based on mesh density 1.45 × 10 6 tetrahedra/cm 3 and initial cortical thickness 2.5 mm. . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1 Global Cortical Expansion for Right and Left Hemispheres . . . . . . . . . . . . xiii xiv

  Human brain structure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.1.2 Growth . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 1.1.3 Normative brain growth models . . . . . . . . . . . . . . . . . . . . . . . . 7 1.1.4 Cortical expansion and folding . . . . . . . . . . . . . . . . . . . . . . . . 8 1.1.5 Quantitative description of cortical folding . . . . . . . . . . . . . . . . . . 10 1.2 Biomechanical framework of cortical folding . . . . . . . . . . . . . . . . . . . . . 11 1.2.1 Concepts in biomechanics . . . . . . . . . . . . . . . . . . . . . . . . . .

	Chapter	1
	Human Brain Development and
	Cortical Folding Modeling
	1.1 Human brain development . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .	5
	1.1.1 . 12
	1.2.2 Biomechanical folding hypotheses and models . . . . . . . . . . . . . . . . 16
	1.3 Biomechanical modeling by Finite Element Method . . . . . . . . . . . . . . . . . 24
	1.3.1 Finite Element Method (FEM) . . . . . . . . . . . . . . . . . . . . . . . . 24
	1.3.2 Biomechanical modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
	1.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

Linear and 1st Gompertz Linear and 2nd Gompertz Linear and Logistic 1st Gompertz and 2nd Gompertz 1st Gompertz and Logistic 2nd Gompertz and Logistic

		1.05	Correlation of distances to convex hull
	Pearson correlation coefficient	0.70 0.75 0.80 0.85 0.90 0.95 1.00	3D GI 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0

  . In order to compare quantitatively the uniformity Comparison of the 3D GI in the pink surface area between different initial geometries. Simulations are generated based on mesh density 10 6 tetrahedra/cm 3 , α t = 1.829t and

		2.4 2.6	ellipsoid with a fissure complete ellipsoid Impact of initial geometry on GI
		2.2
	3D GI	1.6 1.8 2.0
		1.4
		1.2
		1.0
			0.0	0.2	0.4 Time of model 0.6	0.8	1.0
	Figure 2.26:

Table 2 .

 2

	2: Kullback-Leibler divergence of angular uniformity
		time 0.32 time 0.55 time 0.77 time 1.0
	Complete ellipsoid	0.0023	0.0021	0.0022	0.0025
	Ellipsoid with a fissure	0.0276	0.0044	0.0017	0.0023

Table 3 .

 3

	1: Total folding power of B4, B5 and B6
		GW 22 GW 29 GW 32 GW 34
	H i = 0.74 mm	3.76	73.11	325.59	854.95
	H i = 1.48 mm	3.76	59.26	177.48	366.06
	H i = 1.98 mm	3.76	50.89	125.71	272.16
	H i = 2.48 mm	3.76	39.06	94.39	194.27
	H i = 2.98 mm	3.76	33.39	73.80	152.26
	H i = 5.96 mm	3.76	15.48	19.49	29.01

  parametrizes time of model, t m ∈ [0, 1] corresponds to t GA ∈ [22weeks, adult]. α(t m ) controls the magnitude of expansion at t m relative to the initial state. y is the distance from the top surface in material coordinates. H is the cortical thickness at t m . The cortical layer has an areal growth by a factor of g 2 = 8[START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF] relative to the victitious stress-free state, thus the global linear cortical growth is defined as α(t m ) = ( √ 8 -1)t m = 1.829t m in this model. Since the time in the brain growth model

	reads as
	g(y, t m ) = 1 +	α(t m ) 1 + e 10( y H -1) ,	(4.4)
	where t m

Table 4 .

 4

	1: Global Cortical Expansion for Right and Left Hemispheres
	Noninjured1 28-30 wk 30-33 wk 33-38 wk 28-38 wk
	RH	1.30	1.35	1.67	2.93
	LH	1.26	1.37	1.70	2.93
	Noninjured2 27-31 wk 31-33 wk 33-37 wk 27-37 wk
	RH	1.56	1.24	1.42	2.75
	LH	1.51	1.22	1.48	2.73
	Noninjured3 27-29 wk 29-33 wk 33-36 wk 27-36 wk
	RH	1.22	1.47	1.44	2.58
	LH	1.21	1.44	1.50	2.61
	Noninjured4 27-30 wk 30-34 wk 34-36 wk 27-36 wk
	RH	1.46	1.58	1.22	2.81
	LH	1.38	1.57	1.19	2.58

  Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 5.2 Perspective . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 5.2.1 Geometric smoothing and mesh correction . . . . . . . . . . . . . . . . . . 98 5.2.2 Initial geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.2.3 Cortical thickness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.2.4 Model improvement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

	Chapter	5
	Conclusion and Perspective
	5.1

Acknowledgements

This work, being part of the project Biophysically-based computational modeling to explore impact of neonatal stroke on brain growth, was supported by the Fondation pour la Recherche Médicale (grant DIC20161236453), INSERM and Institut thanks to my supervisor, Prof. Bohi, a postdoctor supervised by Julien, always discussed and solved problems with us; Benjamin Fouquet, a research engineer supervised by François, is trying to help us optimize the model; Mireia Alenyà, a PhD student of Oscar, provided us help during the last two years of my PhD,

A Appendix

Implementation of Models

The biomechanical human cerebral cortical folding model [START_REF] Tallinen | Gyrification from constrained cortical expansion[END_REF][START_REF] Tallinen | On the growth and form of cortical convolutions[END_REF], the human brain longitudinal length growth model and the spatio-temporal cortical expansion model have been implemented in Python which are available at https://github.com/rousseau/BrainGrowth/. All scripts were coded in python 2.7, but they are compatible to python 3.7 (tested on the version 3.7.3).

Author Publications

Journal articles

Conference abstracts

1. Wang, X., Bohi, A., Al Harrach, M., Dinomais, M., Lefèvre, J., and Rousseau, F. (2019, June). On early brain folding patterns using biomechanical growth modeling. In OHBM 2019-Annual meeting of the Organization for Human Brain Mapping.