Dr Sylvain Brochard

Rodolphe Bailly

Chloé, Sarah Nathan Guillaume

Anne Zhengyang

Emma Zakaria

Prior to anything else, I would like to express my special

List of Tables

Résumé Étendu

La paralysie cérébrale (PC) est définie comme un ensemble de troubles permanents du développement du mouvement et de la posture, responsables de limitations d'activité, causés par des lésions non progressives au cours du développement du cerveau chez le foetus ou le nourrisson. [START_REF] Peter | A Report: The Definition and Classification of Cerebral Palsy April 2006[END_REF] Elle entraîne un ensemble de déformations musculo-squelettiques, comme l'équin, l'une des conséquences les plus courantes chez les patients atteints de PC. La déformation est causée par la spasticité des muscles et modifie la morphométrie musculaire, et peut entraîner une déformation osseuse par l'interaction os-muscle. La pratique clinique a signalé la déviation par rapport à la norme des torsions fémorales et tibiales ou de l'angle cervico-diaphysaire pendant la croissance des patients. Cependant, pour l'articulation de la cheville, les changements morphométriques osseux liés à l'équin sont moins étudiés. Peu d'informations sont disponibles. Pour une meilleure compréhension de la déformation en équin, l'analyse morphologique de l'articulation de la cheville, incluant le calcaneus, le talus et l'extrémité distale du tibia, est essentielle. L'objectif de cette thèse est d'étudier la déformation osseuse morphologique de l'articulation de la cheville dans la population PC. Avant de caractériser la déformation anormale de la forme de l'os, il faut répondre à une question : comment définir la forme normale de l'os ? La modélisation de la normalité fournira les informations anatomiques de la population de TD (typique en développement). La modélisation statistique de la forme (SSM) est un outil pour comprendre les formes en utilisant des méthodes statistiques, y compris la forme moyenne et la variation de la forme. Il a été largement utilisé dans le domaine médical pour analyser les structures anatomiques. Dans [START_REF] Lenz | Statistical Shape Modeling of the Talocrural Joint Using a Hybrid Multi-Articulation Joint Approach[END_REF] et [START_REF] Krähenbühl | Morphologic Analysis of the Subtalar Joint Using Statistical Shape Modeling[END_REF], les auteurs ont analysé la morphométrie de l'articulation tibiotalaire et de l'articulation subtalaire de l'adulte avec l'outil SSM. Dans le cas des données pédiatriques, l'analyse morphologique reste un domaine ouvert. La première contribution de cette thèse est d'analyser la forme de l'os de la cheville d'enfants avec une approche SSM basée sur les voxels. Grâce à l'analyse de la forme, un modèle moyen de l'IRM de l'articulation de la cheville des enfants, également appelé atlas, est construit. La variation de la forme chez les enfants en bonne santé est analysée, et la corrélation avec les mesures cliniques est étudiée, ce qui fournit des informations sur la variation anatomique des os d'intérêt et une explication qualitative du modèle.

Une fois la normalité définie, les questions suivantes se posent : existe-t-il une différence morphologique entre les os de la cheville du TD et du CP ? Si oui, où se situent ces différences ? Une anomalie est une observation qui s'écarte considérablement de la normalité [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]. Dans le contexte de l'analyse des formes, une anomalie peut être définie par des paramètres morphologiques locaux qui s'écartent xvi considérablement de la distribution normale de la forme moyenne. Les trois principales classes d'algorithmes permettant de détecter les anomalies sont les modèles probabilistes, la classification à une classe et les modèles de reconstruction [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]. Dans le contexte considéré dans ce travail, un élément clé réside dans la petite taille de l'échantillon qui détermine le choix des méthodes de détection appropriées. Deux techniques complémentaires de détection d'anomalies, l'une basée sur une approche statistique et l'autre sur un critère de reconstruction, sont agrégées pour combiner les informations statistiques et biomécaniques. Le pipeline de la méthode proposée est illustré dans la Figure 1. En raison de la variabilité des réponses individuelles au traitement conçu sur la base de l'analyse de groupe, une approche personnalisée des soins médicaux est actuellement nécessaire [START_REF] Damiano | Meaningfulness of Mean Group Results for Determining the Optimal Motor Rehabilitation Program for an Individual Child with Cerebral Palsy[END_REF]. Une approche personnalisée pour comprendre la variation inter-sujet devient nécessaire pour une réadaptation orientée vers le patient. Par conséquent, le troisième objectif de cette thèse est d'analyser la forme des os à la fois au niveau du groupe et au niveau personnel, afin de déterminer si le modèle de déformation personnel correspond au modèle de déformation de la population.

Une étude récente [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF] a montré l'importance de l'évaluation et de la validation de ces outils dans les applications cliniques. Les résultats quantitatifs et qualitatifs montrent que les outils SSM ont différents niveaux de cohérence et différentes capacités à capturer la variabilité au niveau de la population. Ce qui ressort de cette étude est la nécessité de comparer les résultats obtenus à l'aide de plusieurs méthodes d'analyse de forme. Nous adoptons une approche basée sur les surfaces et une approche basée sur les voxels. Les deux méthodes capturent des régions de déformation similaires au niveau du groupe (voir la Figure 2).

Cette thèse est organisée comme suit :

Le chapitre 1 présente le contexte clinique, notamment la PC et la déformation en équin spastique, l'anatomie de l'articulation de la cheville et la déformation musculo-squelettique dans la PC. La technologie pertinente, notamment le SSM et la détection des anomalies, est également présentée. Le chapitre 2 modélise la forme du calcaneus, du talus et du tibia. Tout d'abord, l'acquisition et le prétraitement des données utilisées dans cette thèse sont présentés. Ensuite, la forme des os de la population TD est modélisée avec une méthode SSM basée sur les voxels et la variation de la forme parmi la population TD est présentée. En outre, la corrélation entre les paramètres de forme du modèle et les coefficients cliniques est analysée pour évaluer et comprendre le modèle construit.

Le chapitre 3 présente l'analyse des anomalies morphologiques de l'articulation de la cheville. Dans ce chapitre, nous analysons les déformations osseuses anormales de manière multi-résolution : à l'échelle globale, par région et par voxel. Diverses métriques quantifiables pour mesurer la morphométrie sont introduites pour l'analyse globale et régionale. Nous avons proposé un pipeline d'analyse des anomalies adapté à la taille limitée de l'ensemble de données pour l'analyse basée sur les voxels. Des expériences sur des données IRM synthétiques et réelles sont présentées afin d'évaluer la méthode proposée et d'étudier les différences de forme des os entre les populations.

Le chapitre 4 présente une analyse personnalisée complète de la forme des os des articulations de la cheville. Tout d'abord, l'analyse de groupe est réalisée avec deux méthodes SSM différentes afin de comparer les résultats obtenus avec différentes méthodes. Ensuite, l'analyse personnalisée de la forme est exécutée avec une méthode SSM basée sur les voxels, qui permet de comprendre au niveau du patient.

Le chapitre 5 conclut cette thèse et présente les perspectives de travaux futurs, notamment l'analyse des anomalies cinématiques basée sur l'IRM dynamique, l'interprétation des résultats de l'analyse des anomalies et le recalage d'images basé sur l'apprentissage profond.

Motivation

Cerebral palsy (CP) is a group of disorders that affect a persons ability to move and maintain balance and posture. It causes a set of musculoskeletal deformities, such as the equinus, one of the most common consequences that occurs to patients with CP. The deformity is caused by the spasticity of muscles and changes the muscle morphometry, and may lead to bone deformation through the bone-muscle interaction. The clinical practice has reported the deviation from the norm of the femoral and tibial [START_REF] King | Torsional Problems in Cerebral Palsy[END_REF], [START_REF] Hyeong | Tibial Torsion in Cerebral Palsy: Validity and Reliability of Measurement[END_REF] or the cervico-diaphysial angle [START_REF] Fernandez | Adductor Myotomy in Cerebral Palsy: Uni or Bilateral[END_REF] during the growth of the patients. However for the ankle joint, the bone morphometry changes related to the equinus are less studied. Little information is available. For a better understanding of the equinus deformity, the morphological analysis of the ankle joint, including the calcaneus, talus and distal end of the tibia, is essential.

Before the characterization of abnormal bone shape deformation, one question needs to be answered: how to define the normal bone shape. Modeling the normality with Statistical Shape Modeling (SSM) can provide the anatomical information of the typical developing (TD) population. SSM is a powerful tool for understanding shapes using statistical methods, including the mean shape and the shape variation. It has been widely used in the medical field to analyze anatomical structures. With the help of SSM, the normality of ankle joint bones can be accessed.

Once the normality is defined, then the following questions will be asked: Is there any morphological difference existing between TD's and CP's ankle joint bones? If yes, where are the differences located? An anomaly is an observation that deviates considerably from the normality [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]. In the context of shape analysis, an anomaly can be defined by local morphological parameters that deviate significantly from the normal distribution to the mean shape. Anomaly detection is a technology that aims to localize abnormal samples based on probability, which enables us to detect the possible shape difference between two groups.

Chapter 1. General introduction

Due to the variability of individual responses to the treatment designed based on the group analysis, a personalized approach to medical care is currently required [START_REF] Damiano | Meaningfulness of Mean Group Results for Determining the Optimal Motor Rehabilitation Program for an Individual Child with Cerebral Palsy[END_REF]. A personalized approach provides the possibility for patient-oriented rehabilitation. Therefore the third aim of this thesis is to analyze the bone shape at both group level and personal level, to figure out if the personal deformation pattern corresponds to the pattern of the population.

Thesis overview

The ultimate objective of this thesis is to investigate the morphological bone deformation of ankle joint in children with CP. The main contributions of this thesis are:

• A shape analysis of the ankle joint bone of children with a voxel-based SSM approach and the construction of atlas of children ankle joint.

• A multi-resolution children ankle joint morphometry analysis to compare the morphological difference between the TD population and the CP population.

A small-dataset-adapted pipeline for voxel-based anomaly analysis is developped to localize the abnormal bone deformation in voxel scale.

• A comprehensive personalized shape analysis allowing a morphometry analysis in the subject level.

• A comparison between different SSM tools in order to confirm the robustness and the reliability of the obtained results.

Thesis organisation

This thesis is organized as follows:

Chapter 2 briefly introduces the clinical background, including the CP and the spastic equinus deformity, the ankle joint anatomy and the musculoskeletal deformation in CP. The relevant technology, including SSM and anomaly detection, is also presented.

Chapter 3 models the shape of the calcaneus, talus and tibia. Firstly, we introduce the acquisition and the preprocessing of the data used in this thesis. Then the bone shape of the TD population is modeled with a voxel-based SSM method and the shape variation analysis among the TD population is presented. In addition, to evaluate and understand the constructed model, the correlation analysis between the shape parameters of the model and the clinical coefficients is performed.

Chapter 4 presents the ankle joint morphological anomaly analysis. In this chapter, we analyze abnormal bone deformation in a multi-resolution way: global-scale, region-based, and voxel-based. For global and regional analysis, some quantificational metrics to measure morphometry are introduced. For voxel-based analysis, we propose an anomaly analysis pipeline adapted to the limited-size dataset. Experiments on our MRI data reveals the possible bone shape differences between populations.

Chapter 5 presents a comprehensive personalized bone shape analysis of ankle joints. Firstly, the group analysis is performed with two different SSM methods to compare the results obtained with different approaches. Secondly, the personalized shape analysis is executed with a voxel-based SSM method, which provides the understanding of the pathology at the patient level.

Chapter 6 concludes this thesis and presents the perspectives for future works, including the In vivo kinematic anomaly analysis based on dynamic MRI, interpretation of anomaly analysis results, and the deep-learning-based registration.

Clinical context

Cerebral palsy and spastic equinus deformity

Cerebral palsy (CP) is the most common physical disability in childhood, affecting 1-4 out of every 1000 newborns and 17 million people worldwide are affected by CP [START_REF] Wimalasundera | Cerebral Palsy[END_REF], [START_REF] Graham | Cerebral Palsy[END_REF]. Cerebral palsy is defined as a set of "permanent disorders in the development of movement and posture, responsible for activity limitations, caused by non-progressive damage during brain development in the fetus or infant." [START_REF] Peter | A Report: The Definition and Classification of Cerebral Palsy April 2006[END_REF] This non-progressive neurological damage results in a set of movement and posture disorders.

According to the type of motor dysfunction and activity restriction, three classes of CP are distinguished: spastic, ataxic and dyskinetic [START_REF] Graham | Cerebral Palsy[END_REF]. The mixture of three types also exists. These classifications reflect the areas of the brain that are damaged. The spasticity is the most common symptom of CP that appears in 70% of patients. Foot equinus is the most common deformity in children with spastic cerebral palsy [START_REF] Metaxiotis | Repeated Botulinum Toxin A Injections in the Treatment of Spastic Equinus Foot[END_REF], which is typically defined as the inability to dorsiflex the foot above plantigrade, with the hindfoot in neutral position and the knee in extended position [START_REF] Kay | Outcome of Gastrocnemius Recession and Tendo-Achilles Lengthening in Ambulatory Children with Cerebral Palsy[END_REF], [START_REF] Shore | Surgical Correction of Equinus Deformity in Children with Cerebral Palsy: A Systematic Review[END_REF]. In children with bilateral cerebral palsy, equinus prevalence is 83.3% and tends to increase with age [START_REF] Horsch | Prevalence and Classification of Equinus Foot in Bilateral Spastic Cerebral Palsy[END_REF]. It manifests poor muscle control and weakness around the ankle and foot, leading to abnormal gait patterns and bone deformations during growth. The ankle joint consists of 4 bones: calcaneus, talus, tibia and fibula, as illustrated in Fig. 2.1. The articulation between calcaneus and talus is subtalar joint, also called talocalcaneal joint. This joint enables the inversion and exversion of the foot. The talocrural joint is the articulation between bone of the leg (tibia and fibula) and talus, permitting dorsiflexion and plantarflexion of the foot. It consists of tibiotalar joint and tibiofibula joint. In this thesis, we focus on calcaneus, talus and tibia, also the subtalar joint and tibiotalar joint.

Anatomy of ankle

Bones of interest

The calcaneus, as known as the heel bone, is the largest bone of human foot. The anatomy of calcaneus is illustrated in Figure 2.2.

The talus connects with two bones of lower leg, the tibia and the fibula. At the foot end, it articulates with the calcaneus to form the subtalar joint. The anatomy of talus is illustrated in Figure 2.3. It can be divided into three parts: the head, the neck and the body.

The distal end of tibia is the end contributing to the ankle joint. The distal end has five surfaces: anterior, posterior, medial, lateral and distal. The lateral surface exhibits a triangular notch which attaches to the fibula. The distal articular surface connects to the talar dome. The medial side is a projection called medial malleolus. The anatomy of tibia is illustrated in Figure 2.4.

Tendons and ligaments of ankle joint

Tendons are fibrous tissues and play a role of connection between bone and muscle. It is able to transmit the mechanical forces of muscle contraction to the skeletal system and withstand significant amounts of tension. An illustration of foot tendons is shown in Figure 2.5.

Ligaments are fibrous tissues that are similar to tendons but they connect bone to bone. They either attach to the surface of the bone or fuse with the outer layer of the joint capsule to strengthen the stability of the joint. The lateral, medial and posterior view of ankle and foot ligaments are presented in Figure 2.6.

Musculoskeletal deformation in CP

A common cause of equinus deformity is the imbalance between the plantar flexor muscles and the ankle's dorsal flexor muscles. It is possible that this imbalance is caused by spasticity, which can affect the entire triceps surae muscle complex or perhaps the isolated gastrocnemius muscle. Literature reports that fixed equinus deformity is not initially present but develops during adolescence. In [2], experiments revealed the slower muscles and bones growth in CP than their TD siblings. At first, the equinus is presented as dynamic contracture of muscles and impact to gait pattern. Over the time dynamic tightness can develop into fixed deformity occurs in soft tissue, articulars and bones. In [START_REF] Morrell | Progressive Bone and Joint Abnormalities of the Spine and Lower Extremities in Cerebral Palsy[END_REF], authors reported the correlation of rocker-bottom deformity and subluxation of the talonavicular joint with equinus deformity. In [START_REF] Schless | Combining Muscle Morphology and Neuromotor Symptoms to Explain Abnormal Gait at the Ankle Joint Level in Cerebral Palsy[END_REF], authors reported the morphological muscle changes related to CP. The muscles of children and adolescents with spastic CP were significantly smaller than those with typically developing muscles. Hence little information is available to the bone morphological changes regarding to the ankle joint.

The muscle quality has influence to the bone growth. It has been demonstrated that several myokines are involved in bone metabolism, whose secretion appears to be related to exercise and muscle contraction [START_REF] Bente | Adolph Distinguished Lecture: Muscle as an Endocrine Organ: IL-6 and Other Myokines[END_REF], [START_REF] Walsh | Adipokines, Myokines and Cardiovascular Disease[END_REF]. For instance, IL-6 acts on osteoclastic activity and activates bone resorption. Muscle also seems to have a paracrine action on bone metabolism. Studies have shown that muscle hypertrophy and bone anabolism are coupled. A powerful bone anabolic stimulus is also provided by exercise-induced alterations in the plasma membrane of myofibers via FGF2.

Therefore, one hypothesis can be proposed that the morphological differences may exist between children with equinus deformity and the typical developing children. Investigating this difference will enable a better understanding of CP and a suitable rehabilitation program.

Magnetic resonance imaging (MRI) in the study of ankle joint

MRI is a medical imaging technique used in radiology to visualize anatomy. Unlike the technologies working with ionizing radiation, such as CT, MRI scanners use strong magnetic fields, magnetic field gradients, and radio waves to generate images of the organs in the body. The basic principles of MRI are introduced in [START_REF] Vijay | Magnetic Resonance Imaging: Principles and Techniques: Lessons for Clinicians[END_REF]. The hydrogen nucleus possesses the property of spin that results in the rotation or process of the nucleus when applying a strong external static magnetic field B 0 . The velocity of the rotation and the process is the Larmor frequency ω 0 , and proportional to the field strength ω 0 = γ 0 B 0 . A second radio frequency (RF) magnetic field B 1 is applied to the nuclei under the B 0 . The RF energy is usually applied in short pulses at the resonant frequency which equals ω 0 . The RF pulses excite the nucleus and cause an energy transition. The absorption and the relaxation of the energy induces an electromagnetic signal that can be detected by a suitably tuned coil of wire. Different tissues produce different electromagnetic signals. The location information is provided by using the magnetic field gradients. The non-invasive and ionizing radiation-free nature of MRI plays an irreplaceable role in pediatric studies, as in this thesis.

Static MRI

Static MRI techniques have been used for accurate diagnosis of ankle joint disorders and ankle joint morphological analysis. In [START_REF] Elstob | MRI Ankle and Subtalar Characteristics in Haemochromatosis Arthropathy: A CaseControl Study[END_REF], Elstob et al. examined the MRI features of the ankle related to genetic hemochromatosis, which are numerous and large cysts on ankle MRI, large osteophytes and the presence of extensive fullthickness cartilage loss in the ankle joint and middle subtalar articulation. In [START_REF] Schütz | Biochemical Cartilage Alteration and Unexpected Signal Recovery in T2* Mapping Observed in Ankle Joints with Mobile MRI during a Transcontinental Multistage Footrace over 4486 Km[END_REF], the effect of long-distance running on the ankle cartilage has been studied with the help of T2w MR images scanned by a 1.5T MRI on a mobile unit. Also, MRI is used to diagnose and grade the osteoarthritis of the ankle [START_REF] Aboelmagd | The Norwich Osteoarthritis of the Ankle MRI Score (NOAMS): A Reliability Study[END_REF].

MRI is not only applied in bone analysis but also the soft tissues such as the ligaments and tendons in ankle joint studies. In [START_REF] Chan | Medial Joint Line Bone Bruising at MRI Complicating Acute Ankle Inversion Injury: What Is Its Clinical Significance?[END_REF], authors studied the ligament damage and the medial joint line bone bruising following an acute ankle inversion injury with MRI examinations. In [START_REF] Siegler | Mechanics of the Ankle and Subtalar Joints Revealed through a 3D Quasi-Static Stress MRI Technique[END_REF] and [START_REF] Ringleb | The Effect of Ankle Ligament Damage and Surgical Reconstructions on the Mechanics of the Ankle and Subtalar Joints Revealed by Three-Dimensional Stress MRI[END_REF], the mechanics of the ankle joints have been studied with 3D stress MRI, and the effect of the ligament damage. MRI can be used to evaluate the quality of ankle ligaments, such as anterior talofibular ligament, lateral fibulotalocalcaneal ligament complex, and calcaneofibular ligament in chronic lateral ankle instability [START_REF] Wiewiorski | Delayed Gadolinium-Enhanced MRI of Cartilage of the Ankle Joint: Results after Autologous Matrix-Induced Chondrogenesis (AMIC)-Aided Reconstruction of Osteochondral Lesions of the Talus[END_REF], [START_REF] Nakasa | MRI Signal Intensity Ratio Reflects the Quality of the Anterior Talofibular and Calcaneofibular Ligaments in Patients with Chronic Lateral Ankle Instability[END_REF], [START_REF] Nakasa | MRI Appearance of the Lateral Fibulotalocalcaneal Ligament Complex Injury in the Patients with Chronic Lateral Ankle Instability[END_REF], [START_REF] Choi | The Influence of Suture-Tape Augmentation on Biological Healing of the Anterior Talofibular Ligament in Chronic Ankle Instability: A Quantitative Analysis Using MRI[END_REF]. It is also used for the detection of early architectural changes in cartilage and subchondral bone for chronic lateral ankle instability patients [START_REF] Tao | Using Radiomics to Detect Subtle Architecture Changes of Cartilage and Subchondral Bone in Chronic Lateral Ankle Instability Patients Based on MRI PD-FS Images[END_REF]. In [START_REF] Zhang | Risk Factors for Chronic Ankle Instability after First Episode of Lateral Ankle Sprain: A Retrospective Analysis of 362 Cases[END_REF], authors identified the risk factor for chronic ankle instability following an acute lateral ankle sprain with MRI, including age, BMI, posterior talofibular ligament injury, large bone marrow lesion of the talus and Grade 2 effusion of the tibiotalar joint. In [START_REF] Ishimatsu | Ankle Bone Morphology Affects the Size of Non-Trauma Related Osteochondral Lesions of the Talus in Skeletally Immature Children[END_REF], Ishimatus et al. reported the relationship between ankle morphology and the size of non-trauma-related osteochondral lesions of the talus.

Static MRI is a mature technique for both clinic and research proposes in bone analysis. However, one challenge in the acquisition is that it requires the subject to be examined to remain stationary during the scan in order to obtain an image with high resolution. The acquisition lasts normally for minutes. It may be a challenge for subjects.

Dynamic MRI

For physiologically better represent the dynamic musculoskeletal system, researchers have developed the dynamic MRI technique, which enables the in vivo analysis of joint mechanics. It has been successfully applied in the quantification of patellofemoral cartilage contact kinematics [START_REF] Bhushan | A Methodology to Accurately Quantify Patellofemoral Cartilage Contact Kinematics by Combining 3D Image Shape Registration and Cine-PC MRI Velocity Data[END_REF], in vivo quantification of ankle joint space width [START_REF] Makki | 4D in Vivo Quantification of Ankle Joint Space Width Using Dynamic MRI[END_REF] and ankle joint kinematic analysis and motion reconstruction [START_REF] Makki | In Vivo Ankle Joint Kinematics from Dynamic Magnetic Resonance Imaging Using a Registration-Based Framework[END_REF]. Some technologies are based on the motion-triggered principle, i.e. data are continuously collected over many cycles and then retrospectively sorted with a synchronization trigger [START_REF] Sheehan | Using Cine Phase Contrast Magnetic Resonance Imaging to Non-Invasively Study in Vivo Knee Dynamics[END_REF]. The cine MRI [START_REF] Langner | Cine MRI: A New Approach to the Diagnosis of Scapholunate Dissociation[END_REF] and cine-PC (cine Phase Contrast) MRI [START_REF] Sheehan | 3D Musculoskeletal Kinematics Using Dynamic MRI[END_REF], or fast-PC MRI, belong to this category. The limitation of this technique is the requirement of a long acquisition time. In [START_REF] Sheehan | In Vivo Talocrural and Subtalar Kinematics: A Non-Invasive 3D Dynamic MRI Study[END_REF], authors analyzed the foot joint kinematics, including talocrural and subtalar joints, with cine-PC MRI. Also in [START_REF] Rebmann | Precise 3D Skeletal Kinematics Using Fast Phase Contrast Magnetic Resonance Imaging[END_REF], knee joint kinematics has been studied with cine-PC MRI.

Another category of dynamic MR imaging is real-time sequences. Contrary to the motion-triggered technology, it achieves the imaging of examined anatomy trajectory in one motion cycle by exploiting spatio-temporal redundancy. In [START_REF] Clarke | A Non-Invasive, 3D, Dynamic MRI Method for Measuring Muscle Moment Arms in Vivo: Demonstration in the Human Ankle Joint and Achilles Tendon[END_REF], the ultra-fast turbo gradient echo sequences have been used in the muscle moment arm measurement method, which is demonstrated in the human ankle joint and Achilles tendon. The dynamic images with 20 frames are acquired in 50 seconds. Real-time Fast Field Echo (FFE) sequence is another technique of real-time MRI. As presented in [START_REF] Makki | In Vivo Ankle Joint Kinematics from Dynamic Magnetic Resonance Imaging Using a Registration-Based Framework[END_REF], the dorsi-plantar flexion of the ankle joint is scanned in 18 seconds to acquire a dynamic image of 15 frames and the resolution is 0.56 × 0.56 × 9mm. In

Statistical shape modeling (SSM)

SSM is a tool to analyze the geometric properties of a given set of shapes. It uses statistical methods to compute the quantitative representation of shapes. SSM describes the mean shape and the principle modes of variation of the shapes, therefore it facilitates the population-wise shape comparison and structure reconstruction and has been widely applied in the medical field. In neuroscience, it enables the analysis of the relation shape between brain morphometry and the brain disorders [START_REF] Joseph | Three-Dimensional Surface Deformation-Based Shape Analysis of Hippocampus and Caudate Nucleus in Children with Fetal Alcohol Spectrum Disorders: Analysis of Hippocampus and Caudate Nucleus in Children[END_REF], [START_REF] Younes | Regionally Selective Atrophy of Subcortical Structures in Prodromal HD as Revealed by Statistical Shape Analysis: Atrophy in HD Revealed by Shape Analysis[END_REF], [START_REF] Jikke | Structural Assessment of Thalamus Morphology in Brain Disorders: A Review and Recommendation of Thalamic Nucleus Segmentation and Shape Analysis[END_REF]. In skeletal applications, SSM helps the 3D bone reconstruction, pathology diagnosis, implant design, and surgery planning. A review of SSM and its application in bone shape analysis is available in [START_REF] Sarkalkan | Statistical Shape and Appearance Models of Bones[END_REF].

A shape is the external boundary, outline, or external surface of an object. Let A 1 , A 2 , . . . , A N be the set of N shapes. The correspondence between these shapes can be computed by registration technique and a deformation field D i is obtained which represents the morphological information of A i . The key technology to build the correspondence between shapes is therefore the registration, a technology to transform images from their original space to another space.

Image registration

Image registration is a process to transform several different data into one coordinate system and is widely applied in different domains including medical image processing [START_REF] Noblet | Recalage d'images médicales[END_REF]. The image to which the other images are transformed is considered as the reference I re f and the other images are called as moving images I mov . The objective of registration is therefore to find the transformation T that will maximize the similarity between the reference image and the transformed moving image I mov • T ∼ I re f . To this end, the registration procedure consists of 4 elements: features of the image, similarity or distance between images, transformation, and optimization.

According to the application context, medical image registration can be divided into intra-individual registration and inter-individual registration. In intra-individual registration, the reference image and the moving image belong to the same subject, and in inter-individual registration, where 2 images are from different subjects.

Image features and inter-image distance

The first step of the registration is to extract the features of images, which is used to measure the distance between 2 images during the optimization and to guide the registration. Based on the features extracted, it can be distinguished into 2 categories: geometric approaches and iconic approaches.

The geometric approaches work with some geometric features such as points (also called landmarks) [START_REF] Thirion | New Feature Points Based on Geometric Invariants for 3D Image Registration[END_REF], curves [START_REF] Declerck | Automatic Retrieval of Anatomical Structures in 3D Medical Images[END_REF] and surfaces [START_REF] Malandain | Topological Segmentation of Discrete Surfaces[END_REF]. The features can be defined by experts according to anatomy or automatically detected. Denote the G n (n = 1, 2, . . . , N) is the geometric features, the estimation of transformation is performed by minimizing the distance map of G n (n = 1, 2, . . . , N) between I re f and I mov :

T = arg min N ∑ n=1 ∥G n (I re f ) -G n (I mov )∥ 2 (2.1)
such as ICP (Iterative Closest Point) [START_REF] Besl | Method for Registration of 3-D Shapes[END_REF] algorithm.

The geometric methods are often based on the anatomy which enables the medical explanation and discrimination. However, the accuracy of automatic feature extraction needs to be taken into consideration and the error of feature detection may lead to wrong correspondence. Also, the local registration accuracy is related to the distance to the pre-defined geometric feature. It only guarantees accuracy near the region of features.

Instead of using geometric features, the iconic approaches focus entirely on the features based on the voxel intensity of the image and optimize the transformation by maximizing the similarity criterion between I re f and I mov .

T = arg max sim(I re f , I mov • T) (2.2)
Based on the hypothesis of relationship between I re f and I mov , the similarity measures can be divided into 4 categories [START_REF] Roche | Unifying Maximum Likelihood Approaches in Medical Image Registration[END_REF]:

• Identical relationship: this category assumes that the density of voxel in correspondence between I re f and I mov is identical. Therefore the similarity can be defined as the mean squares difference. This category of methods is only applied in mono-modality registration.

• Affine relationship: since the identical relationship is rarely met in practice, the affine relationship is often used to measure the linear similarity, such as the linear correlation coefficient.

• Functional relationship: this category of relationship is less strict. It assumes that the relationship can be described with a continuous function and not limited to a linear relationship, thus suiting the multi-modality cases. Correlation ratio [START_REF] Roche | The Correlation Ratio as a New Similarity Measure for Multimodal Image Registration[END_REF] and its normalized version belong to this category.

• Statistical relationship: such as mutual information [START_REF] Wells | Multi-Modal Volume Registration by Maximization of Mutual Information[END_REF] and normalized mutual information, it is based on information theory and calculates the statistical dependence with the joint histogram.

Linear registration

In linear registration, the transformation can be described with a linear model:

I re f = Rx + t (2.3)
For 3D application, R is a 3 × 3 rotation matrix, and t is a 3D translation vector. Thus the transformation matrix T can be defined as a 4 × 4 matrix:

T = R t 0 1 (2.4)
The linear transformation consists mainly of rigid transformation, similarity, affine transformation, and projective transformation.

• Rigid transformation: the rigid transformation contains only rotation and translation. The number of DOFs is 6 for 3D images. This type of transformation only suits the intra-subject registration, i.e. I re f and I mov are from the same subject.

• Similarity: in addition to translations and rotations, the similarity possesses an isotropic scale factor. Thus the number of DOFs for 3D images is 7. It preserves angles, parallelism, and length ratio. This type of transformation is also applied in intra-subject registration of different modalities.

• Affine transformation: the affine transformation allows the rotation, translation, reflection, scaling, and shearing. It is able to preserve the parallelism.

The number of DOFs is 12 in 3D. It is often used in simple inter-subject situations but for more complicated cases it is not able to capture the anatomical variability.

• Projective transformation: the projective transformation does not preserve parallelism, length, and angle. But it still preserves collinearity and incidence. Among the applications for this type of transformation are medical imaging for which a point source generates radiation that interacts with the object being imaged in a plane to produce a projection image of the object, as well as virtual reality applications.

projective FIGURE 2.9: Linear transformation applied on 2D square.

Non-linear registration

Due to the large and delicate anatomical variability among population, using linear transformations that apply the same transform matrix to all voxel is not flexible enough to match the anatomy of a single person whose anatomy has changed locally. In this case, the transform matrix can be replaced by a deformation field D that associates a 3D displacement vector v(x, y, z) = [d x , d y , d z ] with voxel located at (x, y, z). In practice, the registration begins normally with a linear stage to roughly align the image and then is followed by a non-linear registration for refined matching.

For biomedical applications, it is crucial to apply constraints to the transformation in order to exhibit special properties including inverse consistency, symmetry, topology preservation, diffeomorphism, etc. [START_REF] Sotiras | Deformable Medical Image Registration: A Survey[END_REF]. The constraints are achieved by adding a regularisation term during the optimization process. Two reviews of nonlinear registration are available in [START_REF] Holden | A Review of Geometric Transformations for Nonrigid Body Registration[END_REF] and [START_REF] Sotiras | Deformable Medical Image Registration: A Survey[END_REF]. According to Sotiras et al. [START_REF] Sotiras | Deformable Medical Image Registration: A Survey[END_REF], the deformable transformation model can be classified into 4 groups.

Transformations derived from physical models [START_REF] Modersitzki | Numerical Methods for Image Registration[END_REF]: five sub-categories are contained in this group: elastic body models, viscous fluid flow models, diffusion models, curvature registration, and flows of diffeomorphisms.

The elastic body algorithms model the images as an elastic body and use the Navier-Cauchy Partial Differential Equation to describe the deformation [START_REF] Christensen | Consistent Image Registration[END_REF]. In viscous fluid flow models, the image under deformation is modeled as a viscous fluid and the transformation is governed by the Navier-Stokes equation [START_REF] Christensen | Deformable Templates Using Large Deformation Kinematics[END_REF]. The diffusion models execute the regularisation by a diffusion equation with Gaussian kernel [START_REF] Thirion | Image Matching as a Diffusion Process: An Analogy with Maxwell's Demons[END_REF]. Curvature-based registration models the deformation by an equilibrium equation [START_REF] Fischer | Curvature Based Image Registration[END_REF]. The flows of diffeomorphisms model the deformation as the integration of velocity of voxel over time [START_REF] Trouvé | Diffeomorphisms Groups and Pattern Matching in Image Analysis[END_REF].

Given two manifolds M and N, a differentiable map f : M → N is called diffeomorphism if f is bijective and f -1 is also differentiable. The properties of diffeomorphism such as smoothness, invertibility and topology-preservation fit the constrains for the registration task and thus some frameworks based on diffeomorphism have been proposed.

Under the Large Deformation Diffeomorphic Metric Mapping (LDDMM) framework [START_REF] Christensen | Deformable Templates Using Large Deformation Kinematics[END_REF], a diffeomorphism ϕ ∈ Di f f (M), M ∈ R 3 can be parameterized with a time-varying velocity field υ(t), therefore the diffeomorphisms are the flow of vector fields γ(t) initiated at the identity given by the following ordinary differential equation (ODE):

∂γ(x, t) ∂t = υ(t, γ(x, t)), γ(x, 0) = 0 (2.5)
The transformations is generated by the integration of this ODE: Transformations derived from interpolation theory: Interpolation theory or approximation theory are used to build models of this class. The approximation theory assumes that displacement estimates are inaccurate. Thus, instead of taking exact values of known displacements, the transformation approximates them smoothly. The methods in this category include radial basis functions [START_REF] Yang | Topology Preservation Evaluation of Compact-Support Radial Basis Functions for Image Registration[END_REF], elastic body splines [START_REF] Davis | A Physics-Based Coordinate Transformation for 3-D Image Matching[END_REF], free form deformations [START_REF] Sederberg | Free-Form Deformation of Solid Geometric Models[END_REF] and signal representation [START_REF] Christensen | Consistent Image Registration[END_REF].

ϕ(x) = γ(x, 1) = x + 1 0 υ(t, γ(t))dt (2.6)

Knowledge-based geometric transformations:

In medical image analysis, the registration may be performed to a specific target image or may involve a specific anatomical organ. Thus it is possible to introduce knowledge about the deformations. There are two ways to introduce the knowledge of deformations. One is performed by learning the statistic based on training data and serves as prior to constrain the registration, such as applying principle component analysis (PCA) [START_REF] Cootes | Active Shape Models-Their Training and Application[END_REF] or Gaussian mixture model (GMM) [START_REF] Glocker | Dense Registration with Deformation Priors[END_REF]. The other method is based on biomechanical or biophysical knowledge and adds the biomechanical or biophysical model to regularize the transformation [START_REF] Konukoglu | Image Guided Personalization of Reaction-Diffusion Type Tumor Growth Models Using Modified Anisotropic Eikonal Equations[END_REF].

Task-specific constraints: Some methods were proposed with a specific goal, such as topology preservation [START_REF] Noblet | 3-D Deformable Image Registration: A Topology Preservation Scheme Based on Hierarchical Deformation Models and Interval Analysis Optimization[END_REF], volume preservation [START_REF] Greene | Constrained Non-Rigid Registration for Use in Image-Guided Adaptive Radiotherapy[END_REF] and rigidity constraints [START_REF] Staring | A Rigidity Penalty Term for Nonrigid Registration: A Rigidity Penalty Term for Nonrigid Registration[END_REF]. The constraints for specific tasks often are used in conjunction with physics-based models and interpolation-based models.

Deep-learning-based registration

The use of deep learning techniques in medical image analysis has increased significantly over the past few years. Deep-learning-based methods have changed the landscape of medical image processing research in a number of fields and attained state-of-the-art performance in a variety of areas, including medical image registration tasks. A review of deep-learning-based medical image registration has been reported in [START_REF] Haskins | Deep Learning in Medical Image Registration: A Survey[END_REF]. The authors proposed to divide the deep-learning-based methods into three categories: deep iterative registration, supervised transformation estimation and unsupervised transformation estimation.

The deep iterative registration consists of deep-similarity-based registration and reinforcement-learning-based registration. The deep-similarity-based registration inserts the similarity into a conventional registration framework with pre-processing, transformation estimation and optimization. The deep network aims to learn a similarity metric with the backward propagation [START_REF] Wright | LSTM Spatial Co-Transformer Networks for Registration of 3D Fetal US and MR Brain Images[END_REF]. Instead of learning similarity metrics, the reinforcement-learning-based approaches use a trained agent to replace the pre-defined optimization algorithm [START_REF] Krebs | Robust Non-Rigid Registration Through Agent-Based Action Learning[END_REF].

Recent works have confirmed the performance of deep iterative registration algorithms. However, the iteration requires a large resource of computation and leads to the difficulty of real-time implementation, especially in high-dimensional deformable registration. Therefore the deep transformation estimation approaches are proposed.

The fully supervised registration uses a neural network to estimate the transformation and optimize with a loss function determined by ground truth data. The ground truth data can be real transformation generated by classical methods [START_REF] Yang | Fast Predictive Image Registration[END_REF], synthetic transformation [START_REF] Sokooti | Nonrigid Image Registration Using Multi-Scale 3D Convolutional Neural Networks[END_REF]. For cases lacking the ground truth transformation, Fan et al. [START_REF] Fan | BIRNet: Brain Image Registration Using Dual-Supervised Fully Convolutional Networks[END_REF] proposed a dual supervision strategy where the loss function is determined by the similarity and the ground truth. Hu et al. applied the GAN structure and optimized the network by simultaneously maximizing the similarity term and minimizing an adversarial loss term [START_REF] Hu | Adversarial Deformation Regularization for Training Image Registration Neural Networks[END_REF].

The supervised methods require a high number of manually annotated ground truths. Also, the reliability of ground truth significantly influences the performance of the trained model. This motivated the appearance of unsupervised models. Some works proposed to utilize the similarity metric to train the models. In [START_REF] Dalca | Unsupervised Learning for Fast Probabilistic Diffeomorphic Registration[END_REF], Dalca et al. combined transformer with diffeomorphism to obtain the velocity and the model is optimized with mean squared error. In [START_REF] Fan | Adversarial Similarity Network for Evaluating Image Alignment in Deep Learning Based Registration[END_REF], Fan et al. proposed a model based on GAN and used the discriminator to qualify the registration performance, instead of the utilization of manually crafted metrics. On the other hand, some works applied a feature-based strategy. Such as in [START_REF] Liu | Tensor-Based Descriptor for Image Registration via Unsupervised Network[END_REF], Liu et al. proposed a tensor-based MIND method using a principle component analysis-based network.

It has been reported the outperformance of the deep-learning-based registration methods to the conventional registration models. However, the limitations of this type of approach may limit the application in the real world. Firstly, the large dataset is obligatory for the training of deep neural networks which is challenging for medical applications. Secondly, the training of deep neural networks has a high demand on computational resources.

Construction of atlas

The construction of atlas is one special application of inter-subject registration where all subjects are registered to a common reference, which can be an initial reference, the other subjects or a common template space. Thus, the atlas is an average model of all subjects used in the construction and therefore a powerful tool to understand shape variability and population-wise difference by computing statistics of [START_REF] Whitcher | Statistical Group Comparison of Diffusion Tensors via Multivariate Hypothesis Testing[END_REF].

To produce an atlas, two main categories of algorithms exist: template-based approach and template-free approach.

The template-based algorithms require an image as the initial template and the other images are registered to this image. The workflow of the template-based method (see Figure 2.11), as presented in [START_REF] Guimond | Average Brain Models: A Convergence Study[END_REF], begins with affine and non-linear registration to the initial. The transformation of image I i to template I re f denotes as T i . Then all transformations T i , i = 1, 2, . . . , N are averaged and the mean transformation T is produced. The T is inversed and denotes as T-1 . The inversed mean transformation T-1 is applied to each transformed image I i,R = I i • T i • T-1 and the new reference I R is produced by averaging the transformed images. By iterating this procedure, the bias towards the initial template will be eliminated, as reported in [START_REF] Guimond | Average Brain Models: A Convergence Study[END_REF].

The template-free algorithms do not rely on an initial reference and therefore avoid the bias towards the reference. The approach in this category can be divided into 2 sub-groups: pairwise method and groupwise method.

The pairwise method has been proposed by [START_REF] Loeckx | Nonrigid Image Registration Using Free-Form Deformations with a Local Rigidity Constraint[END_REF]. The principle idea of this method is to exhaustively calculate the transformation T i,j between all the pairs of images I i , i = 1, 2, . . . , N and I j , j = 1, 2, . . . , N. For each image I i , i = 1, 2, . . . , N, we calculate the average transformation Ti by averaging the T i,j , j = 1, 2, . . . , N.

Then the average transformation Ti

-1 is applied on image I i to produce images I i,R

and the final atlas will be produced by averaging I i,R . The workflow of pairwise method is illustrated in Figure 2.12. The pairwise method is unbiased but requires a quadratic amount of registrations with respect to the number of subjects.

For groupwise method, all images are simultaneously registered to a common space instead of on a specific image space (see Figure 2.13). Therefore this unbiased approach requires only one single groupwise registration and has less demand of computation resource. One state-of-the-art implementation of this method is Advanced Normalization Tools (ANTS) software [START_REF] Avants | Geodesic Estimation for Large Deformation Anatomical Shape Averaging and Interpolation[END_REF], [START_REF] Brian | Advanced Normalization Tools (ANTS)[END_REF]. FIGURE 2.13: Overview of groupwise atlasing method.

Voxel-based shape analysis

In biomedical applications, the shapes are originally extracted by volumetric images. In voxel-based approaches, the correspondence is established by aligning the volumetric representation, such as MRI, to an atlas. The morphometry can be spatially quantified with the help of the deformation fields (deformation-based morphometry, DBM) or their Jacobian/log-Jacobian determinant (tensor-based morphometry, TBM) generated from the registration procedure, which represents the local volume change between the source image and atlas [START_REF] Ashburner | Identifying Global Anatomical Differences: Deformation-Based Morphometry[END_REF].

The Jacobian determinant of a deformation field encodes the local volume change between the source and target image, i.e., a value of 0.9 refers to 10% of volume reduction, and 1.1% denotes 10% of volume augment. Since no folding appears in medical image registration, the values in the Jacobian matrix are asymmetric that range from [0, +∞]. The logarithmic transformation of the Jacobian matrix makes the distribution symmetric and is often applied in practice [START_REF] Leow | Statistical Properties of Jacobian Maps and the Realization of Unbiased Large-Deformation Nonlinear Image Registration[END_REF]. In the following sections, we will refer to the log-Jacobian determinant as log-Jacobian.

The goal of DBM and TBM is then to detect significant differences between groups, which can be achieved with voxel-by-voxel statistical models. The pathological area is then identified by the p-value map. In neuroscience, Studholme et al. analyzed the relation between semantic dementia and the tissue contraction in the left temporal pole, the hippocampus, occipitotemporal gyrus and parahippocampal gyrus with B-spline-based method [START_REF] Studholme | Deformation Tensor Morphometry of Semantic Dementia with Quantitative Validation[END_REF]. In [START_REF] Song | Deformation-Based Morphometry Identifies Deep Brain Structures Protected by Ocrelizumab[END_REF], authors analyzed the deep brain structure based on brain MRI and TBM. In [START_REF] Cao | The Detection of Local Shape Changes via the Geometry of Hotelling's $T 2$ Fields[END_REF], Cao et al. applied the framework of DBM to analyze brain shape differences.

Another voxel-based shape analysis approach is voxel-based morphometry (VBM).

This technique is mainly used in neuroscience. In VBM [START_REF] Ashburner | Voxel-Based MorphometryThe Methods[END_REF], all images are spatially normalized into an atlas with a registration procedure. Then the segmentation of tissues (gray matter, white matter and cerebrospinal fluid) is smoothed and the statistical analysis is performed on the masks to localize the volume change. However, VBM is more sensitive to artifacts [START_REF] Ashburner | Computational Anatomy with the SPM Software[END_REF].

Surface-based shape analysis

The shapes can also be represented as surfaces with a set of points distributed across the surface extracted from volume data, i.e.

A i = {α k i ∥ α k i ∈ R 3 , k = 1, 2, .
. . , n}. α i are n landmarks of the shape A i . The shape A i is therefore described as a vector α, where the coordinates of each landmarks

α k i = (x k i , y k i , z k i ) are concatenated: α = (x 1 i , y 1 i , z 1 i , . . . , x n i , y n i , z n i ).
In order to reconstruct the surface, connectivity information between the points is stored. A point set with connectivity information constitutes a surface mesh. For surface-based SSM, various open-source tools based on pairwise or groupwise methods have been developed. In this section, we introduce three state-of-the-art tools that are adapted to all anatomical structures, called SPHARM-PDM, Deformetrica and ShapeWorks. SPHARM-PDM: SPHARM-PDM is an SSM tool with a parameterization-based correspondence construction method adapted to the anatomies with spherical topologies. Each shape is registered to a unit sphere by using a spherical harmonic (SPHARM) basis function to obtain the spherical parameterization. This procedure is optimized with area preservation and distortion minimization. The SPHARM basis functions Y m l are defined with degree l and order m:

Y m l (θ, ϕ) = 2l+1 4π
(l-m)!

(l+m)! P m l (cosθ)e imϕ , where θ ∈ [0, π], ϕ ∈ [0, 2π], and P m l the associated Legendre polynomials. The SPHARM expression of a surface is obtained by decomposing the three coordinate functions and the surface is transformed from

v(θ, ϕ) = (x(θ, ϕ), y(θ, ϕ), z(θ, ϕ)) T to the form of v(θ, ϕ) = ∑ ∞ l=0 ∑ l m=-l c l,m Y m l (θ, ϕ)
, where c l,m is 3D coefficient vector obtained by using a least-square optimization. Then the correspondence is established by using a first-order ellipsoid from the coefficient obtained through the spherical parameterization.

ShapeWorks:

ShapeWorks is an SSM tool based on a groupwise particle-based shape modeling algorithm. In particle-based shape modeling, each shape is treated as a collection of interacting dynamic particles with mutually repelling forces to optimally cover. The correspondence is constructed by optimizing the shape space and the particle position. For 3 , the correspondence between shapes is estimated by a rigid or similarity transformation to transform the particles to the common coordinate system. Two types of random variables are involved in the optimization of correspondence: a shape space variable S ∈ R 3M and a particle position variable A i ∈ R 3 . The objective of optimization is therefore to minimize the function

N shapes A 1 , A 2 , . . . , A N , each shape has M parti- cles A i = α 1 i , α 2 i , . . . , A M i ∈ R
Q = H(S) -∑ N i=1 H(A i )
where H is an entropy estimation under the assumption of Gaussian distribution in the shape space and Euclidean inter-particle repulsion in the configuration space.

Deformetrica: Deformetrica is an SSM tool based on LDDMM framework and requires an initial template to estimate the final atlas and the correspondence between shapes and the atlas. The construction of atlas of Deformetrica is performed with a Bayesian framework [START_REF] Gori | A Bayesian Framework for Joint Morphometry of Surface and Curve Meshes in Multi-Object Complexes[END_REF]. The shape instance A n is modeled as deformation and a residual, which is defined as A n = T n + ϵ n , where T n is the deformation between n-th shape and the atlas and ϵ n is the residual. The deformation of each shape is characterized by a set of parameters β n . We assume that the parameters follow a Gaussian distribution with a covariance matrix Γ β . Therefore the objective is maximizing the joint posterior distribution of shapes. For the i-th iteration, the objective function is

M i , Γ β i = arg max M i-1 ,Γ β i-1 p(M i-1 , Γ β i-1 ∥ {A n } N n=1
). The optimization is constrained by the fact that the residual should be small and the deformed template should match the shapes.

Benchmark of surface-based SSM tools

SSM tool assessment and understanding of the consequences of SSM tools is critical for the biomedical scenario, since different tools may give different results on synthetic nonclinical data [START_REF] Gao | Shape Analysis, a Field in Need of Careful Validation[END_REF]. In [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF], authors evaluated the SSM tools previously introduced in the context of clinical application. More specifically, the SSM is evaluated with the classical quantitative and qualitative evaluation metrics, as well as the validation with the landmarks/measurements inference and the lesion screening.

Quantitative evaluation

The quantitative metrics to evaluate the model are compactness, generalization and specificity.

The compactness measures the ability to cover as much as possible of variation with as few as possible of parameters. The compactness of m mode is defined as

C(m) = ∑ m j=1 λ j ∑ n-1 i=1 λ i
, where n is the number of shapes used for modeling. The faster the compactness converges to 1, the more compact the model is.

The generalization quantifies the ability to generate the valid instance under the same distribution as training set but not in it. The generalization of m mode is defined as

G(m) = 1 n ∑ n i=1 |z ′ i (m) -z i |,
where s ′ i is the reconstruction of the unseen sample s i , i.e. the generalization is the reconstruction error of an unseen valid instance. Lower this error is, more general the model is.

The specificity is to evaluate the ability of generate new valid instance by randomly generating the parameters. The specificity of m mode is defined as

S(m) = 1 k ∑ k i=1 |z" i (m) -z ′ i |, where z" i (m)
is the new generated instance with m principal components and z ′ i is the nearest shape in the training set. k is the number of new generated samples which need to be large. Smaller the distance between generated sample and its nearest sample is, more specific the model is.

Results: For left atrial appendage (LAA) and Scapula, ShapeWorks outperformed Deformetrica and SPHARM-PDM in compactness. For Humerus, SPHARM-PDM produced a compact model with one split of the data and ShapeWorks obtained a compact model with the split of the data. For femur, SPHARM-PDM performed better in compactness. In the aspect of generalization, ShapeWorks performed better in the LAA, humerus and femur dataset. For scapula, Deformetrica and Shape-Works were comparable. For specificity, ShapeWorks outperformed in the LAA, the humerus and the femur dataset. Deformetrica and ShapeWorks achieved comparable results in the scapula dataset.

Qualitative evaluation

The authors assessed qualitatively the model with the modes of variation of PCA, and cluster analysis. The modes of variation can reflect clinically relevant patterns. Therefore the model can be assessed based on their ability to discover the relevant modes. The clustering is an approach to find natural groups in populations based on the similarity. In clinical context, it can be used for diagnosis such as distinguishing the different levels of illness.

Results: For the LAA, the scapula and the femur shapes, Deformetrica and ShapeWorks was able to capture the clinically relevant variation modes, and correctly identify the clinical differences between control and pathology groups. However, SPHARM-PDM could neither produce the clinically relevant modes of variation, nor detect the group-wise differences. For humerus dataset, all three SSM tools discovered clinically relevant modes of variation, and were able to encode the differences between groups. But the differences detected by SPHARM-PDM has a visible bias.

Landmarks/measurements inference SSM is able to execute the personalized anatomical landmarks and measurements inference by mapping the mean shape to patient space with the correspondence. For ShapeWorks and Deformetrica, the patient-tomean correspondence is constructed by using thin plate splines (TPS) to generate the mean landmarks. Then the patient-specific predicted landmarks are obtained by using a TPS warp built with the mean-to-patient correspondence. For SPHARM-PDM, the tool does not provide the patient-to-mean correspondence. Therefore, the landmarks on the mean shape are manually annotated. The patient-specific prediction from SPHARM-PDM is performed with Procrustes fit to align the annotated landmarks to patient space. Paired sample t-test is used to compare the distance between the ground truth and the prediction. The landmarks/measurements inference validation was performed with the LAA, the scapula and the humerus dataset Results: For all three dataset, Deformetrica and ShapeWorks achieved a comparable performance and better than that of the SPHARM-PDM.

Lesion screening Lesion screening aims to locate the abnormal changes in subjectspecific anatomy and classifies the subject's anatomy to two groups control and pathology, based on the extent of the lesion. The lesion screening can be performed with the deviation between the original pathological shape and its reconstruction based on the model of the control group, i.e., the pathological shape is projected to the PCA subspace of the control population and reconstructed from this projection.

More details of reconstruction-based lesion detection will be introduced in section 2.3.3. The lesion screening is performed on the femur and humerus datasets.

Results: For the lesion identification, ShapeWorks and Deformetrica correctly identified the differences of the femur dataset, while SPHARM-PDM failed it. For humerus, all models succeeded in the identification but SPHARM-PDM detected false positives. For lesion classification, ShapeWorks and Deformetrica were comparable on femur dataset and better than SPHARM-PDM. On humerus dataset, the performance of ShapeWorks is better than Deformetrica and SPHARM-PDM.

Anomaly detection

An anomaly (also called outlier or novelty) is an observation that deviates considerably from some concept of normality [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]. The normality is defined as the distribution P + on X which is the ground truth law of normal behavior in the given task. An observation that deviates significantly from P + is a data point x ∈ X that lies in a low probability region under P + . Denoting the probability distribution function of normality as p + (x), the set of anomalies is defined as A = {x ∈ X | p + (x) ≤ τ} , τ ≥ 0, where τ is a detection threshold that defines the low probability region under P + .

The fundamental assumption of anomaly detection is the concentration assumption [START_REF] Schölkopf | Learning with Kernels: Support Vector Machines, Regularization, Optimization, and Beyond[END_REF], [29]. The data space X is unbounded in most application situation, while the distribution of normal data is bounded by a threshold τ. However, the law of normality P + is known in certain applications. In most situations, the P + remains unknown due to the complex underlying process. Therefore the P + needs to be estimated from the given data. And under the distribution estimated, we could predict whether a new instance x ∈ A. The anomaly detection objective is now to estimate the low-density regions (or equivalently high-density regions), which can be achieved by the statistical hypothesis test [START_REF] Tsybakov | On Nonparametric Estimation of Density Level Sets[END_REF], [148] in data space X . Given a density level α, a corresponding threshold τ α will be determined to decide whether the instance x is an anomaly or not.

As presented in Figure 2.15, according to the methods of normality estimation, three classes of algorithms to perform anomaly detection have been introduced [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]: probabilistic models, one-class classification and reconstruction models. Also, with the development of deep learning, many deep anomaly detection methods appeared in recent years. Due to the limitations we have presented in section 2.2.1.4, in this work, we focus more on classical shallow methods.

Probabilistic models

The target of probability-based methods is to estimate the distribution of normal data, which meets the need of anomaly detection problems. Methods in this category include statistical hypothesis tests, density estimation and deep statistical models. A statistical hypothesis test is a method of statistical inference that associates the value of probability with the decision-making in order to quantify confidence in a hypothesis. In a hypothesis test, two types of hypotheses exist:

-H 0 is the null hypothesis and is for being accepted, which means the difference is caused by sampling error or inter-subjects variability -H 1 is the alternative hypothesis and is for being rejected which indicates the significant difference.

To decide whether H 0 is accepted or rejected, the significance level α which is a probability threshold for rejecting the null hypothesis. As shown in Figure 2.14, if p( x) ≤ τ α , i.e. it is located in rejection region (red zone), it will be regarded as an anomaly. The selection of density level α requires the specialization to application. The common choice of α is 0.05 or 0.01. To make the decision, the p-value of the test hypothesis is calculated. If p < α, then H 0 is accepted, else reject H 0 and accept Two types of errors exist in statistical hypothesis tests: type uppercasei error and type II error. As presented in table 2.1, the type I error is the false positive findings and the type II error is the false negative findings. To eliminate the type I error, multi-test correction is required, such as with the false discovery rate (FDR) method. For the type II error, one way to correct it is to augment the number of samples.

H 1 . truth decision H 0 H 1 H 0 1 -α (true positive) β (type II error) H 1 α (type I error) 1 -β (true negative)
The common used tests include Hotelling's T 2 test [START_REF] Fisher | Statistical Methods for Research Workers[END_REF], analysis of variance (ANOVA) [START_REF] Keselman | Statistical Practices of Educational Researchers: An Analysis of Their ANOVA, MANOVA, and ANCOVA Analyses[END_REF], analysis of covariance (ANCOVA) [START_REF] Keselman | Statistical Practices of Educational Researchers: An Analysis of Their ANOVA, MANOVA, and ANCOVA Analyses[END_REF] and their multivariate version MANOVA [START_REF] Keselman | Statistical Practices of Educational Researchers: An Analysis of Their ANOVA, MANOVA, and ANCOVA Analyses[END_REF], [START_REF] Carl | Multivariate Analysis of Variance and Covariance[END_REF], MANCOVA [START_REF] Carl | Multivariate Analysis of Variance and Covariance[END_REF], etc. Hotelling T 2 distribution is a multivariate probability distribution, a generalization of the Student's T-distribution and tightly related to the F-distribution. For 2-sample Hotelling T 2 statistic, we assume that 2 groups of data follow multivariate normal distributions with the same mean and covariance. Suppose the deformation vector W ij (x) matches the structure at the position x in atlas of subject i in group j, the T 2 statistic is defined as [START_REF] Paul | Detection and Mapping of Abnormal Brain Structure with a Probabilistic Atlas of Cortical Surfaces[END_REF], [START_REF] Cao | The Detection of Local Shape Changes via the Geometry of Hotelling's $T 2$ Fields[END_REF]:

T 2 (x) = N 1 N 2 (N 1 + N 2 )(N 1 + N 2 -2) W µ 2 (x) -W µ 1 (x) T Ψ -1 (x) W µ 2 (x) -W µ 1 (x) (2.7
) N 1 and N 2 are the number of subjects of group 1 and group 2. µ 1 and µ 2 are the respective mean of two groups. Ψ is the covariance matrix between two groups. At each point, let v = (N 1 + N 2 -2) and d = dimension so d = 3 in our case, then the T 2 statistic can be transformed to a F statistic to calculate the p-values: (v-d+1) . Two examples of anomaly detection with Hotelling's T 2 test are illustrated in Figure 2.16, subfigure (a).

F(x) = ((v -d + 1)/d)T 2 (x) ∼ F d,
To model more complex distributions, non-parametric density estimators have been introduced, such as kernel density estimators (KDE) [START_REF] Parzen | On Estimation of a Probability Density Function and Mode[END_REF], [START_REF] Härdle | Chapter 38 Applied Nonparametric Methods[END_REF], histogram estimators, and Gaussian mixture models (GMMs) [START_REF] Roberts | A Probabilistic Resource Allocating Network for Novelty Detection[END_REF], [START_REF] Bishop | Novelty Detection and Neural Network Validation[END_REF]. GMM is an extension of the Gaussian model, with a finite number of K mixtures to describe the distribution of data. It can also be viewed as a clustering method that assumes K prototypical modes. KDE is the most widely used non-parametric density estimator. For (x 1 , x 2 , ..., x n ) which are independent and identically distributed samples, f is their density function. f can be estimated as follows:

f = 1 nh n ∑ i=1 K( x -x i h ),
where K is a non-negative kernel function whose integral is 1. Various kernel functions exist such as Gaussian kernel, Epanechnikov kernel, etc. h is a positive parameter called bandwidth. The early deep statistical models are based on energy-based models (EBMs) [START_REF] Fahlman | Massively Parallel Architectures for AI: NETL, Thistle, and Boltzmann Machines[END_REF], [START_REF] Hopfield | Neural Networks and Physical Systems with Emergent Collective Computational Abilities[END_REF], [START_REF] Lecun | A Tutorial on Energy-Based Learning[END_REF], such as Deep Belief Networks [START_REF] Hinton | A Fast Learning Algorithm for Deep Belief Nets[END_REF] and Deep Boltzmann Machines [START_REF] Salakhutdinov | Efficient Learning of Deep Boltzmann Machines[END_REF]. EBMs characterize the density by an energy function E θ (x) with

p θ (x) = 1 exp(-E θ (x))dx exp(-E θ (x))
Recently, Generative Adversarial Network (GAN) [START_REF] Goodfellow | Generative Adversarial Networks[END_REF] is proposed to deal with anomaly detection problems. GAN consists of 2 neural networks: a generator network ϕ : Z → X is trained in competition with an adversary that challenges it to generate samples whose distribution is similar to the training samples. A discriminator network ψ : X → (0, 1) is trained to discriminate between the generated samples and real data. Using the discriminator directly has been suggested as one approach to using GANs for anomaly detection which has a similar principle as one-class classification [START_REF] Sabokrou | Adversarially Learned One-Class Classifier for Novelty Detection[END_REF]. The other approach is to use the generator to perform a reconstruction-based anomaly detection [START_REF] Deecke | Image Anomaly Detection with Generative Adversarial Networks[END_REF]. Also in [START_REF] Schlegl | Unsupervised Anomaly Detection with Generative Adversarial Networks to Guide Marker Discovery[END_REF], authors defined the anomaly score as the combination of the reconstruction loss and the discrimination loss.

One-class classification

One-class classification is based on the same principle of probabilistic methods which is to study a statistical threshold between normal data and anomaly but oneclass classification learns the boundary directly instead of learning the distribution as an intermediate step. The most well-known methods are support vector data description (SVDD) [START_REF] David | Support Vector Domain Description[END_REF] and one-class support vector machine (OCSVM) [START_REF] Schölkopf | Support Vector Method for Novelty Detection[END_REF]. The objective of SVDD is to learn a hypersphere of minimum volume to separate normal data from outliers. Therefore, the optimization of SVDD can be defined as

min a,R,ξ R + C n ∑ i=1 ξ i s.t. ∥Φ(x i ) -a∥ 2 < R 2 + ξ i , ξ i ≥ 0 = 1, 2, . . . , n (2.8)
where R and a are the radius and center of the hypersphere respectively. ξ is the vector of slack variables and C is the regularization parameter to set a boundary to the fraction of allowed outliers. Φ(x) is a kernel function mapping the data from the original space to the feature space.

In comparison, OCSVM uses a hyperplane to separate data in feature space with a maximum margin from the origin:

min ω,ρ,ξ 1 2 ∥ω∥ 2 + 1 νn n ∑ i=1 ξ i -ρ s.t. ω • Φ(x i ) > ρ -ξ i , ξ i ≥ 0 = 1, 2, . . . , n (2.9)
ω and ρ define the character of hyperplane and ν plays the same role as C in SVDD to set the fraction of allowed outliers. In subfigure (c) of Figure 2. [START_REF] Bell | Talus[END_REF], examples with ball dataset and double moon dataset solved by OCSVM is presented.

Also, for more complex data, there exist deep one-class classification methods such as Deep SVDD [START_REF] Ruff | Deep One-Class Classification[END_REF] and Deep OCSVM [START_REF] Erfani | High-Dimensional and Large-Scale Anomaly Detection Using a Linear One-Class SVM with Deep Learning[END_REF], etc.

Reconstruction-based methods

Reconstruction-based methods learn a model that is optimized to well reconstruct normal data, thereby aiming to detect anomalies by failing to accurately reconstruct them under the model. The pipeline of reconstruction consists of 2 parts: decomposition and reconstruction. The decomposition is realized by an encoding function ϕ e which maps data from data space to a latent space. The reconstruction is the procedure that reconstructs data from the latent representation obtained in decomposition. It is realized by a decoding function ϕ d . Therefore, for data (x 1 , x 2 , . . . , x n ), the objective of the model is

min 1 n n ∑ i=1 ∥x i -(ϕ d • ϕ e )(x i )∥ 2 + R (2.10)
where R is the regularization term. Principle component analysis (PCA) [START_REF] Tipping | Probabilistic Principal Component Analysis[END_REF], Autoencoders (AEs) [START_REF] Pawlowski | Unsupervised Lesion Detection in Brain CT Using Bayesian Convolutional Autoencoders[END_REF], [START_REF] Chen | Unsupervised Detection of Lesions in Brain MRI Using Constrained Adversarial Auto-Encoders[END_REF], VAEs [START_REF] Diederik | An Introduction to Variational Autoencoders[END_REF], [START_REF] You | Unsupervised Lesion Detection via Image Restoration with a Normative Prior[END_REF] and GANs [3], [START_REF] Schlegl | F-AnoGAN: Fast Unsupervised Anomaly Detection with Generative Adversarial Networks[END_REF] are common used methods in this category. The reconstruction error can be simply defined as Euclidean distance, Mahalanobis distance [START_REF] Mclachlan | Mahalanobis Distance[END_REF] or using some Machine Learning classification algorithm such as clustering [5] or OCSVM [4].

PCA [START_REF] Tipping | Probabilistic Principal Component Analysis[END_REF] is one of the most-used classic dimension reduction algorithms. Thanks to its natural out-of-sample extension and inherent reconstruction, it can be adapted for reconstruction-based anomaly detection. The objective of PCA is to seek an orthogonal principal basis W in data space over which data can be linearly projected with the maximized empirical variance:

max n ∑ i=1 ∥Wx i ∥ 2 s.t. WW T = I (2.11)
For reconstruction, the objective is to find an orthogonal projection W T W to a highdimension linear space with a minimized reconstruction error:

min n ∑ i=1 ∥x i -W T Wx i ∥ 2 s.t. WW T = I (2.12)
Linear PCA is only able to represent linear features. The disability of encoding the non-linear features will result in the mis-decomposition and reconstruction. Therefore, the kernel PCA [START_REF] Schölkopf | Nonlinear Component Analysis as a Kernel Eigenvalue Problem[END_REF],is introduced. As a non-linear extension of PCA, it introduced a non-linear kernel function ϕ(x) which maps data from original space to feature space, such as RBF kernel, etc. Thus, the objective of decomposition and reconstruction becomes:

max n ∑ i=1 ∥Wϕ(x i )∥ 2 s.t. WW T = I (2.13) min n ∑ i=1 ∥ϕ(x i ) -W T Wϕ(x i )∥ 2 s.t. WW T = I (2.14)
Examples for unimodal and bimodal anomaly detection with the ball dataset and double moon dataset are presented in Figure 2.16, subfigure (d) and (f).

For this type of methods, we can also identify their probabilistic interpretation. For linear PCA, the data distribution follows from the linear transformation of a d-dimensional latent Gaussian distribution with noise. And for KPCA with RBF kernel it can be considered as equivalent to kernel density estimation.

Autoencoder is a non-linear reconstruction model based on a neural network, consisting of two parts: encoder and decoder. The encoder is aimed at learning a low-dimensional representation of original data and the decoder is for the reconstruction. Different ways to regularize autoencoders have been introduced such as sparse autoencoders [START_REF] Makhzani | K-Sparse Autoencoders[END_REF], [START_REF] Zeng | Facial Expression Recognition via Learning Deep Sparse Autoencoders[END_REF], and denoising autoencoders [START_REF] Vincent | Extracting and Composing Robust Features with Denoising Autoencoders[END_REF], [START_REF] Vincent | Stacked Denoising Autoencoders: Learning Useful Representations in a Deep Network with a Local Denoising Criterion[END_REF], [START_REF] Marchi | A Novel Approach for Automatic Acoustic Novelty Detection Using a Denoising Autoencoder with Bidirectional LSTM Neural Networks[END_REF], etc.

Variational autoencoders (VAEs) is a generative model whose architecture is based on autoencoder. VAEs constrain the latent space by leveraging the encoder and decoder of AEs to parameterize an approximation to the posterior distribution. The input vector x are parameterized on latent space Q as an latent sample z by encodeing the posterior q ϕ (z|x) ∼ N (µ x , σ x ). The decoder parameterizes the mean and the variance of an isotropic Gaussian distribution, therefore p θ (x|z) ∼ N (x; µ(z), σ 2 (z)I). The loss function is the combination of reconstruction error and the KL-divergence: L ϕ,θ = L 1 (x, x) + D KL (q(z), p(z)). This method combines the principle of the probabilistic approach and the reconstruction-based model. OCSVM with RBF kernel has a similar effect of KDE. PCA finds a linear orthogonal subspace that is able to deal with unimodal data yet is less efficient for bimodal double-moon data. Kernel PCA introduced a nonlinear kernel which enables an optimal reconstruction from non-linear components in input space.

Discussion

In this chapter, we first introduce the definition of CP and the spastic equinus deformity in the patient's ankle joint, and the anatomy of the human ankle, in order to better understand the background of this thesis. As described in section 2.1.3, CP causes abnormal muscle development and leads to bone deformation. Yet little information is available on this deformation. The investigation of bone morphological deformation will enable a better comprehension of pathology and rehabilitation of patients. To this end, a morphological analysis needs to be performed.

SSM is a powerful tool to analyze the morphological properties of given shapes. It statistically describes the mean shape and the variation between the given shapes.

One key element of SSM is to build the correspondence between shapes, which can be processed by the image registration technique. Many algorithms have been developed to align two shapes. Under the medical application context, for the registration between the images or meshes of the same subjects, linear registration is often adopted. For cases more complex such the inter-subject registration, the non-linear transformation will provide more flexibility.

The SSM can be performed with voxel-based approaches and surface-based approaches. In voxel-based approaches, the shapes are in a volumetric representation such as images. It aligns the images into the mean template image and the morphological analysis is performed voxel-wisely with the deformation fields or their log-Jacobian. As consequence, the voxel-based analysis can be performed theoretically without the manual segmentation of target anatomy. And this category is able to examine not only the surface but also the inner tissue of the anatomy to be studied. However, the required computation resources are augmented with the increase of the image resolution for the registration and the following analysis step.

In surface-based methods, the shapes are represented as the surface with a set of points. The surface is generated from the segmentation of the target anatomy. Then the correspondence is estimated based on the registration between surface meshed. In ShapeWorks, the transformation is estimated with a global linear step and a nonlinear step with a physical model. In Deformetrica, the transformation is obtained by using the LDDMM framework. Compared to the voxel-based method, in bone shape analysis, the surface-based method is able to focus on the shape contour and ignore the osseous substance inside the bone, which can reduce the computation resource requirement. However the precision of the manual segmentation influences the precision of the surface representation. The generation and the processing of the surface mesh also play an important role in the following analysis. Moreover, the effect of different shape representations in shape analysis is under-studied and remains to be investigated.

To investigate the abnormal deformation of bone shape that occurs in children with CP, the anomaly detection technique will be applied for the population-wise comparison. Anomaly detection is based on statistic theory which aims to find the set that has a low probability in the normal sample distribution. According to [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]

Introduction

As presented in Chapter 2, the equinus deformity caused by CP probably leads to muscle weakness, spasticity and passive amplitude limitation of dorsiflexion, and furthermore, the change of bone morphology, which is still unreported in literature and requires to be explored. Pediatric MRI studies remain challenging to conduct due to complex acquisition settings [START_REF] Mitchell | Residual Equinus After the Ponseti Method: An MRI-Based 3-Dimensional Analysis[END_REF], [START_REF] Makki | In Vivo Ankle Joint Kinematics from Dynamic Magnetic Resonance Imaging Using a Registration-Based Framework[END_REF]. However, visualizing and analyzing morphological changes in ankle bones could provide a better understanding of CP pathology and therefore a more adapted setting up of rehabilitation for patients.

For this purpose, understanding the normality is a pre-requisite procedure, thus requiring us to study the shape of the ankle bone of TD population. Without the definition of standard bone shape as a reference, abnormal deformation would not be defined and detected. Statistical shape modeling (SSM) is a popular tool that provides the mean representation of the population as well as the shape variation [START_REF] Luthi | Gaussian Process Morphable Models[END_REF].

One way to perform the SSM is to compute the mapping of 3D images of each subject into a common reference image space. A series of deformation fields describes then the spatial transformations required to match the different shapes to the same atlas [START_REF] Ashburner | Identifying Global Anatomical Differences: Deformation-Based Morphometry[END_REF]. Information about the individual shapes is encoded in the deformation fields. Such a voxel-based approach is performed directly from scanned images and avoids the extraction and reconstruction of surface mesh [START_REF] Gao | Shape Analysis, a Field in Need of Careful Validation[END_REF].

The shape variation of ankle bones can be related to demographic characteristics, such as age, weight, height, and length of the lower limb. In addition, muscle quality influences bone morphology. Thus the shape variation may also be explained by the muscle strength and the passive dorsiflexion. A correlation analysis between these clinical coefficients and the model parameters is able to create a connection with clinical measurement and interpret the shape from a clinical view. For the further morphological comparison between the two populations, such a correlation analysis will facilitate the localization of disorders, understand the biomechanical interaction between bone morphology and muscle development, and result in a better treatment.

In this chapter, we study the ankle bone shape of a healthy population using voxel-based statistical shape modeling methods. Firstly, the data for this thesis, as well as the acquisition, are introduced. Secondly, the shape model of typical developing children is constructed, and the variation is analyzed with PCA. Thirdly, the correlation between clinical measurement and shape model parameters is analyzed, and thus the clinical interpretation of the model is derived.

Materials and methods

Populations

In this thesis, the prime objective is to investigate the bone morphometry of the ankle joint (calcaneus, talus and tibia) relevant to equinus caused by cerebral palsy. To this end, eleven TD children and nine children with CP with age ranging from 6 to 14 years old is involved in this study which was approved by the regional ethics committee. The CP group includes seven males and two females, and the TD group consists of seven males and four females with no history of the pathology of the lower limbs. The demographic characteristics, including age, weight, height and BMI of the two groups, are demonstrated in Table 3.1. The T-test is performed on these characteristics and no significant inter-group difference is noticed (p > 0.05). All children were selected with no MRI contraindications and no history of lower limb musculoskeletal injury or surgery in the past six months.

Clinical coefficient measurement

In order to determine the factors linked with the bone shape, the muscular strength, the spasticity of CP subjects and the passive range of motion are measured according to the measurement in [START_REF] Bailly | Relation Entre La Morphologie Osseuse Des Membres Inférieurs, La Fonction de La Marche, Les Limitations d'activité et Les Restrictions de Participation Des Enfants Avec Paralysie Cérébrale Marchants[END_REF]. The muscles accessed relevant to ankle joints are tibialis anterior, triceps sural, peroneus long and short and tibialis posterior. Two experienced examiners assessed the muscle strength with a dynamometer (Powertrack II Commander MMT, JTECH Medical). The assessment is performed with "make test": the dynamometer is kept in a static position, and the subject has to push against this resistance [START_REF] Nyström Eek | Isometric Muscle Torque in Children 5 to 15 Years of Age: Normative Data[END_REF]. Three trials are conducted sequentially, separated by periods of rest. The variation coefficient of three values given by the dynamometer has to be less than 10% and the maximum value is kept [START_REF] Nyström Eek | Isometric Muscle Torque in Children 5 to 15 Years of Age: Normative Data[END_REF]. The strength value is normalized by dividing the lower limb's length [START_REF] Ayalon | Reliability of Isokinetic Strength Measurements of the Knee in Children with Cerebral Palsy[END_REF].

Spasticity:

The spasticity of tibialis anterior, triceps sural, peroneus long and short and tibialis posterior is assessed by an experienced examiner with modified Ashworth scale. The 1+ rating was replaced by 1.5 for statistical analysis [START_REF] Mutlu | Reliability of Ashworth and Modified Ashworth Scales in Children with Spastic Cerebral Palsy[END_REF].

Passive range of motion:

The passive joint examination of the lower limb was performed at the ankle by two experienced examiners. Dorsiflexion was measured with the knee flexed and then with the knee extended. The subject was installed in the supine position [START_REF] Viehweger | Multidimensional Outcome Assessment in Cerebral Palsy: Is It Feasible and Relevant?[END_REF] and the measurement was performed using a manual goniometer. The goniometer branches were aligned with the axis of the tibial segment (marked by the lateral malleolus and the head of the fibula) and with the axis of the 5th metatarsal.

The results of measurement of two groups are presented in Table 3.2 as well as the T-test results for each coefficient.

3D image acquisition and bone-of-interest segmentation

3D MRI data were acquired in a single visit after parents signed informed consent forms, using a 3T MR scanner (Achieva dStream, Philips Medical Systems, Best, Netherlands) with a resolution of 0.26 × 0.26 × 0.8 mm 3 and resampled to 0.5 × 0.5 × 0.5 mm 3 for the purpose of adaptation to the clinic (T1-weighted gradient-echo, flip angle 10, matrix 576×576, FOV 150 mm × 150 mm, TR/TE 7.81/2.75 ms, mean acquisition duration: 424.32s). Images of the ankle were taken on the CP group's paretic lower limb and the non-dominant lower limb for the TD group. The acquisition protocol is detailed in [START_REF] Makki | In Vivo Ankle Joint Kinematics from Dynamic Magnetic Resonance Imaging Using a Registration-Based Framework[END_REF] and [START_REF] Garetier | Dynamic MRI for Articulating Joint Evaluation on 1.5 T and 3.0 T Scanners: Setup, Protocols, and Real-Time Sequences[END_REF].

In order to extract the shape of the three bones of interest, we make use of a semi-automatic segmentation approach. Age variability induces developmental variability of bones and cartilages. To overcome such shape variability, the considered regions of interest include bones and growth cartilages. First, three subjects of

3D voxel-based shape modeling

Atlas construction and deformation fields computation

To achieve the 3D shape modeling with a voxel-based approach, we make use of a registration-based framework to compute a mean image model of the TD population, called atlas. To deal with the high variability in shape and appearance of the bones, we propose to estimate the TD atlas using a group-wise diffeomorphic multivariate algorithm [START_REF] Avants | Geodesic Estimation for Large Deformation Anatomical Shape Averaging and Interpolation[END_REF] by not only considering image intensity but also bone segmentation maps and Chamfer distance maps. Such a multivariate approach ensures realistic atlas estimation with sharp details of cortical bones and cartilage. Then, for each subject of the dataset, deformation fields are estimated by non-linear multivariate registration (using segmentation and Chamfer distance maps) onto the mean TD atlas previously computed. The atlas estimation and the patient-to-atlas registration stages are performed with ANTS1 .

Statistical deformation variation analysis with PCA

PCA is one of the most-used dimension reduction algorithm. The objective of PCA is to seek an orthogonal principal basis in data space over which the data can be linearly and decorrelatedly projected with the maximized empirical variance, thus eliminating information redundancy.

In our hypothesis, we assume that a deformation field D follows a multivariate normal distribution D ∼ N (D, Σ), where D is the mean and the Σ is the covariance matrix. During the decomposition procedure, the shape is decomposed into n -1 principal components, also called eigenmode. Therefore, the new representation of a deformation field D is described as:

D = D + n-1 ∑ i=1 b i ϕ i (3.1)
where ϕ i are eigenvectors of covariance matrix Σ, and b i is a shape parameter, whose value is often defined as -3

√ λ i ≤ b i ≤ 3 √ λ i .
λ is the eigenvalues of covariance matrix thus it represents the variance. Under the hypothesis of D ∼ N (D, Σ), we have P(|b i | ≤ 3 √ λ i ) = 99.7%. Each principal component represents one linearly uncorrelated shape variation mode, and the principal components are arranged in the order of the importance, which means the first mode represents the largest variation of bone shape. With PCA, it is possible to visualize the principle variation exists in the dataset from a statistical view.

Clinical coefficient correlation analysis

To better understand the constructed shape model, as well as each eigenmode derived from PCA, we analyse the correlation between the shape parameters λ and the clinical coefficients, including demographic characteristics (age, height, weight, leg length), muscular strength (tibialis anterior, triceps sural, peroneus and tibialis posterior) and the dorsiflexion range (knee extended and flexed). The correlation is analysed with Pearson correlation coefficient, defined as ρ = cov(X TD ,X CP ) σ TD σ CP .

Experiments and results

Evaluation of atlas construction and deformation field generation

The estimated atlas, presented in Figure 3.3, provides visually realistic anatomy of the ankle in the region of the bones of interest (calcaneus, talus and tibia), including crisp details for the contour of bones and cartilages. To quantitatively evaluate the performance of atlas building and registration between each subject and the estimated atlas, the Dice coefficient which measures the overlap between regions of interest is used:

Dice(A, B) = 2(|A|∩|B|)
|A|+|B| . Experiments report that the 3D Dice coefficient for every bones of each subject is greater than 0.97 (Table 3.3), demonstrating high accuracy of shape matching. Examples of obtained deformation field of three bones are illustrated in

Quantitative evaluation of model

To evaluate the robustness of the shape model, some quantitative metrics are applied, including compactness, generalization and specificity. The introduction of these metrics can be found in Section 2.2.3.1.

The compactness of calcaneus, talus and tibia is reported in Figure 3.5-(A). The first mode captures most shape variations. For calcaneus and talus, it contains approximately 30% of shape variance (29.48% for calcaneus and 33.66% for talus). For tibia, the first mode includes 53.42% variations. The first four modes captured over 70% of variation for all three bones (76.08% for calcaneus, 72.95% for talus and 84.51% for tibia). To explain 90% of shape variation, calcaneus and talus require 7 principal components and tibia needs 5 eigenmodes. The evaluation of compactness indicates that the obtained model has the ability to represent the shapes in a compact way, especially for tibia, which is the most simple shape among the three bones. For calcaneus and talus, based on their anatomical complexity, they require more complex model than tibia.

To evaluate the generalization of our model, the leave-one-out evaluation is performed. The results of generalization of three bones are presented in Figure 3.5-(B). For all three bones, the curve of G(m) starts from a higher value (2.37mm for calcaneus, 1.56mm for talus and 1.60mm for tibia), and then the error gradually decreases as the number of modes increases. When n reaches to 9, the generalization is 1.85mm for calcaneus, 1.32 for talus and 1.17 for tibia. However, it is possible to notice that the curves do not converge, indicating that the model requires to re-train on a large dataset to improve the effectiveness.

The specificity is calculated, in our experiment, with 50 generated samples. The specificity of our models are in Figure 3.5-(C). The specificity of calcaneus, the most complex in anatomy, is found average, ranges from 1.46mm to 1.82mm. For talus and tibia, the specificity values are similar to each other and better than calcaneus, ranges from 1.00mm to 1.23mm. However, the curves are still fluctuating at the end of the experiment (number of modes equals to 10), which indicates the curves do not reach to convergence.

Shape variation analysis

Calcaneus

In Figure 3.6-3.11, we present the anterior (Figure 3.6), lateral (Figure 3.7), superior (Figure 3.8), inferior (Figure 3.9), medial (Figure 3.10) and posterior (Figure 3.11) view of first four eigenmode of calcaneus. The first eigenmode represents the bone size globally. At -3 √ λ 1 , the reconstructed calcaneus has a larger volume with an edge more rounded and less obvious. At 3 √ λ 1 , the bone appears flat, with a smaller volume, more sharp angles at the inflection points and clearer edges. In the visualization of deformation magnitude, the highlighted area with a magnitude over 2mm is extensively observed on the superior, lateral and inferior surfaces. On the lateral margin of the inferior surface, the magnitude is over 3mm. It confirms the representation of volume variation of the first mode. On the anterior view, the sustentaculum tali is highlighted with a magnitude over 4mm, illustrating that the first mode also includes the prominence of sustentaculum tali. The second mode mainly shows the size difference between the anterior surface (facet for cuboid) and the posterior surface (calcaneal tuberosity). At -3 √ λ 2 , the reconstruction tends to have a large cuboid facet and the tuberosity is small and flat, while at 3 √ λ 2 , the cuboid facet is thinner, with inward-recessed superior and inferior margin, and the tuberosity is thicker and more rounded. It can be visualized in the deformation magnitude of the inferior surface. The anterior and posterior areas of this surface are highlighted with a magnitude over 4mm, which signifies the change in thickness. As in the first mode, the second mode also represents the prominence of sustentaculum tali. The third mode includes the variation of the back tilt of the posterior facet. The -3 √ λ 2 demonstrates a large tilt while the 3 √ λ 2 shows a erect appearance. The deformation magnitude of the posterior margin of the posterior facet of the subtalar joint is over 4mm. Also, the third mode represents the change in cuboid facet thickness change. The fourth mode contains the variation of the depression of the tarsal sinus, which also influences the shape of the cuboid facet, and the curvature of the superior facet of the tuberosity. These deformations are marked with a magnitude over 2mm. √ λ i

Talus

The visualization of six views of the first four modes of the talus is in Figure 3.12 -3.17 (lateral view: Figure 3.12; medial view: The second mode consists of the variation of the curvature of the medial process. At -3 √ λ 2 , the medial process tend to depress and at 3 √ λ 2 , the medial process is processing. And the prominence of the lateral process and the location of the posterior process are included in the second mode. In mode 3, the variation is mainly located at the medial margin of the talar neck (flatter at 3 √ λ 3 and more curved at -3 √ λ 3 ), and the prominence of the posterior process (processing at -3 √ λ 3 and flat at 3 √ λ 3 ). The fourth mode represents mainly the width of the talar neck and head, and it influences the shape of the talar dome, while the deformation is less visible (magnitude around 2mm).

Tibia

The five view of the first four modes of the tibia is presented in Figure 3.18-3.22 (anterior view: -3 √ λ 1 , the reconstructed bone has a large volume, rounded processes on the lateral surface, a short and rounded medial malleolus and a triangular distal surface. At In each row from left to right: reconstructed shape at -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , reconstructed shape at 3

√ λ i 3 √ λ 1
, the reconstructed bone is small, with pointed processes on the lateral surface, large, squared and pointed medial malleolus and respectively squared distal surface. It can be observed that the first and the second modes represent the length of the tibia. Hence in this thesis, we only focus on the distal end instead of the whole tibia. Therefore the size, as well as the proximal end deformation, is not taken into consideration, and the manual segmentation causes the variation in length. The second mode contains similar features as the first mode and mainly focuses on the shape of the process and malleolus. The third mode primarily represents the shape of the medial malleolus, as in the first and second modes. And the fourth mode contained less variation and the visible shape change is similar to the first three modes. One highlighted deformation is located at the posterior margin of the distal surface, which changes the curvature of the respective bone edge.

Correlation analysis between shape parameters and clinical coefficients

In the previous section, we visualize the shape variations captured by each eigenmode of the shape model. But how to explain these variations remains to be discovered. The first eigenmode of each bone is assumed to relate to age, yet not proved. To access the factors influencing the bone shape, we calculate the Pearson correlation coefficient between the eigenvalues of covariance matrix λ and the clinical coefficients. The clinical coefficients include age, height, weight, leg length, tibialis anterior strength, peroneus and tricep sural, and the passive movement range with the knee extended and flexion.

In the first step, the correlation inter-clinical-coefficient is analyzed to avoid the influence of this correlation in the shape-clinical correlation analysis. The correlation After figuring out the correlation between the clinical coefficients, in the next step, the correlation analysis is performed between the clinical coefficients and the shape parameters λ. The correlation analysis results (see Figure 3.24, 3.25, 3.26) confirm our assumption in shape variation visualization that the first mode of all three bones is strongly related to age (calcaneus: -0.65; talus: -0.68; tibia: -0.72). As a consequence, the age-correlated coefficients, including demographic characteristics and muscle strength, are also statistically relevant to the first mode. For all three bones, the correlation with demographic characteristics is moderate to strong. The muscle strength values are weakly to strongly related to the first mode, except the tibialis posterior strength for calcaneus, where the ρ is nearly 0. No obvious correlation between the first mode and the dorsiflexion range.

The second mode of the calcaneus is strongly related to tibialis posterior strength, with a correlation coefficient of 0.95. The correlation with the strength of the other muscles and demographic characteristics are moderate to high, i.e., the second mode is similar and somehow complementary to the first mode. In the talus model, similar results are visible. The dorsiflexion range is mainly explained by the seventh mode of calcaneus and tibia, with a high correlation coefficient of over 0.5 and 0.64, respectively. For talus, it is related to the sixth mode (negative and moderate correlation) and the tenth mode (positive and strong correlation). FIGURE 3.24: Correlation matrix between clinical coefficients and shape parameters of calcaneus.

Discussion

SSM provides an access to study the morphometry of anatomical structures. As reported in the literature, it is widely applied to bony structures, including but not limited to femur [START_REF] Ebert | Reconstruction of Full Femora from Partial Bone Fragments for Anthropological Analyses Using Statistical Shape Modeling[END_REF], [START_REF] Asvadi | Bone Surface Reconstruction and Clinical Features Estimation from Sparse Landmarks and Statistical Shape Models: A Feasibility Study on the Femur[END_REF], [START_REF] Boutillon | Anatomically Parameterized Statistical Shape Model: Explaining Morphometry Through Statistical Learning[END_REF], [START_REF] Shi | Development of Predictive Statistical Shape Models for Paediatric Lower Limb Bones[END_REF], pelvis [START_REF] Vanden | Virtual Anatomical Reconstruction of Large Acetabular Bone Defects Using a Statistical Shape Model[END_REF], [START_REF] Shi | Development of Predictive Statistical Shape Models for Paediatric Lower Limb Bones[END_REF], tibia [START_REF] Shi | Development of Predictive Statistical Shape Models for Paediatric Lower Limb Bones[END_REF], [START_REF] Schmutz | Anatomical Fitting of a Plate Shape Directly Derived from a 3D Statistical Bone Model of the Tibia[END_REF], fibula [START_REF] Shi | Development of Predictive Statistical Shape Models for Paediatric Lower Limb Bones[END_REF] and scapula [START_REF] Plessers | Virtual Reconstruction of Glenoid Bone Defects Using a Statistical Shape Model[END_REF], [START_REF] Boutillon | Anatomically Parameterized Statistical Shape Model: Explaining Morphometry Through Statistical Learning[END_REF], [START_REF] Salhi | Statistical Shape Modeling Approach to Predict Missing Scapular Bone[END_REF]. In [START_REF] Lenz | Statistical Shape Modeling of the Talocrural Joint Using a Hybrid Multi-Articulation Joint Approach[END_REF] and [START_REF] Krähenbühl | Morphologic Analysis of the Subtalar Joint Using Statistical Shape Modeling[END_REF], authors analysed the adult tibiotalar joint and subtalar joint morphometry with SSM tool. In the case of pediatric data, morphological analysis is still an open area that needs to be explored. An ankle bone morphometry study of typical developing children is presented in this chapter using a voxel-based SSM approach, which is a first step in understanding how cerebral palsy influences ankle joint morphometry. SSM correspondences are constructed using a voxel-based approach. By completing this step, an ankle joint atlas for children is created. Atlas depicts visually realistic structures of the calcaneus, talus, and distal end of the tibia. The quantitative evaluation of deformation fields with DICE indicated that images of all subjects had a high degree of accuracy in the deformation field generation procedure.

The analysis of shape variation reveals the anatomical variation in the calcaneus, talus and tibia of typically developing children. The quantitative evaluation of model compactness indicates that the model achieves a compact representation of ankle bone shape. 70% of variation is encoded in the first four principal components, and the first mode of all three bones contains the most variation and represents the overall volume. First four modes capture the shapes of subtalar joint surfaces on calcaneus and talus, tibiotalar joint surfaces on tibia and talus, cuboid facet on anterior calcaneus, and calcaneal tuberosity. The bone shape can be affected by joint stability or pathology. A flatter talar dome has been identified to result in a lower level of tibiotalar stability [START_REF] Krähenbühl | The Subtalar Joint: A Complex Mechanism[END_REF], [START_REF] Frigg | The Biomechanical Influence of Tibio-Talar Containment on Stability of the Ankle Joint[END_REF], and the chronic ankle instability flatter talar joint surfaces and a flattened calcaneal groundcontact surface [START_REF] Tümer | Typical Shape Differences in the Subtalar Joint Bones Between Subjects with Chronic Ankle Instability and Controls[END_REF]. The shape model with normative parameters can help clinicians locate the significant morphological changes for diagnosis.

Certain studies have accessed the relationship between shape parameters and anatomical coefficients [START_REF] Cerveri | Representative 3D Shape of the Distal Femur, Modes of Variation and Relationship with Abnormality of the Trochlear Region[END_REF], [START_REF] Sintini | Investigating Gender and Ethnicity Differences in Proximal Humeral Morphology Using a Statistical Shape Model: PROXIMAL HUMERAL MORPHOLOGY[END_REF]. Accordingly, under the hypothesis that demographic characteristics and muscle quality influence bone shape, the correlation between demographic characteristics, muscle strength and dorsiflexion is investigated. The first mode of all three bones is strongly correlated with age. The second mode of calcaneus and talus also has a correlation with age. A strong correlation between the second mode and tibialis posterior strength is visible for calcaneus. Visually, the second mode represents the shape of the facet for cuboid, where the plantar portion of the tibialis posterior inserts into, and the sustentaculum tali, at which the recurrent portion of the muscle inserts into. The model has interpretability from a clinical perspective, which offers the possibility to detect the location of functional disability and facilitate the rehabilitation setting.

Limitations of this study are three folds. Firstly, the generalization and specificity do not converge with all principal components. It indicates the requirement for a larger dataset to improve the model's robustness. A large number of data will statistically facilitate the model's learning of more accurate distributions, enhancing the model specificity and introducing more shape variation, improving generalization. Models with sufficient specificity and generalization can handle generative tasks such as predicting the missing part of a shape based on partial observations and augmenting data with the model. Secondly, due to children's growth, the bone shape is variable longitudinally. Even if we take growing cartilage into account, it does not accurately represent the accurate picture of bone development. With real 3D pediatric scanning of a large population, longitudinal modeling of ankle bone shape will provide information on ankle joint development. Thirdly, only the distal end of the tibia is analyzed in this thesis. By taking the entire tibia into consideration, the shape will be better modeled.

The purpose of this chapter is to present a voxel-based model of the typical developing child's ankle joint. An atlas is constructed, which can be used as an anatomical reference in future studies. The shape variation among healthy children is analyzed and the correlation with clinical measurements is studied. It can provide information for clinical routines, such as the diagnosis of some musculoskeletal disorders. In the following chapters of this thesis, we will continue our investigation in ankle bone morphometry based on the data and results obtained in this chapter. The inter-population morphological difference will be accessed.

Introduction

After the normal bone shape of the ankle joint is studied and the atlas of children's ankle joint is constructed, these deformities need to be characterized and visualized through the anomaly analysis of the shape of the ankle bones. Therefore, this work aims to provide informative bone shape maps.

To perform the morphometry analysis with the volumetric representation of shapes, two types of methods based on the images have been developed, according to the scale to study: region-scale analysis and voxel-scale morphometry. For region-scale analysis, some anatomically defined regions of interest are manually labeled, therefore it has the strength of anatomical validity. However, the time consumption and the reliability of manual labeling need to be treated with care. The voxel-scale morphometry is an automated method that examines the whole anatomical structure with a voxel-wise measurement to localize the morphological change. Researchers have performed comparative studies of these two methods in the detection of gray matter abnormalities [START_REF] Seyedi | Comparing VBM and ROI Analyses for Detection of Gray Matter Abnormalities in Patients with Bipolar Disorder Using MRI[END_REF], [START_REF] Giuliani | Voxel-Based Morphometry versus Region of Interest: A Comparison of Two Methods for Analyzing Gray Matter Differences in Schizophrenia[END_REF]. In [START_REF] Seyedi | Comparing VBM and ROI Analyses for Detection of Gray Matter Abnormalities in Patients with Bipolar Disorder Using MRI[END_REF], authors report that the voxel-scale method is more sensitive to morphological change and it is highly recommended to use voxel-scale morphometry. However in [START_REF] Giuliani | Voxel-Based Morphometry versus Region of Interest: A Comparison of Two Methods for Analyzing Gray Matter Differences in Schizophrenia[END_REF], authors suggest that the two methods provide different types of information and should be used in tandem. Therefore in this study, both region-scale analysis and voxel-scale method will be explored, in order to see the consistency and the difference between these two methods.

In the context of shape analysis, an anomaly can be defined by local values of deformation fields that deviate significantly from the normal distribution to the template space. The main classes of algorithms to perform anomaly detection are probabilistic models (such as density estimation or neural generative models), one-class classification, and reconstruction models [START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]. One key element in the context considered in this work lies in the small sample size driving the choice of appropriate detection methods.

In this chapter, we investigated the bone morphometry difference between TD and CP populations from a multi-scale view. From a global scale, the analysis is performed with bone volume comparison. For region-scale analysis, some anatomical regions is defined and the local distribution of deformation characters is analysed. The voxel-wise analysis inspects the morphological difference with a voxel-scale anomaly detection from voxel scale, by using a comprehensive voxel-scale ankle joint bone shape analysis pipeline between TD population and CP population with equinus. Through this study, we initially analyzed abnormal changes in bone morphometry that occur in children with equinus. This study allows us to understand further the pathological mechanism providing a targeted reference for individualized clinical decision-making.

Methods

The clinical data and its pre-processing, the template estimation and the deformation field computation are presented in Chapter 3. Thus in this section, we skip these steps and focus on the morphological anomaly analysis in global, regional and voxel scales.

Global scale analysis -volumetric measurement

The bone volume is an anatomical indice that is related to muscle quality and thus the function [START_REF] Vicente-Rodríguez | How Does Exercise Affect Bone Development during Growth?[END_REF], [START_REF] Kindler | Skeletal Muscle and Pediatric Bone Development[END_REF]. In this study, we compare the volume of calcaneus and talus, including their growth cartilage, according to age. The quantification is performed with the software itk-SNAP1 .

Region-scale morphometry analysis

To better understand the bone deformity on a regional scale, and obtain the morphometry information based on the anatomically valid regions, we divide calcaneus in 7 regions: lateral, medial, superior, inferior, cuboid facet, subtalar facet and posterior region; talus in 5 regions: lateral, medial, superior, inferior and head region; tibia in 3 regions: fibula notch, malleolus, distal surface regions (see Figure 4.2). The regions are defined based on the anatomy of the bones of interest [START_REF] Dufour | Ostéologie[END_REF]. In order to find out the difference in each region from a statistical perspective and based on TBM, for each region, we analyze the regional log-Jacobian determinant (called log-Jacobian) difference by measuring the Wasserstein distance of the log-Jacobian distribution of TD and CP groups. Wasserstein distances are metrics to measure the distance between probability distributions that are inspired by the problem of optimal transportation [START_REF] Victor | Statistical Aspects of Wasserstein Distances[END_REF]. It represents the minimum energy cost to transform a pile of soil from the shape of one probability distribution to the shape of the other distribution, and it is also called Earth Movers distance. Compared to the other distance to measure two distributions, such as Euclidean distance and Hellinger distance, the Wasserstein distance can measure the distance not only of continuous distributions but also of discrete. Moreover, it is able to preserve the geodesic structures when transforming one distribution to the other one.

For probability measures µ and ν, x and y are d-dimensional vectors distributed as µ and ν, the p-Wasserstein distance is defined as W p (µ, ν) = inf(E∥x, y∥ p ) 1/p . If the sample vectors are 1-dimensional and p = 1, the Wasserstein distance can be calculated as W 1 (x, y) = R ∥F x (t), F y (t)∥dt, where F x (t) and F y (t) are cumulative distribution functions of x and y.

Voxel-scale anomaly detection

To focus on bone shapes, the deformation field analysis is restricted on bone edges (computed using Sobel filters). The overall objective is to provide a visually interpretable anomaly map robust to the limited size of the data sets. Given two deformation fields (denoted as {D TD i } i=1,2,...n for TD group and {D CP i } i=1,2,...m for CP group), we investigate 2 complementary anomaly detection techniques, one based on a statistical approach and another based on a reconstruction criterion.

Statistical map

The statistical map is based on Hotelling's T 2 test. Hotelling's T 2 distribution is the multivariate analogue of Student's t-distribution. The T 2 score of a subject i at location x is given as follows:

T 2 (x) = [d i (x) -µ TD (x)] T Σ -1 TD (x)[d i (x) -µ TD (x)] (4.1)
where µ TD and Σ TD are the mean and covariance matrix of the deformation fields based on the TD samples, and d i is the deformation field of the subject i. To provide statistical maps highlighting significant differences in deformation fields, a p-value map is calculated. Using the following link between the Hotelling's T distribution and the F-distribution ( n TD -3 3(n TD -1) T 2 (x) ∼ F 3,n TD -3 ), we can compute the p-value corresponding to a given location, which comes from the F-distribution.

Reconstruction-based map

The reconstruction-based map is computed using kernel PCA [START_REF] Hoffmann | Kernel PCA for Novelty Detection[END_REF], which is a non-linear extension of PCA algorithm introducing a non-linear kernel function ϕ(x), to enable the decomposition and reconstruction of non-linear data such as non-linear deformation fields. The kernel PCA model is trained with deformation fields of TD subjects {D TD i } i=1,2,...n . Firstly, D TD i is decomposed from original space to feature space with an orthogonal principal basis W with the following objective:

max n ∑ i=1 ∥Wϕ(D TD i )∥ 2 s.t. WW T = I (4.2)
Then the reconstruction of D TD i is performed as follows:

min n ∑ i=1 ∥ϕ(D TD i ) -W T Wϕ(D TD i )∥ 2 s.t. WW T = I (4.3)
In this work, ϕ is a Gaussian kernel. Once a model is trained, it can be applied on the deformation field D i of each subject to obtain a model-based reconstructed field Re(D i ). The reconstruction error of subject i at location x is defined as:

E i (x) = 3 ∑ j=1 ∥D i (x j ) -Re(D i (x j ))∥ 2 (4.4)
From a biomechanical perspective, such an error map provides physical information (reconstruction error of the deformation field in mm), complementary to the statistical approach previously described. From an algorithmic point of view, kernel PCA is implemented with scikit-learn library 2 . In order to combine this information with the statistical detection based on the Hotelling T 2 test, the reconstruction error is transformed into p-values (i.e. by converting the error into a z-score using the deformation fields of the TD subjects and then by transforming the z-score into a p-value).

Combined anomaly map computation

To provide an anomaly map of significant bone deformations related to cerebral palsy from this limited-size dataset, we propose to combine the statistical maps with the reconstruction errors calculated from the kernel PCA decomposition. More specifically, we propose to combine the two p-value maps (denoted as p stat and p recon ) using Fisher's method.

We compute a voxel-wise anomaly map of every TD subject and the population anomaly map denoted as AnoMap is then derived by the mean of individual maps, as follows:

AnoMap(x) = 1 n CP n CP ∑ i=1 F (p stat , p recon )(x) (4.5)
where F is the Fisher's method which combines the p-value maps into one χ 2 test statistic with the formula χ 2 2 sin -2log(p stat p recon ). The p-values map is corrected for multiple tests using FDR3 .

Filtering of deformation fields

Deformation fields provide the local transformation map describing the voxelwise volumetric and positional change of subjects. However, due to the geometric constraints, the alignment process is challenging to give an utterly optimal mapping. There may still be a residual error between the aligned source and reference images or the regional discontinuity in the deformation fields. On the other hand, the disease may impact one region but different voxel on a subject scale. Therefore as mentioned in [START_REF] Studholme | An Intensity Consistent Filtering Approach to the Analysis of Deformation Tensor Derived Maps of Brain Shape[END_REF], the analysis scale is one key point to reveal a consistent pattern among the population while it is hard to decide. In this section, we investigate the application of filters on deformation fields and their log-Jacobian determinant.

Gaussian filter

Gaussian filter is a filter widely used for signal smoothing and denoising. Its filter function is a Gaussian kernel f

(x) = 1 √ 2πσ exp( -x 2 2σ 2
). The advantage of the Gaussian filter in statistical morphology analysis is that it can force the deformation fields towards Gaussian random fields, as described in [START_REF] Gaser | Deformation-Based Morphometry and Its Relation to Conventional Volumetry of Brain Lateral Ventricles in MRI[END_REF], [START_REF] Studholme | An Intensity Consistent Filtering Approach to the Analysis of Deformation Tensor Derived Maps of Brain Shape[END_REF].

Statistics filter

One disadvantage of the conventional Gaussian filter is the tendency to lose edge information, which may cause trouble in distinguishing different tissues and identifying different anatomical structures. The statistics filter considers the neighborhood statistics by weighting the original image with the local mean and variance [START_REF] Lee | Speckle Analysis and Smoothing of Synthetic Aperture Radar Images[END_REF]. The estimated despeckled value f of neighborhood (x, y, z) can be described as:

f (x, y, z) = f (x, y, z) + k( f (x, y, z) -f (x, y, z)) (4.6)
where f (x, y, z) is the mean of the voxel in the window and k is a weighting function:

k = σ 2 (x, y, z) σ 2 (x, y, z) + f (x, y, z)σ 2 (4.7)
σ 2 (x, y, z) is the local variance and σ 2 is the variance of whole image. With the help of local statistics, statistic filter is able to enhance the edge effect.

Experiments and results

Global scale analysis -volumetric quantification

In the first step, we perform a global bone volumetric quantification of calcaneus and talus, from a global perspective, to provide a global view of bone morphological differences caused by equinus. The mean calcaneus volume of TD and CP group are 48172.9mm 3 and 41001.6mm 3 , and for talus are 29404.1mm 3 and 25138.3mm 3 . Both calcaneus and talus volume are decreased by 15% in CP group compared to TD group. The gap between the two groups is increased with age, as presented in

Region-scale analysis

For region-scale analysis, the Wasserstein distance is employed to measure the log-Jacobian distribution distance between the two groups. As reported in Table 4.1, the distance between TD and CP in regions ranges from 0.010 to 0.083. The cuboid facet region of calcaneus possesses the most discrepancy, with a distance equal to 0.083. The maximum dissimilarity of the talus is observed in the lateral region with a distance of 0.060. The distance in the inferior region of the talus is also more notable than the other regions (0.041). In the tibia, the distance in all regions is less than 0.018.

Voxel-scale anomaly analysis

Finally, to locate morphological abnormalities more accurately, we performed morphological anomaly analysis from the voxel scale using our pipeline. The optimization of the registration is based on not only the image intensity but also the segmentation map and the distance map. The different segmentation map and Chamfer distance map may result in a different final deformation fields. In the first step, to concentrate on the shape of one bone, the independent bone shape analysis is performed. The independent bone segmentation and distance map is able to well perserve the single bone shape, and may leads to a more precise anomaly analysis. However, if we use the joint bone shape analysis, i.e. the construction of atlas and the deformation field computation is executed with the mask of the three bones together, instead of one bone at one time, the anatomical information of the subtalar joint and the tibiotalar joint can be taken into consideration. Therefore, we also test the joint analysis and compare the results with the independent bone analysis.

Independent bone shape analysis

In the first step, the independent bone shape analysis is performed. In this section, the atlas is created separately for each bone. The registration and the generation of deformation fields, as well as the anomaly analysis, is only for one specific bone at one time.

Figure 4.4 shows our pipeline's voxel-wise anomaly analysis results. The voxelwise analysis is executed with all CP subjects and the atlas of the CP population, i.e., the results are given by the mean of anomaly analysis between TD subjects and each CP subject and between TD subjects and the CP atlas. The analysis based on the subjects' mean and the atlas demonstrates consistency. The areas with a low pvalue are similar in the two types of analysis. The possible impacted areas reported by the pipeline in the lateral and the anterior view of the calcaneus are the lateral margin of the facet for cuboid, and the anterior, middle and posterior facet for the subtalar joint. In the medial and the posterior view, the medial side of calcaneal tuberosity is mildly marked, and some scattered areas close to the sustentaculum tali. In the talus, the lateral process, the head and the anterior facet for calcaneus is the most marked region, as well as the head, neck and the anterior facet for calcaneus. Besides, some scattered regions on the dome region are observed only in the subject analysis, which was the only difference between the results of the subject and the atlas analysis. In the tibia, few anomalies are reported. The voxel-wise analysis corresponded to the regional analysis that the possible anomalies are located in the anterior calcaneus, and the lateral talus, and there is approximately no anomaly in the tibia.

Joint analysis

For joint analysis, the atlas is constructed with the segmentation and distance map of three bones instead of 3 bones independently. The generation of deformation fields is also performed with the whole joint and then the anomaly analysis is independently executed for each bone.

It can be observed in Figure 4.5 that the joint analysis tends to detect fewer possible anomalies than the independent analysis. The morphological changes reported by joint analysis are the lateral margin of the facet for cuboid and the lateral margin of the posterior facet of the subtalar joint surface in the calcaneus. In the talus, the highlighted region is the lateral process. All possible anomalies are less significant than independent analysis. The joint registration optimizes the entire region of interest, i.e., the three bones. Thus the inner anatomical margin, such as the subtalar surfaces for calcaneus and talus and the tibiotalar surfaces for talus and tibia, is not the optimization priority, which may cause the mixing of anatomically distinct but geometrically close voxels. It may indicate that the joint scale is too coarse to analyze the bone morphological changes in our case.

(A) Independent Analysis (B) Joint Analysis

Effect of filtering

To examine the effect of the filtering, we apply the two filters on the deformation field of all subjects and the CP atlas to the TD atlas and its log-Jacobian map. The full width at half maximum (FWHM) of the gaussian kernel is ranged from 4mm, 8mm and 12mm.

Visualization of filtered log-Jacobian map and deformation field magnitude

To examine the smoothing effect and contrast provided by filtering, we visualize the log-Jacobian map and the magnitude of the deformation field between the CP atlas and the TD atlas of the talus. The results are illustrated in Figure 4.6. Both Gaussian filtering and statistical filtering achieved the smoothing of images. Gaussian filtering has a more substantial smoothing effect and is able to effectively remove noise, while the loss of edge information is also visible. The bone edge is difficult to be identified after Gaussian filtering, and the loss of high-frequency information is more significant as the number of neighborhood voxels increases. The statistics filtering is less efficient in removing noise, as the discontinuity of the deformation field still exists after smoothing. However, the edge of the bone is better perserved. In addition, increasing the neighborhood area does not change the statistics filtering results significantly.

Effect on anomaly detection of filtering

We applied the anomaly detection pipeline on the original deformation fields and the filtered version to explore how the filtering influences the detection results.

Discussion

The aim of this study was to investigate the bone morphology relevant to fixed equinus caused by cerebral palsy. To this end, a global volumetric quantification, a region-scale and a voxel-scale analysis were performed. Experiments revealed the abnormal deformation pattern of ankle joint bone which is the first illustration to our best knowledge. Firstly, the results of global volumetric measurement reported that the volume of calcaneus and talus in CP population is 15% smaller than TD population. Secondly, in the results of region-scale and voxel-scale analysis, the morphometry in region scale and voxel scale is analyzed.

Global volumetric analysis

The growth of the rearfoot bone appears to be limited in children with cerebral palsy and fixed equinus ankles. A similar result was reported by Wren et al. [START_REF] Al | Bone Density and Size in Ambulatory Children with Cerebral Palsy: Bone Density and Size in CP[END_REF], who reported a smaller vertebral size in children with cerebral palsy compared to TD population. The volumetric reduction of ankle short bone may indicate the growth problem of children with cerebral palsy. One point of further study is to confirm this hypothesis and figure out the factors relevant to volumetric reduction, which will provide a guide for prevention and treatment.

Region-scale and voxel-scale anomaly analysis

The region-scale and voxel-scale morphometry provided the results with consistency. In region-scale analysis, the Wasserstein distance reports the cuboid facet of calcaneus and lateral talus are the regions with more morphological changes. In voxel-scale analysis, the possible abnormal areas are localized at the anterior-lateral calcaneus, posterior facet of subtalar joint of calcaneus, the talar neck and lateral process of talus. Two methods can provide mutually complementary information. The region-scale analysis provides a view from a more coarse scale and has an anatomical explanation. The voxel-scale analysis provides more precise information in an accurate and efficient way.

Equinus deformity is primarily caused by contracture of the gastrocsoleus muscle. These muscles are connected to the calcaneus through the Achilles tendon at calcaneal tuberosity and are involved in plantar flexion of the ankle [START_REF] Horsch | Prevalence and Classification of Equinus Foot in Bilateral Spastic Cerebral Palsy[END_REF]. Thus the calcaneal tuberosity may be impacted by pathology and corresponds to the voxel-wise analysis results on the medial and posterior view of calcaneus. The tight gastrocsoleus muscle also causes the compensation of the foot or leg during gait. During this compensation, the hindfoot severely pronates through the subtalar joint, allowing the midtarsal joint to unlock and dorsiflexion to increase through the oblique axis of the midtarsal joint [START_REF] Gourdine-Shaw | Equinus Deformity in the Pediatric Patient: Causes, Evaluation, and Management[END_REF]. The pronation of subtalar joint may finally lead to the bone deformation with growth, as detected in calcaneus and talus. Researchers have reported that the compensation in the sagittal plate through dorsiflexion at the naviculocuneiform joint, leading to pes planus [6], [START_REF] Johnson | Biomechanics of the First Ray Part V: The Effect of Equinus Deformity[END_REF], [START_REF] Perry | Gait Analysis: Normal and Pathological Function[END_REF], which may cause the deformation at the talar lateral process and the calcaneal cuboid facet. The literature also reports that the equinus led to a change in tibiotalar angle and resulted in a small-dome and flatter talus deformation [START_REF] Shalaby | Correction of Complex Foot Deformities Using the V-Osteotomy and the Ilizarov Technique[END_REF], [START_REF] Kolb | The Impact of Different Types of Talus Deformation after Treatment of Clubfeet[END_REF]. This deformation was potentially detected on the tibiotalar surface of talus. In the future study, it would be interesting to explore the correlation between our results and the gait parameters such as muscle strength and spasticity, as well as dorxiflexion range, to figure out the interaction of muscle and bone morphology and how they affect gait from biomechanical aspect.

Filtering of deformation fields

The deformation fields describe the volumetric change in voxel scale. However, the scale of analysis is still difficult to decide, even if the regional analysis has been performed. Equinus may change the bone morphometry in a scale finer than voxel but coarser than the defined regions. In order to adjust the analysis scale, as well as to eliminate the registration error, two types of filter have been tested on the log-Jacobian maps and the deformation fields. Experiments report that Gaussian filter is efficient in noise removal while loses the edge information. The effect of statistic filter is non-significant. The filtering achieves the enhancement of deformation tensor map in brain shape analysis [START_REF] Studholme | An Intensity Consistent Filtering Approach to the Analysis of Deformation Tensor Derived Maps of Brain Shape[END_REF]. Hence, the brain is a more complex anatomical structure with different tissues. The filtering is able to take into account information from the neighborhood so that the analysis can be adaptively scaled to the size of a particular anatomical structure by adjusting the FWHM without providing a segmentation map. In our case, the bones are simpler and can be considered as a whole anatomy structure, i.e. few anatomical margin inside the bone. Moreover, we performed separate registration and independent shape analysis for each bone, and the registration is with segmentation map and distance map. It ensures the maintenance of bone shape edge to distinguish the different anatomical structure. Thus the filtering in our case obtained less effectiveness in our skeletal study.

Limitations and perspectives

Despite the fact that this study offers exciting new information concerning the bone morphology of children with cerebral and equinus, there are some limitations that must be considered in order to make a meaningful interpretation of the findings. First of all, children in TD group tends to have larger BMI than those in CP group. In the context of bone volume comparisons, these non-significant differences may have an impact on the results but it is hard to explain a large gap of 15% bone volume reduction. Lower BMI may be another consequence of the growth problem caused by cerebral palsy. Secondly, this study was performed on a limitedsize dataset due to the scarcity of pediatric datasets. The lack of data can cause the mis-estimation of distribution which will lead to an inaccurate detection result. Although we revealed some deformation pattern from our dataset, it is better to confirm our results by increasing the number of subjects included. Due to the limitation of acquisition settings of pediatric MRI, an alternative way to evaluate the analysis results is to perform the comparative study with SOTA methods. Thus in the next chapter, we will introduce the shape analysis with one SOTA tool and our method. Thirdly, the group mean difference is reported in this chapter. However, in clinical practice, intervention designed based on population analysis may not be adapted to each individual. Therefore, in the next chapter, we will present a comprehensive personalized shape analysis for each individual.

Introduction

To determine whether an intervention is effective, studies report mean differences between groups in clinical trial settings. However, even in the case of statistically significant mean group effects, the intervention may not be effective for every participant in the study [START_REF] Damiano | Meaningfulness of Mean Group Results for Determining the Optimal Motor Rehabilitation Program for an Individual Child with Cerebral Palsy[END_REF]. As pointed out in the review [START_REF] Damiano | Meaningfulness of Mean Group Results for Determining the Optimal Motor Rehabilitation Program for an Individual Child with Cerebral Palsy[END_REF], CP is a very heterogeneous population with a large variation in the disorder distribution, character and severity. It may also be the case that even for interventions deemed highly effective in CP, a range of individual responses may be observed, from a negative response to a strong positive effect. Studies [START_REF] Damiano | Can Strength Training Predictably Improve Gait Kinematics? A Pilot Study on the Effects of Hip and Knee Extensor Strengthening on Lower-Extremity Alignment in Cerebral Palsy[END_REF], [START_REF] Steele | Characteristics Associated with Improved Knee Extension after Strength Training for Individuals with Cerebral Palsy and Crouch Gait[END_REF] have reported that lower extremity strengthening can both alleviate abnormal gait posture and exacerbate it, depending on the individual degree of muscle spasticity. A mean group response may not be informative in individual medical decision-making.

The concept of personalized medicine is developed based on the measurement of biomarker -a characteristic that is objectively measured and evaluated as an indicator of normal biological processes, pathogenic processes, or pharmacologic responses to a therapeutic intervention [START_REF]Biomarkers and Surrogate Endpoints: Preferred Definitions and Conceptual Framework[END_REF]. In pharmacogenomics, it has been recognized that differences in response from benefit to serious side-effect exists and the importance of personalized medication is pointed out [START_REF] Heller | Genetics/Genomics and Drug Effects[END_REF]. The personalized medicine can be now generalized to the other domains of health care such as precision psychiatry [START_REF] Brisa | The New Field of Precision Psychiatry[END_REF]. However in the rehabilitation of musculoskeletal disorders caused by CP, personalized bone morphometry has not been adequately investigated. The large standard deviation of inter-subject responses to therapy in the clinic is ineligible and a personalized approach to medical care is currently required. In order to better respond to the various needs of patients, it appears necessary to develop tools that allow a personalized approach based on a fine morphological analysis to propose a dedicated, effective follow-up.

For this purpose, it is required to rely on a dedicated analysis of the morphology of the ankle joint for each subject with the help of the SSM. In the last chapter, the population scale shape analysis is presented. Such an analysis of the shape of the ankle bones provides access to identify and visualize the deformities in the view of the whole group. Thus the main objective of this chapter is to develop a personalized approach, which could provide a better understanding of CP at the individual level as well as a more suitable rehabilitation program for those individuals [START_REF] Kedem | Foot Deformities in Children with Cerebral Palsy[END_REF]. However, a recent study [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF] has shown the importance of the evaluation and validation of these tools in clinical applications. Specifically, this study compared three widely used state-of-the-art SSM tools, namely ShapeWorks [START_REF] Cates | ShapeWorks[END_REF], Deformetrica [START_REF] Durrleman | Morphometry of Anatomical Shape Complexes with Dense Deformations and Sparse Parameters[END_REF], [START_REF] Bône | Deformetrica 4: An Open-Source Software for Statistical Shape Analysis[END_REF] and SPHARM-PDM [START_REF] Styner | Framework for the Statistical Shape Analysis of Brain Structures Using SPHARM-PDM[END_REF]. The quantitative and qualitative results show that the SSM tools have different levels of consistency and different abilities to capture variability at the population level. What becomes apparent through this study is the need to compare results obtained using multiple shape analysis methods.

In order to propose a customized approach, we focus on tools to provide informative deformation maps at a population level and also for each CP child based on shape analysis approaches. Morphometry, which is the study of the geometry of shapes, can be performed using voxel-based methods or surface-based approaches. In this chapter, we propose to investigate both types of approaches to provide the most comprehensive analysis possible at the patient level. Shape analysis relies on mapping (also called registration or matching) between subjects (or templates). Using a reference template, information about the individual shapes can be encoded in the deformation fields [START_REF] Ashburner | Identifying Global Anatomical Differences: Deformation-Based Morphometry[END_REF]. In such a context, personalized shape analysis is the study of the deformation fields for each subject.

In this chapter, we study the shape of the ankle joint of children with CP from high-resolution MRI data using two shape analysis approaches. More specifically, our contributions are three-fold: 1) a group analysis using two SSM methods, 2) an analysis of three ankle bones for a complete visualization of the joint, and 3) a subject analysis for a fine study of the deformation patterns for each child.

Methods

The overall flowchart is shown in Figure 5.1. In this work, the analysis is performed both on the surface-based approach and the voxel-based approach.

Surface-based shape analysis

Statistical Shape Modeling is a mathematical approach to quantify 3D shape variation. One approach for SSM is to analyze meshes computed for image segmentation. In this work, we rely on the surface-based SSM software called Shapeworks [START_REF] Cates | ShapeWorks[END_REF] 1 that has been recently used for ankle joint analysis using weightbearing computed tomography [START_REF] Krähenbühl | Morphologic Analysis of the Subtalar Joint Using Statistical Shape Modeling[END_REF], [START_REF] Lenz | Statistical Shape Modeling of the Talocrural Joint Using a Hybrid Multi-Articulation Joint Approach[END_REF]. A recent benchmarking study has shown the potential of Shapeworks with respect to other surface-based SSM [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF], which has been introduced in Section 2.2.3.1. As presented in Section 2.2.3, ShapeWorks is a groupwise particle-based shape modeling method that does not rely on surface parameterization. Shapeworks handles each surface mesh as a set of particles that describe the surface geometry. Such a particle-based representation avoids many of the problems inherent in parametric representations (i.e. limitation to specific topologies for instance). Shapeworks takes as input binary segmentation of each bone of every subject. Correspondences between surfaces (relying on particles) are estimated using signed distance images. Procrustes analysis is used to remove scaling (i.e. size) from the shape modeling analysis. Mean shapes are generated for CP and TD groups and deformation fields between the two groups are used to study the shape differences. In this work, the calcaneus is modeled with 2048 particles and the other two bones are modeled with 1024 particles. The length of the tibia is normalized with a cutting plane perpendicular to the tibial shaft.

Voxel-based shape analysis

The details of our voxel-based approach are introduced in Chapter 2. As a reminder, we briefly present this approach in this section.

For the voxel-based approach, we make use of an image registration-based framework to compute a mean image model of the TD population by using a group-wise diffeomorphic algorithm [START_REF] Avants | Geodesic Estimation for Large Deformation Anatomical Shape Averaging and Interpolation[END_REF] with image intensity, bone segmentation maps and signed distance maps. Then, for each CP subject of the dataset, deformation fields are estimated by non-linear multivariate registration onto the mean TD template previously computed. The template estimation and the patient-to-template registration stages are performed with ANTS2 .

For the group level analysis, a template of CP group is constructed in the same way as the TD template. The deformation field represents the correspondence between two groups, which is computed by registering the CP template to the TD template.

Results

To analyze the similarity and the variation between two different SSM methods, in this section, we first compare the results obtained in each bone at group level, and visualize the whole joint results in order to provide an overall view. To further understand bone morphometry from a personalized view, a subject-level analysis is performed. The subject-level analysis aims also to investigate the difference between group mean results and personalized results.

Analysis at group level

Our first objective is to provide a comparison between the two groups of interest (CP and TD), with the voxel-based and particle-based methods. Figure 5.2 shows the magnitude of the deformation fields for the 3 bones of interest from the CP group toward the TD group. These deformations correspond to the average shape variations between the two populations. This figure indicates that the two SSM methods studied in this work provide similar results regarding the main deformation regions of talus and calcaneus. These results tend to show that both SSM methods capture the same patterns of shape variation at the group level for these two bones.

Analysis at subject level

The group-level analysis shows shape differences in some regions of the bones of interest. However, the group-level analysis only provides average deformation patterns, without allowing a personalized analysis for each child. An analysis at the individual level is necessary to propose a more personalized approach. Figures 5.4, 5.5, 5.6 show the magnitude of the deformation fields between each CP subject (ordered by age) and the TD atlas using the voxel-based method for respectively, the calcaneus, the talus and the tibia in one view. The other views can be found in Appendix A. The main deformation patterns revealed by the group analysis may be observed in some subjects, but not necessarily in all subjects. Age does not seem to be related to the observed deformation patterns.

Discussion

Motivated by the observed volume difference of calcaneus and talus between TDs and CPs, this study is aimed to investigate the ankle joint bone morphology relevant to fixed equinus caused by cerebral palsy. To understand population-wise pathological change and to adjust patient-adapted rehabilitation program, the morphological analysis was performed on both on population level and individual level. To this end, we make use of two SSM methods. Such an analysis provides the possibility to analyse the morphological properties and detect large deformations than can be caused by pathology [START_REF] Ruiter | Automatic Image Matching for Breast Cancer Diagnostics by a 3d Deformation Model of the Mamma[END_REF], [START_REF] James | Detection and Characterisation of Bone Destruction in Murine Rheumatoid Arthritis Using Statistical Shape Models[END_REF].

Cross-approach comparison

In [START_REF] Gao | Shape Analysis, a Field in Need of Careful Validation[END_REF] and [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF], authors described different levels of consistency between SSM tools and stressed the importance of validating these tools in medical applications. Different SSM methods can lead to results with large variations. Under these conditions, a comparison between several shape analysis methods is necessary.

In this study, we adopt a surface-based approach and a voxel-based approach. Both methods capture similar group-level deformation regions for talus and calcaneus. Indeed, both methods reveal deformations of the anterior lateral edge of the calcaneus, which is located near the cuboid facet and the anterior facet of the subtalar joint, as well as the posterior-lateral edge of the posterior facet of the subtalar joint and the lateral process of the calcaneal tuberosity. In the talus, the area showing obvious deformity in both methods is the lateral process, as well as the neck of the talus. Compared to the calcaneus and the talus, the results for the tibia show less similarity between the methods. However, deformities of the infero-lateral part of the tibial malleolus are evident, which in relation to the deformities found in the calcaneus and talus can be explained by the valgus deformity of the hindfoot frequently found in children with CP and equinus gait. This deformity, exacerbated by weight-bearing, causes dislocation of the hindfoot with malalignement of the ankle bones, increasing local mechanical stress [START_REF] Otjen | Foot and Ankle Musculoskeletal Imaging of Pediatric Patients With Cerebral Palsy[END_REF]. The increase in stress is consistent with the deformities seen in the lateral-anterior border of the calcaneus, the lateral process of the talus and the tibial malleolus. However, differences still exist between two methods. In calcaneus, the calcaneus tuberosity and its medial process have large deformation magnitude in surfacebased approach, while this change is not reported by voxel-based approach. In talus, the lateral side of posterior facet is reported by voxel-based approach but not by mesh-based approach. In tibia, the fibula notch region obtained with mesh-based method has a magnitude of 1mm while approximately 0.5mm in voxel-based analysis.

Subject-level analysis

The subject-level analysis reveals the deformation pattern of each CP subject. As presented in section 5.3.2, the subject-level analysis may correspond to the group analysis, but not necessarily.

Interestingly, the talus appears to be the bone with the least variation in intersubject deformation. At the same time, it appears to be the bone of the ankle complex with the most deformation. This may be related to its anatomical-physiological characteristics (absence of muscle insertion, poor vascularisation and its anatomical position at the crossroads of the mechanical constraints of the talus and tibia [START_REF] Dufour | Ostéologie[END_REF]. Just as it is more sensitive to osteochondral lesions in sportsmen and women [START_REF] Barbier | Osteochondral Lesion of the Talus: What Are We Talking About?[END_REF], it seems to be more sensitive to deformation in CP individuals with equinus gait. This suggests that the evolution of the shape of the talus should be monitored more closely.

Despite this trend, no precise deformation pattern is observed for any of the three bones of the ankle complex. Although individual #3 (age 9.87) appears to be closest to the deformation pattern described above at the group level and some individuals are close to it in some deformations, alternative patterns are found for each individual. This does not mean that these individuals do not have deformities that could be explained by their loaded foot position or their gait. Further exploration is required to establish these relationships. The presence of major deformities at the individual level means that individualised diagnostic imaging assessment, in conjunction with gait analysis, should be encouraged in order to provide specific and personalized treatment.

Perspectives

This study offers exciting new information concerning the bone morphology of children with cerebral palsy and equinus. However, it is necessary to continue such a study in order to confirm these results. Firstly, this study is performed on a limitedsize dataset due to the scarcity of pediatric datasets. Although the similar deformation pattern is revealed from our dataset with different tools, it would be better to confirm our results by increasing the number of subjects included. Secondly, as already highlighted in [START_REF] Gao | Shape Analysis, a Field in Need of Careful Validation[END_REF] and [START_REF] Goparaju | Benchmarking Off-the-Shelf Statistical Shape Modeling Tools in Clinical Applications[END_REF], this work confirms the need to analyze the results obtained with SSM methods with care and a future line of research must be the understanding of the potential variability of the results obtained by SSM tools. It is essential to provide simple and efficient tools for shape analysis in a personalized context.

Conclusion and Perspectives

Conclusion

Under the assumption that the fixed equinus deformity may impact the bone morphology of the ankle joint, in this thesis, we investigate the bone shape of the ankle joint of children with CP and fixed equinus deformity. More specifically, the scientific questions are three-fold: 1) how to define the normal bone shape of TD children; 2) is there any shape difference between TD and CP populations; 3) if yes, where is the difference located; 4) whether the results at patient-level are consistent with the results at population-level.

To answer the first question, Chapter 2 shows the shape of the calcaneus, talus and tibia with the help of SSM. We present a voxel-based model of the TD population. Through the model, an atlas of TD children's joints, i.e., the mean image of children's ankle joint MRI, is first constructed. The atlas represents the children's ankle joint anatomy and can be used as a reference in future studies. The shape variation among healthy children is analyzed, and the correlation with clinical measurements is studied, which provides information about the anatomical variation of the bones of interest and a qualitative explanation of the model.

To cope with the second and the third questions, we perform a multi-resolution analysis in Chapter 3, including global volumetric analysis, regional analysis and voxel-scale analysis, which reveals the morphological difference from different levels. The global volumetric measurement reports a volume reduction in CP group, indicating that the equinus deformity may cause the limited growth of the rearfoot bone. The region-based analysis and the voxel-scale analysis are based on a deformation-based approach. In region-based analysis, we compare the regional mean Wasserstein distance of log-Jacobian between TD and CP groups. It demonstrates that the cuboid facet of the calcaneus and the lateral surface of the talus are the regions with more significant differences than the other regions. For voxelscale analysis, we implement an anomaly analysis method adapted to a limited-size dataset that aggregates two different approaches. The aggregation enables the combination of statistics and biomechanics to understand the population-wise morphological difference. The analysis reports that the most possible areas impacted by the pathology are the anterior-lateral calcaneus and the lateral process of the talus.

Regarding the fourth question, we analyze the deformation pattern of equinus at a subject level in Chapter 4. Results have certain consistency, but the variation among subjects is more visible. No precise pattern is observed for any of the three bones, which suggests a personalized treatment in clinical practice.

Another contribution presented in Chapter 4 is cross-SSM-tool comparison. Two different SSM tools are applied to our data: our voxel-based approach and one stateof-the-art SSM tool, ShapeWorks. Inter-approach consistency can be observed, while differences still exist. Both methods report the cuboid facet of the calcaneus and the lateral surface of the talus but ShapeWorks points out the calcaneal tuberosity in addition. Therefore the results obtained with SSM tools need to be treated with care.

Perspectives

In vivo Kinematic Analysis Based On Dynamic MRI

Equinus deformity results in the abnormal deformation of bone and muscle morphometry, as well as ankle kinematics. In this thesis, the 3D shape of the calcaneus, talus and distal end of the tibia is analyzed. Future studies could investigate the inter-population difference in the aspect of kinematics. With the help of 3D+time dynamic MRI sequences, the dorsiflexion movement is captured, and the 3D kinematics is estimated (Fig. 6.1-6.3) [START_REF] Makki | In Vivo Ankle Joint Kinematics from Dynamic Magnetic Resonance Imaging Using a Registration-Based Framework[END_REF]. The analysis of these sequences will enable the understanding of the pathology in a dynamic view.

Interpretation of Anomaly Analysis Results

In this thesis, the possible area of calcaneus and talus impacted by the equinus deformity is reported. But It still lacks interpretation from an anatomical and clinical view. For a better understanding of the causes and the effects of the pathology, future studies should investigate the correlation with the clinical measurements and the gait parameters, in order to find the determinant coefficient affecting the bone shape.

Deep-learning-based Image Registration

The diffeomorphic registration algorithm succeeded in our MR image registration task. However, the execution time is not negligible. A faster registration method is highly required for the real-world application. Many one-step deep-learningbased transformation estimation algorithms have been proposed and successfully applied in the medical context [START_REF] Haskins | Deep Learning in Medical Image Registration: A Survey[END_REF]. In future studies, the deep-based registration methods need to be explored. Abstract: Cerebral palsy (CP), a common physical disability in childhood, often causes abnormal movement patterns and posture, resulting in an abnormality in musculoskeletal morphometry. However, the bone deformation occurring in the CP population remains unclear. To better understand the pathology and improve the rehabilitation of patients, the main objective of this thesis is to comprehensively analyze the morphological change in children with CP based on image studies. In the first step, the normal bone shape of typical developing (TD) children is modeled with statistical shape modeling (SSM). A mean shape template called atlas of TD population is created as a reference. To locate the possible abnormal area, we analyze the bone morphometry of the CP population in a multiresolution manner. The bone shape is analyzed from the global, region-based and voxelbased scales. For voxel-based scale analysis, in order to cope with the limited size of the pediatric dataset, we introduce a voxel-based anomaly detection method adapted to the small dataset. Moreover, due to the variability between individuals and the requirement of personalized medical care, we study bone shape at both populational and personal levels. Besides, we compare the results of our voxel-based approach and one state-of-the-art SSM tool called ShapeWorks. The comparison demonstrates the consistency between methods, confirming the reliability of the obtained results.

2. 1

 1 Possible choices for decision making in anomaly detection............................. 3.1 Demographic characteristics of the children included per group. .................. 3.2 Muscle strength, spasticity and dorsiflexion of TD and CP groups. ............... 3.3 DICE of region of interest between atlas and subject's segmentation ............ 4.1 Regional Wasserstein distance of calcaneus, talus and tibia............................. xv

FIGURE 1 :

 1 FIGURE 1 : Pipeline d'analyse des anomalies de l'articulation de la cheville basé sur la déformation pour les équins.

FIGURE 2 :

 2 FIGURE 2 : Analyse comparative de la forme de l'articulation de la cheville avec une méthode basée sur les surfaces et une méthode basée sur les voxels.

19 1. 2 20 1. 3

 192203 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Thesis overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Thesis organisation . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

1 39 2. 3 46 2. 4

 1393464 Clinical context . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 2.1.1 Cerebral palsy and spastic equinus deformity . . . . . . . . 23 2.1.2 Anatomy of ankle . . . . . . . . . . . . . . . . . . . . . . . . . 24 2.1.3 Musculoskeletal deformation in CP . . . . . . . . . . . . . . 25 2.1.4 Magnetic resonance imaging (MRI) in the study of ankle joint 28 2.2 Statistical shape modeling (SSM) . . . . . . . . . . . . . . . . . . . 30 2.2.1 Image registration . . . . . . . . . . . . . . . . . . . . . . . . 30 2.2.2 Voxel-based shape analysis . . . . . . . . . . . . . . . . . . . 38 2.2.3 Surface-based shape analysis . . . . . . . . . . . . . . . . . . Anomaly detection . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 2.3.1 Probabilistic models . . . . . . . . . . . . . . . . . . . . . . . 42 2.3.2 One-class classification . . . . . . . . . . . . . . . . . . . . . . 45 2.3.3 Reconstruction-based methods . . . . . . . . . . . . . . . . . Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

FIGURE 2 . 1 :

 21 FIGURE 2.1: Bones of ankle joint. The ankle joint consists of 4 bones: calcaneus, talus, tibia and fibula. In this thesis, we focus on the calcaneus, talus, tibia. Extracted from [66].

FIGURE 2 . 2 :

 22 FIGURE 2.2: Anatomy of calcaneus. Figure source: Case courtesy of DrMatt Skalski, Radiopaedia.org, rID: 23709[START_REF] Vajuhudeen | Calcaneus[END_REF]

FIGURE 2 . 3 :FIGURE 2 . 4 :

 2324 FIGURE 2.3: Figure source: Case courtesy of Dr Matt Skalski, Radiopaedia.org, rID: 31891 [16]

FIGURE 2 . 5 :

 25 FIGURE 2.5: Tendons of foot. Reproduced from [66].

FIGURE 2 . 6 :

 26 FIGURE 2.6: Ligaments of ankle joints. Reproduced from [66].

Figure 2 .FIGURE 2 . 7 :

 227 FIGURE 2.7: Example of static MR image: (A) sagittal plane, (B) coronal plane, (C) axial plane.

Figure 2 . 8 ,

 28 one example of real-time T1w FFE sequence of ankle joint dorsi-plantar flexion is presented.The main challenges of dynamic MRI are the low resolution of acquired images and the artifacts induced by motion, which require processing including the superresolution and the removal of artifacts.

FIGURE 2 . 8 :

 28 FIGURE 2.8: 4 frames of a sequence of dynamique MRI of ankle joint.

FIGURE 2 . 10 :

 210 FIGURE 2.10: Illustration of non-linear registration: the moving image is transformed by the deformation field under LDDMM framework to align with the reference.

FIGURE 2 . 11 :FIGURE 2 . 12 :

 211212 FIGURE 2.11: Overview of 2-iteration template-based atlasing methods.

FIGURE 2 .

 2 FIGURE 2.14: An illustration of density level α and its threshold τ α

FIGURE 2 . 15 :

 215 FIGURE 2.15: Anomaly detection methods arranged by depth of feature map and the principle of models. Based on the depth of feature map, the methods can be divided into conventional shallow methods and deep methods. Based on the principle of models, three commonused classes of approaches are identified: one-class classification, probabilistic models and reconstruction-based models. According to the feature map, algorithms can be divided in to traditional shallow methods and novel deep methods. Additionally, distance-based methods are available as a complement to these three groups. Figure extractedfrom[START_REF] Ruff | A Unifying Review of Deep and Shallow Anomaly Detection[END_REF]

FIGURE 2 . 16 :

 216 FIGURE 2.16: Example of anomaly detection with different methods: (a) Hotelling's T 2 test of the probabilistic model; (b) KDE of the probabilistic model; (c) OCSVM of one-class classification; (d) PCA of the reconstruction-based model; (e) Kernel PCA of reconstruction model.The yellow points are detected as normal samples and the black points are detected as abnormal samples. The parametric Hotelling's T 2 statistic is limited to an ellipsoidal unimodal density, which has difficulty in modeling double-moon data. KDE with an RBF kernel is more flexible. OCSVM with RBF kernel has a similar effect of KDE. PCA finds a linear orthogonal subspace that is able to deal with unimodal data yet is less efficient for bimodal double-moon data. Kernel PCA introduced a nonlinear kernel which enables an optimal reconstruction from non-linear components in input space.

51 3. 2 56 3. 3 56 3. 3 . 2 58 3. 3 . 3 58 3. 3 . 4 67 3. 4

 512563563258335834674 , three classes of algorithms exist: probabilistic model, one-class classification and reconstruction-based model. Due to the scarcity of pediatric MRI data, careful selection of the model is required. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Materials and methods . . . . . . . . . . . . . . . . . . . . . . . . . 52 3.2.1 Populations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 3.2.2 Clinical coefficient measurement . . . . . . . . . . . . . . . . 52 3.2.3 3D image acquisition and bone-of-interest segmentation . . 53 3.2.4 3D voxel-based shape modeling . . . . . . . . . . . . . . . . 54 3.2.5 Clinical coefficient correlation analysis . . . . . . . . . . . . Experiments and results . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.3.1 Evaluation of atlas construction and deformation field generation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Quantitative evaluation of model . . . . . . . . . . . . . . . . Shape variation analysis . . . . . . . . . . . . . . . . . . . . . Correlation analysis between shape parameters and clinical coefficients . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71

FIGURE 3 . 1 :

 31 FIGURE 3.1: MRI and bone of interest segmentation of each subject: (A) TD subjects, (B) CP subjects. The title of each sub-figure indicates the number of the subject and the age of the subject. In the segmentation mask: blue: calcaneus; yellow: talus; green: tibia.

FIGURE 3 . 2 :

 32 FIGURE 3.2: 3D rendering of each subjects: (A) TD subjects, (B) CP subjects. The title of each sub-figure indicates the number of subject and its age.

Figure 3 FIGURE 3 . 3 :FIGURE 3 . 4 :

 33334 FIGURE 3.3: Estimated ankle atlas of TD population using the voxelbased approach: (A) sagittal view, (B) segmentation, (C) 3D rendering.

FIGURE 3 . 5 :

 35 FIGURE 3.5: Quantitative evaluation of model: (A) compactness, (B) generalization, (C) specificity.

FIGURE 3 . 6 :

 36 FIGURE 3.6: Anterior view of first four eigenmodes of calcaneus: Each row represents one eigenmode and is arranged in order from top to bottom. In each row from left to right: reconstructed shape at -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , reconstructed shape at 3 √ λ i

Figure 3 .

 3 13; superior view: Figure 3.14; inferior view: Figure 3.15; anterior view: Figure 3.16; posterior view: Figure 3.17). The visualization of the other two views is in Appendix A. The first mode contains

FIGURE 3 . 7 :FIGURE 3 . 8 :FIGURE 3 . 9 :FIGURE 3 . 10 :FIGURE 3 . 11 :

 373839310311 FIGURE 3.7: Lateral view of first four eigenmodes of calcaneus: Each row represents one eigenmode and is arranged in order from top to bottom. In each row from left to right: reconstructed shape at -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , reconstructed shape at 3 √ λ i

Figure 3 .

 3 18; posterior view: Figure 3.19; distal view: Figure 3.20; lateral view: Figure 3.21; medial view: Figure 3.22). The other two views are presented in Appendix A. The first mode consists of numerous anatomical differences, including the global volume, the distal surface shape and the medial malleolus shape. At

FIGURE 3 . 12 :FIGURE 3 . 13 :FIGURE 3 . 14 :FIGURE 3 . 15 :FIGURE 3 . 16 :FIGURE 3 . 17 :

 312313314315316317 FIGURE 3.12: Lateral view of first four eigenmodes of talus: Each row represents one eigenmode and is arranged in order from top to bottom. In each row from left to right: reconstructed shape at -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , reconstructed shape at 3 √ λ i

FIGURE 3 . 18 :

 318 FIGURE 3.18: Anterior view of first four eigenmodes of tibia: Each row represents one eigenmode and is arranged in order from top to bottom. In each row from left to right: -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , 3 √ λ i

FIGURE 3 . 19 :FIGURE 3 . 20 :FIGURE 3 . 21 :FIGURE 3 . 22 :

 319320321322 FIGURE 3.19: Posterior view of first four eigenmodes of tibia: Each row represents one eigenmode and is arranged in order from top to bottom. In each row from left to right: reconstructed shape at -3 √ λ i , mean shape with deformation magnitude of 3 √ λ i , reconstructed shape at 3 √ λ i

FIGURE 3 . 23 :

 323 FIGURE 3.23: Correlation matrix betweens clinical coefficients.

FIGURE 3 . 25 :

 325 FIGURE 3.25: Correlation matrix between clinical coefficients and shape parameters of talus.

FIGURE 4 .

 4 FIGURE 4.1: voxel-scale ankle joint anomaly analysis pipeline for population with equinus.

FIGURE 4 . 2 :

 42 FIGURE 4.2: Illustration of regions of calcaneus, talus and tibia. (A): calcaneus: lateral (dark blue), medial (turquoise), superior (green), inferior (blue), cuboid facet (orange), subtalar facet (dark red) and posterior (yellow). (B): Talus: lateral (turquoise), medial (orange), superior (blue), inferior (dark red), head (yellow). (C): Tibia: fibular notch (red), malleolus (yellow) and distal surface (turquoise).

Figure 4

 4

FIGURE 4 . 3 :

 43 FIGURE 4.3: Bone and growth cartilage volumes according to age: (a): Calcaneus; (b): Talus.

FIGURE 4 . 4 :

 44 FIGURE 4.4: Anomaly detection results: The first row is the mean of each CP individual analysis and the second row is the analysis between TD atlas and CP atlas.

FIGURE 4 . 5 :

 45 FIGURE 4.5: Comparison between independent analysis and joint analysis.

FIGURE 4 . 6 :

 46 FIGURE 4.6: Filtered deformation field magnitude and log-Jacobian Map: The first row is the magnitude of filtered deformation field and the second row is the filtered log-Jacobian map. As displayed in Figure 4.7, the Gaussian filter removes the scatter regions, such as detected red points on the superior face of the talus. The detection of the Gaussian filter emphasizes the lateral process and anterior facet, which are also the most highlighted regions in the initial analysis. The hint on the talus head is also diminished. The statistics filter results are highly similar to the results without filtering. No significant difference is visible. The results of detection correspond to the results of visualization of deformation fields.

FIGURE 4 . 7 :

 47 FIGURE 4.7: Anomaly detection results of talus with original and filtered deformation fields: The first row is the lateral view. The second row is the superior view. The third row is the view of head. In each row, the sub-figures follow this order: (A) the original, (B) 4mm Gaussian filtered, (C) 8mm Gaussian filtered, (D) 12mm Gaussian filtered, (E) 4mm statistics filtered, (F) 8mm statistics filtered, (G) 12mm statistics filtered.

94 5. 4

 944 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 5.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 5.2.1 Surface-based shape analysis . . . . . . . . . . . . . . . . . . 91 5.2.2 Voxel-based shape analysis . . . . . . . . . . . . . . . . . . . 92 5.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 5.3.1 Analysis at group level . . . . . . . . . . . . . . . . . . . . . . 92 5.3.2 Analysis at subject level . . . . . . . . . . . . . . . . . . . . . Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 5.4.1 Cross-approach comparison . . . . . . . . . . . . . . . . . . . 94 5.4.2 Subject-level analysis . . . . . . . . . . . . . . . . . . . . . . . 96 5.4.3 Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97

FIGURE 5 . 1 :

 51 FIGURE 5.1: Comparative ankle joint shape analysis with surface-based method and voxel-based method.

Figure 5 .FIGURE 5 . 2 :FIGURE 5 . 3 :

 55253 Figure 5.3 shows a group-wise deformation zone in the lateral process of the talus and a deformation zone on the anterolateral aspect of the calcaneus. These two areas are anatomically opposite each other. A zone of deformity is also seen on the inferolateral aspect of the tibial malleolus.

FIGURE 5 . 4 :

 54 FIGURE 5.4: Subject-level shape analysis of calcaneus in lateral view of CP population.

FIGURE 5 . 5 :

 55 FIGURE 5.5: Subject-level shape analysis of talus in lateral view of CP population.

FIGURE 5 . 6 :

 56 FIGURE 5.6: Subject-level shape analysis of tibia in lateral view of CP population.

Contents 6 . 1 99 6. 2 100 6. 2 . 1

 6199210021 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . In vivo Kinematic Analysis Based On Dynamic MRI . . . . . 100 6.2.2 Interpretation of Anomaly Analysis Results . . . . . . . . . . 100 6.2.3 Deep-learning-based Image Registration . . . . . . . . . . . 100

FIGURE 6 . 1 :

 61 FIGURE 6.1: 3D normative kinematics of tibio-talar joint.

FIGURE 6 . 2 :

 62 FIGURE 6.2: 3D normative kinematics of calcaneal-tibial complex.

FIGURE 6 . 3 :

 63 FIGURE 6.3: 3D normative kinematics of subtalar joint.

FIGURE A. 2 :

 2 FIGURE A.2: Subject-level shape analysis of calcaneus in inferior view

TABLE 2 .

 2

1: Possible choices for decision making in anomaly detection: α and β correspond to type I error and type II error.

  1 h K( x h ) is scaled kernel function. As for Hotelling's T 2 test, examples of KDE with Gaussian kernel are available in Figure 2.16, subfigure (b).

TABLE 3 .

 3 1: Demographic characteristics of the children included per group.

		Age (years) Weight (kg) Height (cm)	Leg Length (cm)
	TD (mean±std) 10.38±1.95 38.89±13.21 146.55±14.75 77.50±9.11
	CP (mean±std) 10.35±2.39 32.7±6.65	140.39±11.01 74.50±6.78
	p-value	0.97	0.22	0.48	0.37
	Muscular Strength:

TABLE 3 .

 3 2: Muscle strength, spasticity and dorsiflexion of TD and CP groups.

			TD	CP	p-value
	Normalized muscle strength (N/cm)	Tibialis anterior Tricep sural Peroneus 1.40±0.42 0.52±0.42 0.00016 1.76±0.66 0.78±0.54 0.0015 2.82±0.83 1.31±0.86 0.00078
		Tibialis posterior	1.45±0.63 0.74±0.50 0.055
		Tibialis anterior	0±0	0±0
	Spasticity	Tricep sural	0±0	1.44±0.86 5.10 × 10 -5
		Peroneus 0±0	0±0
		Tibialis posterior	0±0	0.31±0.59 0.14
	Dorsiflexion (degree)	Knee extented Knee flexion	12.2±4.78 -3.75±6.94 9.94 × 10 -6 20.5±5.99 7.5±6.55 0.00028

different ages

(7, 10 and 12)

are manually segmented by an expert. A registrationbased label propagation method is then used to segment the subjects of the data sets, with manual correction if necessary. Newly segmented subjects are included at each propagation stage. Finally, smooth segmentation maps are computed by training a standard patch-based U-net network

[START_REF] Ronneberger | U-Net: Convolutional Networks for Biomedical Image Segmentation[END_REF]

trained on the expert-based segmentation maps and manually corrected if required. The MR images and their segmentations of all subjects are illustrated in Figure

3

.1. The 3D renderings of bones of interest are in Figure

3

.2.

TABLE 3 .

 3

		3: DICE of region of interest between atlas and subject's seg-
		mentation
	Subject TD01 TD02 TD03 TD04 TD05 TD06 TD07 TD08 TD09 TD10 TD11
	DICE	0.990 0.993 0.991 0.994 0.993 0.974 0.991 0.992 0.991 0.976 0.974
	Subject CP01 CP02 CP03 CP05 CP06 CP08 CP09 CP10 CP13
	DICE	0.997 0.994 0.992 0.975 0.997 0.996 0.996 0.993 0.995

TABLE 4 .

 4 1: Regional Wasserstein distance of calcaneus, talus and tibia

	Calcaneus	lateral	medial	superior	inferior	cuboid facet	subtalar facet	posterior
		0.023	0.014	0.023	0.019	0.083	0.025	0.011
	Talus	lateral 0.060	medial 0.010	superior 0.014	inferior head 0.041 0.023
	Tibia	distal surface malleolus fibula notch 0.018 0.017 0.017

http://stnava.github.io/ANTs/

https://itksnap.org/

https://scikit-learn.org/stable/

https://www.statsmodels.org

https://sciinstitute.github.io/ShapeWorks

http://stnava.github.io/ANTs/

Acknowledgements

This work is funded by ANR (AI4CHILD ANR-19-CHIA-0015-01), Philips, Fondation de lAvenir, Paris, France, and Fondation Motrice, Paris, France and Institut Mines Télécom Atlantique. thanks to my supervisors, Pr.

A Supplementary Views of Subject-level Shape Analysis

To complete the visualization of subject-level analysis, the other views except the lateral visualization angle are presented here.