Tristan Deruelle

Electromagnetic Compatiblity

Keywords: Magnetic Resonance Imaging, MRI, Elastography, Prostate, Cancer, Porosity, Poro-elasticity, Biot. LNE Reports

Prostate cancer is the second most prevalent cancer in men worldwide. It is suspected when the PSA density is high or/and the superficial prostate feels hard during digital rectal examination. Multiparametric MRI is now recommended prior biopsy when detecting for cancer. However, image interpretation is challenging, even for specialists, and brings many false-positive. Elastography is a technique to assess tissue stiffness by inducing small vibrations. It could provide a 3D map of the stiffness of the prostate. We believe that MR elastography could complement the current multiparametric MRI. Given prostate location and consitution, wave propagation is difficult though. The current work presents the design of a non-invasive wave generation device for the prostate. Then, a new field separation algorithm is presented. This algorithm provides a better estimation of the stiffness, and the correction of artefact generated by common vibrators. Finally, this algorithm can have applications in porous media. Indeed, in poro-elastic materials, a slow compression wave propagates. We observe such a wave in an agar gel, in a foam phantom, and in vivo in human kidney graft. In addition to the classic shear wave velocity estimation, it is now possible to estimate the compression wave velocity. This is an additional piece of information that the operator can use in its diagnostic. In the future, more porous parameters could be derived.

Remerciements

Cette thèse marque la fin d'un cursus scolaire, mais peut-être le début de l'apprentissage. Pour en arriver ici aujourd'hui, j'ai été épaulé par différentes personnes que je souhaiterais remercier ici. Tout d'abord, je voudrais remercier Aline Bel-Brunon, Emmanuelle Canet-Soulas, Sabine Bensamoun, Ingolf Sack, ainsi que Jean-Luc Genisson, d'avoir accepté de constituer mon jury et d'évaluer mes travaux. Ce fut un honneur d'échanger avec vous et d'avoir vos remarques et questions sur ce travail.

Je tiens ensuite à remercier Olivier Rouvière et Rémi Souchon, de m'avoir encadré au cours de ces dernières années. Olivier, merci pour tes explications toujours claires et tes remarques pertinentes. Merci d'avoir pu dégager un peu de temps lorsque nécessaire, malgré ton agenda surbooké. Rémi, merci de m'avoir guidé chaque jour tout en me laissant une grande liberté dans mes recherches. J'ai pu apprendre à remettre tous mes résultats en question et à douter de chaque paramètre, parfois pour le meilleur, ou bien le pire. Je ne commence maintenant plus de phrase par « normalement » ou « en théorie ». Mais cela m'a instauré une grande rigueur. Je retiendrai également nos voyages à Paris, tantôt pour les biscotos en train, tantôt en kangoo les pieds dans l'eau, mais aussi nos sessions à l'IRM le soir à finir à pas d'heure, en se demandant toujours si une mise à jour aurait de nouveau effacé toutes nos séquences.

En parlant d'IRM, je tiens à remercier Scott Kruse pour les séquences élasto que nous avons pu utiliser à l'Hôpital Edouard Herriot. Merci également à Ingolf Sack pour la séquence élasto que nous utilisons au CERMEP. Par ailleurs, merci à Danielle de nous avoir aidé, et surtout d'avoir supporté nos longues, répétitives et bruyantes sessions.

Merci à toutes les personnes du laboratoire qui ont apporté leur contribution de près ou de loin. Par exemple, merci Stefan et Johannes pour votre expertise en physique des ondes et mousses. Merci également Catherine d'avoir réussi à gérer l'université malgré mes commandes à répétition ou mes nombreux petits morceaux de contrat.

Un grand merci à tous mes collègues du LabTAU qui ont su rendre ces dernières années agréables, que ce soit lors des joyeuses pauses ou pour les réunions « mousse ». Cela passe aussi par des petits succès comme réussir à organiser et faire tenir le repas de noël dans notre petite salle de réu. A propos de festivités, je m'excuse au nom du groupe de vous avoir fait courir sous une canicule par plus de 40°C lors de la journée du labo. Le bon côté des choses, c'est que ça fait des souvenirs, certains s'en souviendront. Mention spéciale à mes camarades Elodie, Jade, Myléva, et Tristan avec qui j'ai suivi quelques formations, mais aussi beaucoup râlé. Merci encore à Myléva d'avoir été là et d'être la seule à savoir réanimer mon pc lorsque je perdais ma connexion à distance.

Merci à tous mes amis, proches ou à l'étranger, que je ne pourrais tous citer. Merci d'avoir été là. Nous nous sommes un peu moins vus ces derniers temps, mais je compte pouvoir rattraper ça maintenant.

Merci à toute ma famille qui a été présente tout au long de cette aventure, qui m'a soutenu, et permis de prendre quelques bouffées d'oxygène. Merci à ma compagne, Pauline, qui a été à mes côtés au quotidien. Merci à vous de m'avoir parfois forcé à sortir de mes équations ou de mon code, et de m'avoir changé les idées quand plus rien ne fonctionnait.

Enfin, merci à ceux que j'ai oublié.

Merci à tous, et merci lecteur.

Résumé

Le cancer de la prostate est le deuxième cancer le plus prévalent chez l'homme dans le monde. Il est suspecté par un test de PSA et/ou un ressenti plutôt dur lors d'un toucher rectal. Lors du dépistage, une IRM multiparamétrique est recommandée pre-biopsie. Malheureusement, l'interprétation des images n'est pas aisée, même pour des spécialistes, et fait apparaître des faux-positifs. L'élastographie est une technique permettant d'estimer la rigidité des tissus lors de l'induction de petites vibrations. Cette technique permettrait d'avoir une cartographie 3D de la dureté de la prostate. Nous pensons que l'élastographie par IRM peut aider l'IRM multi-paramétrique actuelle. De par la localisation et la consitution de la prostate, la propagation des vibrations est difficile. Ces travaux présentent la conception d'un dispositif non invasif de génération d'ondes, spécifique pour la prostate. Ensuite, un nouvel algorithme de séparation de champ est présenté. Cet algorithme permet une meilleure estimation de la rigidité et la correction d'artefacts induits par les vibreurs conventionnels. Enfin, cet algorithme peut avoir des applications dans les milieux poreux. En effet, dans les milieux poro-elastiques, une onde de compression lente se propage. Nous montrons la présence d'une telle onde dans un gel d'agar, dans une mousse poreuse, et in vivo dans le greffon renal. En plus de l'estimation classique de la vitesse des ondes de cisaillement, il est maintenant possible d'estimer la vitesse de l'onde de compression lente. C'est une information supplémentaire que peut utiliser le praticien dans son diagnostique. Dans le futur, des paramètres de porosité pourront être évalués.

Mots-clés : Imagerie par Résonance Magnétique, IRM, Elastographie, Prostate, Cancer, Porosité, Poro-élasticité, Biot.

Résumé substantiel

Le cancer de la prostate est le deuxième cancer le plus prévalent chez l'homme dans le monde. Un des premiers tests de dépistage est le toucher rectal. Au cours de l'acte, le praticien va palper la partie surperficielle, accessible, de la prostate. Un ressenti plutôt dur est suspect et mène à effectuer des tests complémentaires (taux de PSA, IRM, biopsies…). Récemment, l'IRM multiparamétrique (pondéré T2, diffusion, dynamique) est recommandée pré-biopsie. Cependant, l'analyse de ces images est très difficile, même pour des experts. Elle fait également ressortir de nombreux faux-positifs. L'élastographie est une technique qui permet d'estimer l'élasticité des tissus (module de cisaillement ou vitesse des ondes de cisaillement), et permet d'obtenir l'équivalent d'une carte 3D de la « rigidité » des tissus. L'ajout de cette nouvelle information aux données actuelles d'IRM multiparamétrique permettrait d'aider les radiologues dans l'interprétation des images, la localisation des foyers cancéreux, ainsi que l'évaluation de leur aggréssivité.

L'élastographie par IRM (ERM) repose sur trois étapes. La première est de générer des ondes mécaniques dans les tissus d'intérêt. Ensuite, l'IRM est utilisée pour capturer des instantanés du champ d'onde afin d'obtenir un film 3D de la propagation des ondes. Ces films peuvent-être acquis pour les déplacements dans les 3 directions de l'espace de manière à obtenir le vecteur complet de déplacement (trois coordonnées spatiales au cours du temps). Enfin, les propriétés mécaniques des milieux sont reconstruites grâce à différents algorithmes.

L'obtention de cartes de vitesse de cisaillement dans la prostate fait face à plusieurs défis. Le premier est d'avoir des ondes d'amplitude suffisante dans toute la région d'intérêt. Or, la prostate est entourée de tissus atténuants. Les ondes qui parviennent jusqu'à l'organe subissent déjà une perte en amplitude significative. De plus, la prostate est un organe relativement rigide, donc difficile à faire vibrer. La propagation d'ondes dans toute la prostate est donc ardue. Une solution pourrait être d'utiliser des fréquences plus basses, qui se propagent mieux. Cependant, la résolution des reconstructions dépend des longueurs d'onde présentes dans la région d'intérêt. Il y a donc un compromis à trouver entre propagation et résolution. Le second challenge est de ne considérer que les ondes de cisaillement. En effet, les ondes sont généralement induites par des dispositifs agissant principalement en compression. Le champ d'onde est la superposition d'ondes de cisaillement (d'intérêt), et d'ondes de compression. Ces dernières ont des longueurs d'onde trop grandes pour pouvoir être mesurées, mais sont néanmoins présentes et biaisent les mesures de vitesse de cisaillement. La méthode classique pour éliminer les contributions de la compression est d'appliquer l'opérateur mathématique rotationnel. Cependant, cette opération fait intervenir le calcul de dérivées. Le champ de déplacement dans la prostate étant déjà de faible amplitude, le SNR y est faible. Le calcul de dérivées sur ce champ bruité va d'autant plus dégrader le SNR, ce qui peut fortement détériorer les reconstructions des propriétés des tissus.

Ces travaux présentent la conception d'un dispositif de génération d'ondes mécaniques pour la prostate. Une attention particulière a été portée sur le confort du patient ainsi que sur la sécurité. De ce fait, une approche non invasive, de manière transpérinéale a été préférée. La source des vibrations est un actionneur piézoélectrique mécaniquement amplifié. Des tests de sécurité électrique et de compatibilité électromagnétique ont été conduits par un laboratoire indépendant. Les performances du dispositif ont été testées sur un fantôme élastique.

Dans un second temps, le problème de la compression est adressé. Un nouvel algorithme de séparation de champ a été développé. Basé sur le théorème de Helmholtz, il consiste à calculer les potentiels vecteurs ou scalaires associés aux composantes de cisaillement et compression en résolvant une équation de Poisson. Ces potentiels permettent ensuite de calculer le champ de cisaillement et le champ de compression. L'algorithme permet donc, à partir du champ de déplacement complet, de le séparer en un champ de cisaillement d'une part, et champ de compression d'autre part. Dans le cas de tissus purement élastiques, le champ de compression n'apporte pas d'information exploitable (vitesse extrêmement rapide, très difficile à mesurer). Le champ de cisaillement est obtenu sans introduire de dérivées supplémentaires. Grâce à cette méthode, les artéfacts liés aux vibreurs compressifs sont corrigés. Les valeurs de vitesse des ondes de cisaillement sont plus justes, et la résolution est légèrement impactée. La résistance au bruit est très bonne, et permet de reconstruire la géométrie des milieux à faible SNR. Ceci peut être très avantageux dans le cadre de l'ERM de prostate. L'algorithme devient d'autant plus intéressant dans les milieux poreux. Un milieu poreux peut être défini comme une matrice solide saturée avec un liquide. Sous l'action d'une contrainte, les deux milieux vont être affectés, et les deux milieux vont réagir l'un à l'autre. La théorie de Biot explique que dans ces milieux poro-elastiques, une seconde onde de compression se propage. Celle-ci a la particularité d'être lente, plus rapide que l'onde de cisaillement, mais du même ordre de grandeur. Ce type d'onde a intéressé les géologues et a été mis en évidence dans des milieux de type roche. Dans le domaine médical, elle a été observée dans les os, milieu pouvant s'apparenter à la roche, mais jamais (à notre connaissance) dans les tissus mous. Or, les tissus mous tels que dans le rein, foie ou cerveau sont des matrices cellulaires dont les membranes sont poreuses et qui sont irriguées par des vaisseaux et capillaires. Il paraît donc raisonnable de les considérer comme poreux. Nous montrons la présence de ces ondes in vivo dans le greffon rénal. En plus des cartes de vitesse des ondes cisaillement, une carte de vitesse des ondes de compression est reconstruite. La seule présence de cette onde indique que le milieu est poreux. Sa vitesse est un nouveau biomarqueur, une nouvelle information pour le radiologue. Ces vitesses sont reliées aux paramètres de porosité intrinsèques du milieu.

Ce sont autant de nouveaux potentiels descripteurs de l'état des tissus, qui, dans le futur, pourraient être évalués.

Table of figures

1

This chapter provides an introduction to prostate and prostate cancer. It gives an insight of possible screening and treatment procedures, and explains the stakes related to prostate cancer.

PROSTATE IN THE SPOTLIGHT

Prostate is a gland that is part of the male reproductive system. Its main role is to produce prostatic fluid for the sperm. The prostate is located below the bladder, above the perineum, and is surrounded by muscles and collagen, forming a capsule. It is traversed by the urethra and ejaculatory duct. Right below the prostate, muscular tissues form the sphincter, controlling the urine flow. The prostate has the shape of a chestnut, with the apex pointing downwards and forwards, and the base pointing upwards and backwards [START_REF] Pietro | Atlas of Imaging Anatomy[END_REF]. Figure 1-1 shows the anatomy of the pelvis and the place of the prostate. In young men, the prostate is around 3cm high, 4cm wide, 2.5cm thick, for a weight of 20g, but get bigger as men get older (the prostate can reach 100g for elderly men) [2]. The prostate is commonly divided into three regions: peripheral (ZP), central (ZC), and transitional (ZT) zones. Figure 1-2 shows the typical layout for an idealized normal prostate. For elderly men, the transitional zone progressively replaces the central zone.

CANCER FIGURES

Prostate cancer is the second most prevalent cancer in men worldwide. In 2020, more than 1 400 000 new cases were counted (Figure 1-3). This represents 14.1% of all new cancers (behind lung cancer with 14.3%). It is the most prevalent in France (25.4%) and in the USA (17.1%). It is the fifth deadly cancer worldwide (in men), and third in France (second in the USA) [5]. The incidence and mortality also varies vary with age [START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF]. Cancer is not uniform within the prostate. The central zone hosts around 5% of cancers, the transitional zone represents 20% of the cancers, and the peripheral zone (accessible during digital rectal examination), has 75% of the cancers.

Another characteristic of prostate cancer is that it appears as multiple foci, of different sizes, shapes, and aggressiveness. Those foci are often stiffer than healthy tissue, making them more detectable during digital rectal examination (DRE). Some authors suggest that small lesions are less dangerous than large ones [START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF].

CANCER DETECTION

When checking for prostate cancer, two first tests are performed:

-Digital Rectal Examination : the practitioner palpates the prostate. A stiff lump or stiff prostate raises suspicion. -Blood test to measure the level of Prostate Specific Antigen (PSA): a high level of PSA is suspicious, but not cancer specific. Note that new markers are being evaluated (e.g. PHI [START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF])

If needed, the doctor may suggest further investigations.

Standard ultrasound imaging (US) is not reliable at detecting prostate cancer [START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF]. Consequently, international guidelines now recommend multiparametric Magnetic Resonance Imaging (mpMRI) [START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF]. It consists in three MRI scans: an anatomic T2-weighted image, a diffusion-weighted image, and a dynamic contrast enhanced (DCE) image. With these new images, the radiologist may spot suspicious regions. Lesions are scored using the Prostate Imaging Reporting and Data System (PIRADS). This PIRADS grading system gives general rules to analyze the 3 images (PI-RADS | American College of Radiology (acr.org)). Based on local anatomical and contrast change, a score can be given from 1 (very low risk) to 5 (very high risk).

Then, after imaging and image grading, biopsies are performed. Small pieces of tissues are extracted with needles. The operation is performed either through the rectal wall (transrectal), or by transperineal approach (however, studies suggest lower risk of infection with transperineal approach [START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF]). The location of the biopsies is set with a sextant template. Other regions can be targeted thanks to the MR images. Biopsies are routinely performed with the guidance of ultrasound. The samples are then analyzed under the microscope.

GRADING OF THE LESIONS

The next step is to assess the aggressiveness of the cancers that are eventually present in the biopsy samples. The key is to discriminate between regions of low risk and those that are potentially life-threatening.

The histopathology analysis will return the ground truth about whether or not the cancer is significant. Indeed, it was shown by Gleason [11] that the architectural pattern of cancer cells reflect the aggressiveness of the tumor. Figure 1-4 shows the different grades, from grade 1 where the cell structure is normal, to grade 5 that presents only cancerous structure. Cancerous tumors are heterogeneous, with multiple foci, some will grow slowly, other will be more aggressive, some will remain confined in the prostate, others will be capable of metastasis [START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF]… To assess the overall risk for the patient, multiple tools have been defined, relying on Gleason's score, PSA level, size and location of the foci… Some of these tools are the d'Amico criteria [START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF], TNM (Tumor, Node, Metastasis) classification [START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF], CAPRA score [START_REF] Cooperberg | Multiinstitutional validation of the UCSF cancer of the prostate risk assessment for prediction of recurrence after radical prostatectomy[END_REF], EAU [START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF]. More information about the different diagnostic tools can be found in a book by Stone [START_REF] Stone | The Prostate Cancer Dilemma[END_REF].

MANAGEMENT OF PATIENTS

Once malignant tumors have been found, different treatment options are available:

-Active surveillance: in this scenario, the patient does not receive any treatment. However, surveillance is increased to follow the evolution of the cancer. This is usually prescribed for low risk tumors. If it becomes dangerous, a treatment will be chosen. The goal is to avoid over-treatment, and keep the same quality of life. Indeed, more patients will die with prostate cancer than of prostate cancer. -Radical prostatectomy: To eradicate the cancer, the prostate is surgically removed. Continence and erectile functions may be affected. A new surgical robot (DaVinci) has made the operation easier [START_REF] Binder | Robot-assisted laparoscopy in urology. Radical prostatectomy and reconstructive retroperitoneal interventions[END_REF].

-Radiotherapy: The radiation can either be external or internal. For external beam radiation therapy, one or several radiation beams are aimed at the cancer so that not the whole body is exposed to radiations. The whole prostate is usually treated. For internal radiation (brachytherapy), a radioactive material is placed inside the cancer so that a higher dose is delivered more locally. -Hormonal therapy: prostate cells rely on the androgen hormone to grow.

Hormonal therapy consists in androgen deprivation. Without the hormone, the cells will stop replicating, and the cancer may shrink or grow more slowly. This technique is however mainly used for highly aggressive or metastatic tumors. [19] -Investigational therapies: new therapeutic methods have emerged and may be considered for the treatment of localized cancer. Those therapies are available in specialized centers and/or in clinical trials.

o Cryotherapy: cryo-needles are inserted in the prostate, which decrease the temperature to freeze the tissues, and destroy the cells. [START_REF] Fahmy | CRYOSURGERY FOR PROSTATE CANCER[END_REF] o High Intensity Focused Ultrasound (HIFU): ultrasound transducers are used to focus ultrasound mechanical energy in a spot. Cells are destroyed by temperature elevation, mechanical forces, and/or cavitation. [START_REF] Madersbacher | High-Energy Shockwaves and Extracorporeal High-Intensity Focused Ultrasound[END_REF] o Other methods include Irreversible Electroporation [START_REF] Onik | Irreversible Electroporation: Implications for Prostate Ablation[END_REF], laser ablation [START_REF] Lindner | Image Guided Photothermal Focal Therapy for Localized Prostate Cancer: Phase I Trial[END_REF],

focal radiation therapy [START_REF] Kovács | Focal radiotherapy as focal therapy of prostate cancer[END_REF] or photodynamic therapy [START_REF] Azzouzi | TOOKAD ® Soluble vascular-targeted photodynamic (VTP) therapy: determination of optimal treatment conditions and assessment of effects in patients with localised prostate cancer[END_REF].

NEED FOR A BETTER IMAGING

The previous paragraph shows that there are treatment options to prostate cancer. However, they usually involve the whole prostate and may cause collateral damage. If touched, surrounding nerves, muscles, and urethra, may be damaged causing possible urinary and/or sexual issues. It is thus important to treat only when necessary, and where necessary. The latter is called focal therapy.

Accurate detection and accurate assessment of cancer aggressiveness are the key to a better management. In the early days, 90% of patients diagnosed with low-risk prostate cancer by PSA and biopsy in the US were treated with definitive therapy [START_REF] Stone | The Prostate Cancer Dilemma[END_REF]. Since then, there has been a controversy about the use of PSA level. That is why more complex staging rules are defined, and the tendency is to combine multiple markers. Biopsies reveal the ground truth, but they provide a discrete vision of the prostate. If under-sampled, cancers might have been missed. If over-sampled, small low-risk cancers might be considered as aggressive tumors. With the wide availability of ultrasound scanners, biopsies are usually performed with ultrasound. However, with the recent recommendation to use mpMRI, MR images can be fused with the US image for special targeting of suspicious lesions. mpMRI has shown great hope at detecting cancers but faces some limits. The first one is specificity. Indeed, even though mpMRI presents a high sensitivity, its specificity remains low. For example, Drost showed that for the detection of grade 2 or higher prostate cancer, MRI had a sensitivity of 91% and a specificity of 0.37% [START_REF] Drost | Prostate MRI, with or without MRI-targeted biopsy, and systematic biopsy for detecting prostate cancer[END_REF]. The other limit is related to reproducibility. Despite the use of the PI-RADS, MR hardware and sequence parametrization can vary, yielding different images, and, qualitative analysis of the signal intensity in the images can vary between operators. It is thus quite troublesome to correctly analyze the images, hence a low reproducibility [START_REF] Dickinson | Magnetic Resonance Imaging for the Detection, Localisation, and Characterisation of Prostate Cancer: Recommendations from a European Consensus Meeting[END_REF], [START_REF] Westphalen | Variability of the Positive Predictive Value of PI-RADS for Prostate MRI across 26 Centers: Experience of the Society of Abdominal Radiology Prostate Cancer Disease-focused Panel[END_REF]. New tools like Computer Aided Diagnostic (CAD) are being developed to help the radiologists [START_REF] Gaur | Can computer-aided diagnosis assist in the identification of prostate cancer on prostate MRI? a multi-center, multi-reader investigation[END_REF]- [START_REF] Saha | End-to-end prostate cancer detection in bpMRI via 3D CNNs: Effects of attention mechanisms, clinical priori and decoupled false positive reduction[END_REF].

THE PLACE FOR MAGNETIC RESONANCE ELASTOGRAPHY

Elastography is a technique to assess tissues stiffness. During DRE the doctor palpates the prostate and feels if it is rather stiff or soft. However, he can only access the superficial part of the prostate. Elastography could provide a 3D map of the stiffness of the whole prostate. The measure is quantitative so should be less prone to variation between centers and operators. It provides additional information to the radiologist and could help to discriminate between healthy tissue, low-risk lesions, and high-risk lesions. Given that MR scans are now performed pre-biopsy, an MR elastography sequence could be added to the track list. We hypothesize that Magnetic Resonance Elastography (MRE) could help image and characterize prostate cancers.

CHAPTER

2 Elastography 2

BASIC PRINCIPLE

Since the dawn of time, manual palpation was the first (if not only) means of diagnostic for accessible organs. When palpating, one is actually applying a certain pressure on the tissue and trying to assess its mechanical response (does it deform easily or not?). Tissue mechanical response to external stimulus (palpation, induced mechanical force…) represents the field of biomechanics and was first formalized by Hooke in 1678 [START_REF] Hooke | Lectures de potentia restitutiva, or, Of spring explaining the power of springing bodies : to which are added some collections viz[END_REF]. It was only in the twentieth century, in the late 80's -early 90's, that the application for medical imaging emerged, and it was shown that tissue biomechanics could reflect pathological change [START_REF] Fung | Biomechanics[END_REF], [START_REF] Sarvazyan | Biophysical Bases of Elasticity Imaging[END_REF]. Elastography is this method of analyzing tissue deformations and trying to infer stiffness parameters. To do so, the first step is to induce a force to deform the tissues of interest. This can be done using a compression device, a vibrator, focused ultrasound, or even natural physiological source (e.g. arterial pulsation). Then, the tissue deformation is measured using a suitable method (ultrasound, MRI, optics…). Finally, mechanical parameters are assessed.

MECHANICAL BACKGROUND

A material is considered linear elastic if after undergoing a deformation, it can return to its original state, and its unloading path is the same as its loading path. Elastic moduli are related to the stress applied to the volume and to the produced strain. A material can deform either with a change of volume, or at constant volume. These are referred as compression and shear deformations (Figure 2-1), and are associated with different moduli:

-Shear deformation is associated to the Shear modulus μ.

-Isotropic pressure on all faces is associated to the Bulk modulus K -Uniaxial, normal stress is associated to the Young's modulus E. -Poisson's coefficient characterizes the relation between normal and shear strain.

It gives an insight of the compressibility of the medium. For an incompressible medium, = 0.5. -Lamé parameters are another couple of material dependent moduli: (first Lamé parameter), and μ (second Lamé parameter, also Shear modulus).

Tristan DERUELLE 31 All of these moduli are interconnected, and the knowledge of two gives all the others. The standard static mechanical method to measure a sample elasticity assesses the Young's modulus. The method consists in applying a certain uniaxial stress to a sample, as shown in Figure 2-2, and measuring the resulting strain = ( )/ . The subsequent Young's modulus is computed with Hooke's law: = / . In dynamic regime, several assumptions are made to simplify the equations. A linear isotropic medium is considered. Under the assumption of small deformation, the wave equation can be written as:

= ( + 2μ) ( ) μ ( ) (2-1)
Where represents the particle displacement, (1000kg/m 3 ) the Lamé parameters. By taking the notation for the nabla operator, as the divergence, and as the curl operator, equation (2-1) can be rewritten as:

= ( + 2μ) μ (2-2)
Using the Helmholtz theorem, the solution can be decomposed as the sum of an irrotational (curl-free) and a solenoidal (divergence-free) fields:

= + (2-3) = = (2-4) = 0 = 0 (2-5)
With a scalar potential, and a vector potential.

Given these properties, equation (2-2) can be rewritten as:

+ = ( + 2μ) μ (2-6)
This equation shows in fact that two waves are propagating: a compressional (curl-free) wave and a shear (divergence-free) wave. These two waves have very different behaviors.

In a compressional wave, the particles move parallel to the direction of wave propagation. In a shear wave, the particles move perpendicular to the direction of propagation (see Figure 2-4). These two waves travel at different speeds, given by:

= = + 2 (2-7) = = (2-8)
With kP and kS the velocity is directly related to the square root of the compression or shear moduli. Elastic moduli can thus be derived from the measure of the wave velocity. The shear modulus being several orders of magnitude smaller than the bulk modulus (Figure 2-3), the shear wave velocity is much smaller than the compression wave velocity. In practice, shear wave velocity is in the order of 1-10 m/s while the compression wave travels at the speed of sound (~1500 m/s). It is thus easier to image shear waves.

However, the materials are usually dissipative, and damping must be taken into account.

In that case, the elastic moduli become functions of the frequency and can be written:

( ) = + ( ) (2-9) μ ( ) = μ + μ ( ) (2-10)
1 and μ1 the time independe 2 and μ2 the relaxation terms associated to viscosity. μ1 and μ2 are sometimes referred to as storage modulus and loss modulus.

The wave numbers become complex:

= + = + 2μ (2-11) = μ (2-12)
The general formulation for a propagating harmonic plane wave is:

( , ) = ( ) = ( ) (2-13)
As the wave must decay, k2 must be negative: = with 0 the attenuation coefficient. k1 is related to the phase velocity: = .

Hence, taking the square of equation (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF] yields to:

= 2 (2-14) = (μ + μ ) μ + μ (2-15)
The real and imaginary parts of equations (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF] and (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF] are equal. Solving the system gives the expressions of the phase velocity and damping:

= μ 1 + μ μ (2-16) = 2μ 1 + μ μ 1 + μ μ (2-17)
From there, one common simplification is to assume μ2 μ1. Equations (2-16) and (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF][START_REF] Cooperberg | Multiinstitutional validation of the UCSF cancer of the prostate risk assessment for prediction of recurrence after radical prostatectomy[END_REF][START_REF] Stone | The Prostate Cancer Dilemma[END_REF] simplify to:

= μ (2-18) = μ 2μ (2-19)
Finally, the Kelvin-Voigt model is generally used, for which

= = μ (2-20) = 2μ μ (2-21)
Under those assumptions, the shear wave speed is roughly the same as the non-viscous speed; and the attenuation is proportional to the square of the frequency.

ULTRASOUND ELASTOGRAPHY

The first investigations in elastography were performed using ultrasound [START_REF] Lerner | Sono-Elasticity: Medical Elasticity Images Derived from Ultrasound Signals in Mechanically Vibrated Targets[END_REF]- [START_REF] Dickinson | Measurement of soft tissue motion using correlation between A-scans[END_REF] to track tissues' displacements. One of the first elastography method, sonoelastography was proposed by Lerner and Parker in 1987 [START_REF] Lerner | Sono-Elasticity: Medical Elasticity Images Derived from Ultrasound Signals in Mechanically Vibrated Targets[END_REF]. It consisted in applying an external vibration, and imaging tissue displacements using Doppler techniques. However, the measure is not quantitative, it produces a relative image of local stiffness.

Later, in 1991, Ophir presented strain elastography (or compression elastography) [START_REF] Ophir | Elastography: A Quantitative Method for Imaging the Elasticity of Biological Tissues[END_REF]. Unlike sonoelastography, strain imaging is performed in quasi-static. Displacements are estimated by comparing B-mode images before and after applying a small compression.

The advantage is that compression can easily be applied with the ultrasound probe held during the scan. This measure of the strain is related to the elasticity by Hook's law = / . For a given stress, a smaller strain ( ) is associated to a higher stiffness (higher Young's modulus E). In practice, the stress that is actually applied is usually not known. Consequently, the produced images are relative. Besides, in order to compare two areas in an image, they must be subject to the same stress. Having a uniform stress across the image is the main difficulty in strain elastography.

Later work focused on transient excitation such as using a single vibration cycle. By doing so, the compression wave is separated from the shear wave because of the large difference in velocity. The propagation of the shear wave solely can thus be tracked. Such method was investigated by Catheline et al. in the late 90s' [START_REF] Catheline | A solution to diffraction biases in sonoelasticity: The acoustic impulse technique[END_REF]. They used a low frequency mechanical excitation and modulated ultrasound to track tissue deformations. This led to the commercial product FibroScan (Echosens, Paris, France) in 2003, routinely used for staging hepatic fibrosis [START_REF] Sandrin | Transient elastography: a new noninvasive method for assessment of hepatic fibrosis[END_REF].

At the same time, the use of ultrasound acoustic radiation force (ARF) was investigated [START_REF] Torr | The acoustic radiation force[END_REF]. Ultrasounds can be focused to generate a remote "push" in the tissues.

This technique for elastography was presented by Sarvazyan in 1998 [START_REF] Sarvazyan | Shear wave elasticity imaging: a new ultrasonic technology of medical diagnostics[END_REF] and developed as an imaging modality by Nightingale [START_REF] Nightingale | On the feasibility of remote palpation using acoustic radiation force[END_REF]. Elasticity maps can be estimated by measuring the induced displacement or tissue relaxation. The development of ultrafast imaging led to "Supersonic" shear imaging [START_REF] Bercoff | Supersonic shear imaging: a new technique for soft tissue elasticity mapping[END_REF] in 2004. Several ARF pushes are rapidly fired along a line, generating a cone of shear waves. The propagation of the shear waves is tracked with the ultrafast scanner capable of 10 000 frames/sec. The velocity of the shear wave can thus be measured. The advantage of using ARF pushes is that the sole ultrasound transducer can create the vibration and measure the displacement. However, the shear waves of small amplitude are rapidly attenuated, resulting in a small region available for the measurements. Also, the intensity of the ultrasound is limited to avoid mechanical and thermal effects in the tissues. Finally, due to attenuation, it may be difficult to image deep organs.

PRINCIPLES OF MAGNETIC RECONANCE ELASTOGRAPHY

Unlike ultrasound, MRI has the ability to image any part of the body, non-invasively, and in 3D. The method for capturing propagating mechanical waves with MRI was first shown by Muthupillai et al in 1995 [START_REF] Muthupillai | Magnetic resonance elastography by direct visualization of propagating acoustic strain waves[END_REF]. The principle is similar to that of ultrasound-based methods, and can be decomposed into the following steps (Figure 2345):

1. Generate waves in the body (typically around 25-300 Hz).

2. Take snapshots of the propagating waves using MRI.

3. Perform (or not) some pre-processing to the images 4. Estimate the wave velocity in each point of the image thanks to various algorithms MRI does not have the ability to take many images per seconds. The sampling frequency is thus much lower than the excitation frequency. The solution is to synchronize the vibration with the acquisition to end up in a controlled stroboscopic state. It is crucial to have a perfect synchronization between the excitation and the MRI scanner. In particular, relative phase between each snapshot must be known.

WAVE GENERATION

The first step of MRE is to generate waves in the body. Although it could be possible to generate waves with acoustic radiation force [START_REF] Souchon | Transient MR elastography (t-MRE) using ultrasound radiation force: Theory, safety, and initial experiments in vitro[END_REF], continuous mechanical vibrators are usually used. Indeed, the organ may not be accessible with ultrasound, and it may be troublesome to correctly position the ultrasound transducer and the patient inside the MR scanner. The ultrasound probe should also be perfectly MR compatible and immune to the electromagnetic fields and impulses. Finally, due to the low penetration of waves and slow framerate of the MRI scanner, ultrasound pushes may need to be repeated a certain number of times, increasing the duration of the acquisition. Mechanical and thermal effects may also require special attention.

Mechanical vibrator can induce larger displacement amplitude than ultrasound. However, depending on the accessibility of the organ of interest, it can still be challenging to transmit waves of sufficient amplitude. For example, the liver is quite accessible and superficial so it is easy to vibrate. General Electric's (GE) commercial hepatic device (GE Resoundant) simply consists in an acoustic drum-like driver. On the other hand, the brain is (fortunately) protected by the skull so the waves are more likely to be guided by the bone and not transmitted to the brain.

The vibrator should be adapted to the organ of interest and account for the surrounding anatomy. In particular, waves can be guided away from the organ of interest, or can be damped before reaching it. In the latter case, vibration frequency may also play a role (following equation (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF][START_REF] Cooperberg | Multiinstitutional validation of the UCSF cancer of the prostate risk assessment for prediction of recurrence after radical prostatectomy[END_REF][START_REF] Stone | The Prostate Cancer Dilemma[END_REF][START_REF] Binder | Robot-assisted laparoscopy in urology. Radical prostatectomy and reconstructive retroperitoneal interventions[END_REF][START_REF] Pagliarulo | Contemporary Role of Androgen Deprivation Therapy for Prostate Cancer[END_REF][START_REF] Fahmy | CRYOSURGERY FOR PROSTATE CANCER[END_REF][START_REF] Madersbacher | High-Energy Shockwaves and Extracorporeal High-Intensity Focused Ultrasound[END_REF]).

Note that some techniques rely on the natural waves present in the body (arterial pulsation, breathing…). Wave velocity can be estimated using noise correlation techniques [START_REF] Catheline | Time Reversal of Elastic Waves in Soft Solids[END_REF]. This technique can be adapted to low frame-rate scanners such as MRI [START_REF] Catheline | Tomography from diffuse waves: Passive shear wave imaging using low frame rate scanners[END_REF], [START_REF] Zorgani | Brain palpation from physiological vibrations using MRI[END_REF], but requires the acquisition of a consequent number of images (~100 times conventional MRE). The resulting acquisition time is not realistic for clinical use, except maybe in specific applications where using external vibrators is not recommended (e.g. cerebral trauma).

MR WAVE IMAGING

For a prior understanding of the MRI acquisition, the reader is referred to several books: [START_REF] Liang | [END_REF]- [START_REF] Hirsch | Magnetic Resonance Elastography[END_REF].

Once the vibrations are in place, the waves are imaged with MRI. The displacement is encoded in the phase of the MR image using a motion-sensitive sequence. During image acquisition, additional magnetic Motion Encoding Gradients (MEG) are used to sense the displacement. However, the phase will also be sensible to magnetic field inhomogeneities and susceptibility effects. The captured phase will thus be the superposition of a static background (from these undesirable effects), and the wave information at a particular time. One solution to retrieve only the wave contribution, is to acquire two images of the wave but with opposite polarity. This can either be done by capturing two images with a re orientation (e.g. +x and -x orientations). The latter enables to increase motion sensitivity and get the effective wave displacement:

= + (2-22) = (2-23) = 1 2 ( ) (2-24)
More specifically, the wave phase is encoded thanks to the synchronization between a particle oscillation and the MEG. Under vibrations, a spin will oscillate around its equilibrium position. With the application of a MEG the spin will accumulate phase while moving:

( ) = ( ) + ( ) ( ) (2-25)
With G the gradient field, and x the trajectory of the spin.

Figure 2-6 illustrates the phenomenon. Let's take the example of a particle without original phase, and with a sinusoidal displacement. Without MEG, its phase at the end of a cycle will be the same as at the beginning. With the application of a MEG, the more the particle moves away from its equilibrium position over time, the more phase it accumulates. It is as if the particle had the original phase of a neighbor. Then, the accumulated phase over a cycle can be written as:

( ) = ( ) + ( )sin ( + ) (2-27)
The accumulated phase thus depends on the original phase offset. Particles which are at different moment in their oscillatory cycle will have different phase. For example, if the particle and the MEG are in phase quadrature, there will be no accumulated phase. That is how a phase wave image is encoded. The design of the MEG function G(t) can also play a role for different properties (gradient moment nulling). More details are given in [START_REF] Hirsch | Magnetic Resonance Elastography[END_REF].

For a full MRE acquisition, such images are acquired multiple times: acquisition is repeated for the different slices, for the different motion encoding directions (usually for x-y-z, the three directions of space, and eventually each one with positive and negative orientations to get a stroboscopic movie of the wave propagation.

PRE-PROCESSING 2.7.1 Unwrapping

The resulted phase wave images are mapped to [-s usually present phase jumps. The first processing to perform is thus unwrapping. Phase unwrapping is not limited to MRE but affects domains on a wide range. MRE has the particularity to have 3D data, so it is preferable to perform the unwrap in 3D. Indeed, waves travel in the three directions of space, so if a wave is coming through the imaging plane, it must be unwrapped correctly and ensure the 3D continuity. Common algorithms for 2D images were presented in a book [START_REF] Ghiglia | Two-Dimensional Phase Unwrapping: Theory, Algorithms, and Software[END_REF], in particular Flynn's minimum discontinuity algorithm, Goldstein branch cut algorithm. Faster 2D and 3D algorithms were developed based on pixels reliability and non-continuous path [START_REF] Herráez | Fast two-dimensional phase-unwrapping algorithm based on sorting by reliability following a noncontinuous path[END_REF]- [START_REF] Marcsous | unwrap[END_REF]. The latter are the ones used in the python scikit-image package. Another interesting algorithm make the use of the Laplacian operator [START_REF] Schofield | Fast phase unwrapping algorithm for interferometric applications[END_REF]- [START_REF] Loecher | Phase unwrapping in 4D MR flow with a 4D single-step laplacian algorithm[END_REF]. We will describe the last two algorithms here.

Laplacian based unwrap

It is important to note that for the wrapped phase and true phase :

= + 2 (2-28)
with k an integer (positive or negative).

So:

= (2-29) = (2-30)
with N the number of dimensions.

Hence, the Laplacian of the phase can be written as:

= (2-31)
Or:

= [cos( ) sin( ) sin( ) cos( )]
(2-32)

The forward and inverse Laplacian can be implemented in various ways, using a finite difference scheme, or more computationally efficient using the Fast Fourier Transform, or other transforms (Discrete Cosine Transform, Discrete Sine Transform…) depending on the boundary conditions.

It should be kept in mind that the Laplacian is a second order differential operator so when applying the inverse, two integration constants cannot be recovered. Those constants are offset and linear trends. In MRE, those may arise from magnetic field inhomogeneities so it might be better to suppress them. This method is extremely fast, but may induce some bias in specific situations.

Quality guided unwrap

Unlike the Laplacian-based algorithm, quality guided methods are based on the measure of the reliability of the pixels. Reliability assesses the possibility of the pixel being wrapped added or subtracted to unwrap the pixel. In the method proposed by Herráez [START_REF] Herráez | Fast two-dimensional phase-unwrapping algorithm based on sorting by reliability following a noncontinuous path[END_REF] and extended in 3D by Abdul-Rahman [START_REF] Abdul-Rahman | Fast three-dimensional phase-unwrapping algorithm based on sorting by reliability following a noncontinuous path[END_REF], the reliability function is also defined on a second order difference. This provides a measurement of the degree of convexity/concavity of the phase which yields to a better detection of inconsistencies. The second order difference is similar to computing a centered difference Laplacian on 9 neighbors in 2D (orthogonal and diagonal neighbors) or 27 in 3D. However, for the computation of this Laplacian, each pair of pixels are unwrapped independently. For example, in 2D, let's consider the pixels of an image as in Figure 234567. The second difference operator D is computed as:

D = 1 + 2 + 1 + 2 (2-33)
And for example:

1 = [ ( 1, ) ( , )] [ ( , ) ( + 1, )] (2-34)
where is a simple unwrapping function. For example, ( ) = ( , 2 ) .

Finally, the reliability R of a pixel is defined as the inverse of D: R=1/D.

The reliability of each pixel is computed. Then, each edge is assigned a reliability. An edge is the connection between two adjacent pixels (vertically, horizontally, or eventually in depth). The reliability of an edge is defined as the sum of the reliability of the two adjacent pixels. All edges are sorted by reliability. Unwrapping is performed from the most to the least reliable. The path is thus not continuous, but based on the edge reliability. At the beginning, the two pixels with the most reliable edge are unwrapped and form a group. Then, for the following edges, if one pixel belongs to a group, and the other not, then, that pixel is unwrapped and added to the group. If both pixels do not belong to a group, they are assigned to a new group. If both pixels belong to two distinct groups, the smallest group is unwrapped and merged into the largest one. Once all edges are processed, the borders of the image are unwrapped.

This yields to a relatively fast and robust algorithm. In some cases, the quality guided algorithm, may also face some issues. In the example shown Figure 2-9, the quality guided algorithm could not unwrap properly one region. The Laplacian based method successfully unwrapped the data, but the full amplitude was not recovered. Quality guided algorithm can recover the exact original field. If data is too noisy or if there is no need to recover the full amplitude, the Laplacian algorithm may be considered. Also, successive unwrapping processes could be performed.

Example of unwrapping

Frequency filter

The wavefield should then be filtered at the desired frequency. Indeed, other frequencies may be present, either naturally from the body (e.g. arterial pulse), or because the source used introduces other harmonics (desired or not). This is done by a Fast Fourier Transform and selection of the harmonic of interest. This operation also removes the phase background.

The issue of compression

MRI senses and images all kind of displacements. It does not differentiate the shear displacement of interest from the compression displacement, but rather detect both contributions. The presence of compression can thus bias the estimation of the shear modulus. Three different methods are commonly used to handle this. The first one is to assume that the wavelength is too large, and the associated displacements can be neglected. This assumption can lead to overestimating the shear modulus, especially close to the source of vibration. The second option is to use a high-pass filter to remove long wavelength contributions. This is usually coupled to a low-pass filter to remove noise. However, with this technique, it is to up the user to choose the type of filter and its characteristics (e.g. cut-off frequency). The output is thus user-dependent. The last option, the most physical way to date, is to apply a curl operator [START_REF] Sinkus | Viscoelastic shear properties of in vivo breast lesions measured by MR elastography[END_REF]. As explained in paragraph 2.2, the compressional part of the wave equation can be removed by applying a curl operator. However, this operation involves the computation of derivatives and thus decreases the Signal to Noise Ratio (SNR). To use this technique, the original SNR should be high enough, and/or advanced method for the computation of noisy derivatives can be used (e.g. using Savitzky-Golay filters [START_REF] Savitzky | Smoothing and Differentiation of Data by Simplified Least Squares Procedures[END_REF], [START_REF] Schafer | What Is a Savitzky-Golay Filter?[END_REF]).

RECONSTRUCTION

There have been many different investigations to reconstruct tissue properties: direct inversion, finite-element-method based inversions, multifrequency dual elasto-visco inversion, just to name a few. One of the most widely used is called Local Frequency Estimation (LFE) and was first proposed by Knutsson et al [START_REF] Knutsson | Local multiscale frequency and bandwidth estimation[END_REF] for image processing. It was later brought to MRE by Manduca et al [START_REF] Manduca | Visualization of tissue elasticity by magnetic resonance elastography[END_REF], [START_REF] Manduca | Image analysis for magnetic resonance elastography[END_REF]. The reconstructions performed in this work will use this technique, and the present paragraph describes its operation. It should be kept in mind that for an accurate estimation, at least one half of the wavelength should be present in the Region of Interest (ROI) [START_REF] Manduca | Magnetic resonance elastography: Non-invasive mapping of tissue elasticity[END_REF].

LFE computes the local spatial frequency in each voxel. Knowing the spatial and temporal frequencies, the wave velocity can then be computed.

The principle of this algorithm is the following:

Given an analytic monochromatic signal s, the local frequency is defined as the derivative of its phase:

= 1 2 arg( ) (2-35)
The local frequency can be estimated as the ratio of the output of two filters. These filters are usually log-normal wavelets and can be written as:

( ) = (2-36)
where is the central frequency of the filter, and = ( ) with B the bandwidth at -6dB.

Taking the ratio of two filters of central frequencies , yields:

( ) ( ) = (2-37)
A particular case occurs when 2 ln = 1.

In that case:

= (2-38)
And equation (2-37) becomes:

( ) = ( ) (2-39)
Multiplying equation (2-39) by the Fourier transform of the signal gives:

( ) ( ) = ( ) ( ) (2-40)
Noticing that ( ) is similar to taking the derivative of in the Fourier domain, by taking the inverse Fourier transform, it can be shown that the local frequency can be written as:

= (2-41)
This works well for monochromatic signals. When multiple frequencies are present, the local frequency can be estimated as a weighted average of equation (2-41) over multiple pairs of filters (equation ). The successive pairs of filters can be determined using equation .

= (2-42)
where are the weighting factors. To apply more weight to useful frequencies, the weighting factors can be defined as the output of the current filter so that = . In that case, equation simplifies to:

= (2-43)
The convolution can be easily computed in the Fourier domain:

= ( ) ( ) (2-44)
The method can be extended from 1D to nD by adding a spatial filter to the original lognormal filter. The new expression can be written:

( ) = ( ) ( ) (2-45)
with the spatial component:

( ) = ( ) , > 0 0, (2-46)
and is the filter directing vector.

Equation is extended by summing also over the k directions:

= , , , , (2-47)
Usually, orthogonal directions and with opposite polarities are chosen (4 directions in 2D, 6 in 3D). Some studies showed that better results could be obtained by filtering using more directions and filtering in time to reduce reflections influence [START_REF] Manduca | Spatio-temporal directional filtering for improved inversion of MR elastography images[END_REF] (an extension to 3D and with more direction is proposed in Appendix A). The basic workflow of LFE algorithm is summarized in Figure 2

CONCLUSION

Elastography technique is both old and new, using palpation and modern imaging capabilities. This methodology stirred great enthusiasm for pathologies assessment. Ultrasound based techniques brought significant contributions to elastography and is even now being used routinely for staging hepatic fibrosis. Although MRE principle differs little from ultrasound elastography, MRE has the ability to image a 3D volume, track displacements in any space direction, and image any part of the body non-invasively. It still faces similar issues to ultrasound, especially in tissue excitation.

THE CHALLENGES OF PROSTATE ELASTOGRAPHY

Since the early days of elastography, interest has been brought to the prostate with the main question: can we make prostate elastography clinically relevant?

Indeed, as explained in the introduction, one of the first screening tests for prostate cancer is palpation of the superficial accessible region. Virtual palpation of the entire prostate could be a valuable asset to radiologists.

Prostate stiffness has been assessed ex vivo using standard mechanical techniques or elastography. However, in vivo assessment is more challenging. Indeed, the prostate is relatively deep, stiff, and surrounded by damping tissues. It is thus difficult to excite and image. To face this problem, three solutions are possible. The first one is to find an anatomical window where the source of excitation can be placed closer to the prostate, or find a region where the propagation is less damped or more guided towards the prostate. The second solution is to increase the amplitude of the source. The third solution is to decrease the frequency of the vibrations as low frequencies are less damped (following equation (2)(3)(4)(5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF][START_REF] Cooperberg | Multiinstitutional validation of the UCSF cancer of the prostate risk assessment for prediction of recurrence after radical prostatectomy[END_REF][START_REF] Stone | The Prostate Cancer Dilemma[END_REF][START_REF] Binder | Robot-assisted laparoscopy in urology. Radical prostatectomy and reconstructive retroperitoneal interventions[END_REF][START_REF] Pagliarulo | Contemporary Role of Androgen Deprivation Therapy for Prostate Cancer[END_REF][START_REF] Fahmy | CRYOSURGERY FOR PROSTATE CANCER[END_REF][START_REF] Madersbacher | High-Energy Shockwaves and Extracorporeal High-Intensity Focused Ultrasound[END_REF]).

The difficulty to perform prostate elastography lies in the tradeoff between patient's acceptance, sufficient excitation amplitude, and high enough frequency.

PROSTATE ELASTOGRAPHY

First work on using elastography for prostate was done by Lee et al. in 1991. They used sonoelastography on in vitro human prostate [START_REF] Lee | Sonoelasticity imaging: results in in vitro tissue specimens[END_REF]. They reported that normal prostates present uniform vibration, while stiff inclusions present a lack of vibration. This preliminary work demonstrates the interest of prostate elastography. It was pursued by Parker et al. in 1993 [70]. With precompression between 1-10%, they measured a linear Young's modulus between 2-4 kPa in healthy in vitro prostates.

Later, in 1998, with the advent of strain elastography, Krouskop et al [START_REF] Krouskop | Elastic Moduli of Breast and Prostate Tissues under Compression[END_REF] reported quantitative measurements, in in vitro human prostate with cancers. They used strain elastography with a mechanical indenter to reach 2 and 4% of precompression with frequencies between 0.1 and 4 Hz. They reported Young's modulus between 55 and 71 kPa in the normal prostate, and between 96 and 241 kPa in the cancer. Prostate cancer was significantly stiffer than normal tissue.

There have then been numerous studies to assess prostate stiffness ex vivo. However, the mechanical properties strongly depend on the method, the frequency, but also tissue condition. In particular, in vivo and in vitro specimens may not present the same elasticity.

McGrath showed that tissue fixation alters the stiffness of the specimens [START_REF] Mcgrath | Quasi-static magnetic resonance elastography at 7 T to measure the effect of pathology before and after fixation on tissue biomechanical properties[END_REF].

Things are different in vivo. The physicist is looking to quantitatively measure the prostate stiffness. As for the medical doctor, he is looking to detect the cancers, so cancers should show a different contrast. Some studies were performed in that way using strain elastography with an endorectal probe. In 2005, König et al. performed B-mode imaging and strain elastography to a cohort of 404 men [START_REF] König | INITIAL EXPERIENCES WITH REAL-TIME ELASTOGRAPHY GUIDED BIOPSIES OF THE PROSTATE[END_REF]. Sextant biopsies revealed 151 cases of cancer (37.4%). Among the 151 cases, 127 were detected using elastography (84.1%). A similar study was performed by Tsutsumi et al in 2007 on 51 patients with prostate cancer [START_REF] Tsutsumi | The impact of real-time tissue elasticity imaging (elastography) on the detection of prostate cancer: clinicopathological analysis[END_REF]. Elastography could detect 94% of anterior tumors, 76% of middle tumor, and 57% of posterior tumors. These studies show that strain elastography can detect prostate cancer. However, it is limited by the penetration depth of ultrasound, and by the lack of uniform compression over the gland. This can yield intra and inter operator differences, which is exacerbated as the measurement is qualitative [START_REF] Correas | Ultrasound elastography of the prostate: State of the art[END_REF].

With the commercial development of the Aixplorer system (Aixplorer system, SuperSonic Imagine, Aix en Provence, France), shear wave elastography took over. With this method, real-time quantitative measurements can be performed. One major study was performed by Correas et al, in 2015, on 184 patients with biopsy reference in the peripheral zone [START_REF] Correas | Prostate Cancer: Diagnostic Performance of Real-time Shear-Wave Elastography[END_REF]. An example of ultrasound-based images is presented on Figure 3-1. He reported a median Young's modulus of 53 kPa for prostate cancer, and 21 kPa for benign regions.

Stiffness increased with Gleason's score. Sensitivity, specificity, Positive and Negative Predicted Values, and area under the receiver operating characteristic curve were respectively 96%, 85%, 48%, 99%, and 95%. This demonstrates the use of elastography for prostate cancer. Some limitations are the small reconstructed area, image stabilization and penetration issues [START_REF] Correas | Ultrasound elastography of the prostate: State of the art[END_REF]. Clinical use of ultrasound elastography for the prostate is significant. To overcome anatomical barrier and image the prostate, endorectal probes are used. The work done with ultrasound shows that elastography can bring valuable piece of information for cancer diagnosis. The World Federation for Ultrasound in Medicine and Biology (WFUMB) even provides guidelines and recommendations for a better understanding and practice [START_REF] Barr | WFUMB Guidelines and Recommendations on the Clinical Use of Ultrasound Elastography: Part 5. Prostate[END_REF].

Ding et al [START_REF] Ding | Development and validation of a nomogram based on multiparametric magnetic resonance imaging and elastography-derived data for the stratification of patients with prostate cancer[END_REF] and Wei et al [START_REF] Wei | Prostate Cancer Gleason Score From Biopsy to Radical Surgery: Can Ultrasound Shear Wave Elastography and Multiparametric Magnetic Resonance Imaging Narrow the Gap?[END_REF] used ultrasound elastography in combination with mpMRI. They showed that this could improve cancer detection and stratification.

On the other hand, MRI is not limited by penetration depth, image acquisition is not limited to a 2D plane, and is not invasive. However, most of the work relies on being able to vibrate the prostate. Common MRE system (e.g. GE Resoundant) were designed for the liver and are therefore not suited for targeting the prostate. It cannot produce sufficient displacement.

To generate waves close to the prostate, several attempts of transurethral probes were designed ( [START_REF] Chopra | In vivo MR elastography of the prostate gland using a transurethral actuator[END_REF], [START_REF] Arani | Transurethral prostate magnetic resonance elastography: Prospective imaging requirements[END_REF]), and tested on phantom and canine subjects. The proximity to the prostate allowed them to reach higher frequencies (up to 600 Hz) for a better resolution in the reconstructions. The clinical acceptance of such an approach is questionable.

Similarly to ultrasound procedure, endorectal setup proved to be feasible by Arani in 2013 [START_REF] Arani | Incorporating endorectal MR elastography into multi-parametric MRI for prostate cancer imaging: Initial feasibility in volunteers[END_REF], [START_REF] Arani | The feasibility of endorectal MR elastography for prostate cancer localization[END_REF]. They used piezo-ceramic actuators coupled to an endorectal MR coil to generate the waves. It allowed them to get close to the prostate and have a strong MR signal. This was tested on 6 healthy volunteers. They reported waves of amplitude from 1 to 50 μm, and measured them in the prostate up to 300 Hz. Depending on the slices, the peripheral zone was either stiffer of softer than the central part. The repeatability was not as high as they expected, mainly due to 2D processing and wave attenuation.

One of the first non-invasive work was done by Kemper et al in 2004 on 7 healthy volunteers [START_REF] Kemper | RöFo -Fortschritte auf dem Gebiet der Röntgenstrahlen und der Bildgeb[END_REF], using an external oscillator placed on the pubic bone. With a frequency of 65 Hz, and an amplitude between 200-400 μm at the source, they could measure a mean displacement of 53±31 μm in the prostate. They reported that the shear modulus in the peripheral region (3.3±0.5 kPa) appears higher than in the central region (2.2±0.3 kPa). This work demonstrated the feasibility of in vivo MRE of the prostate. However, they raised the need to have a strong contact between pelvis and transducer, as well as the extended acquisition time (30mins) which requires patient's full cooperation.

In 2011, Li et al pursued with the pubic excitation [START_REF] Li | A feasibility study of MR elastography in the diagnosis of prostate cancer at 3.0T[END_REF]. The study involved 10 healthy volunteers and 8 patients with prostate cancer. They used an excitation frequency of 100 Hz, but no data was given regarding the generated displacement. They reported shear modulus in prostate cancer and healthy peripheral zone. Cancerous tissues (6.55±0.47 kPa) were stiffer than healthy tissue in the peripheral zone (2.26±0.45 kPa). They also noticed a good correlation between the Gleason score and elasticity in the cancer. However, they suffered from a poor resolution (3.91 mm isotropic voxels) yielding a few voxels per region.

However, Sahebjavaher explains that not only the pelvis is difficult to move with large amplitude without significant discomfort to the patient, but this method may also not be repeatable across subjects of different sizes [START_REF] Sahebjavaher | Transperineal prostate MR elastography: Initial in vivo results[END_REF]. Consequently, they proposed a transperineal approach, where the device comes in contact with the patient's undergarment. The first prototype was presented in 2013 as a hydraulic vibrating system. This was tested on 6 healthy volunteers with a frequency of 45 Hz. They managed to measure vibrations up to 30 μm in the prostate, and a global Young's modulus of 8.6±3.1 kPa (equivalent to a shear modulus of 2.9±1.0 kPa). Compression waves were removed using the curl operator. Measurements were repeatable and comparable with Kemper. In 2014, they upgraded the device to an electromagnetic transducer [START_REF] Sahebjavaher | Prostate MR elastography with transperineal electromagnetic actuation and a fast fractionally encoded steady-state gradient echo sequence[END_REF]. Care was put into making a proper shielding to avoid electromagnetic disturbances, and image artifacts. With that transducer, vibrations of 400 μm can be reached at 70 Hz. Displacement amplitude decreases with frequency. At 200 Hz, 70 μm of displacement can be generated. The new system was tested on 6 healthy volunteers and one patient with prostate cancer. In the prostate, in vivo, displacements of 50 μm could be measured at 70 Hz. Young's modulus was 7.9±1.0 kPa in the peripheral zone, 9.8±1.4 kPa in the central zone, which agrees with the previous study. For the patient, a peak at 20 kPa was measured in a Gleason 4+3 tumor. In 2015, the study was driven further, with 11 patients

Tristan DERUELLE 51 scheduled with radical prostatectomy. Subsequent histopathology was compared with the elastography results [START_REF] Sahebjavaher | MR elastography of prostate cancer: Quantitative comparison with histopathology and repeatability of methods[END_REF]. A frequency of 70 Hz was used and displacement of 55 μm were measured in the prostate. Shear and compression contributions seems to be considered in their reconstruction algorithm. Reported elasticity was 7.5±1.9 kPa in the peripheral zone, 9.7±3.0 kPa in the central gland, and 9.0±3.4 kPa in the transition zone. Cancers with Gleason's score of at least 3+3 were significantly different from normal tissue. However, prostate cancer was usually stiffer but not always. MRE was able to diagnose prostate cancer with an Area Under the Curve of 0.7. Performance and repeatability were limited by patient motion during long acquisitions, the presence of standing waves and the use of low mechanical frequency. It was suggested that some of these limitations could be overcome by using multiple frequencies and different processing algorithms.

More recently, in 2018, Dittman et al investigated a multifrequency approach using multiple air-pressurized drivers positioned on the pelvic region. One anterior and two posteriors drivers were used [START_REF] Dittmann | Tomoelastography of the prostate using multifrequency MR elastography and externally placed pressurized-air drivers[END_REF]. The study involved 12 healthy volunteers and 5 patients with suspected prostate cancer. Frequencies of 60, 70 and 80 Hz were used. A displacement of 15-35 μm could be achieved in the prostate. The curl operator was used to remove the compression contribution. Similar stiffness between the central gland and the peripheral zone was reported (shear wave velocity was respectively 2.27±0.2, and 2.22±0.23 m/s). Suspected regions showed a higher stiffness (shear velocity between 3.29-3.63 m/s). Figure 3-2 shows an example of the images centered on the prostate. At low frequency, the multifrequency approach did not perform better than using a single frequency. This setup was also used in a study conducted by Asbach et al in 2020 [START_REF] Asbach | In Vivo Quantification of Water Diffusion, Stiffness, and Tissue Fluidity in Benign Prostatic Hyperplasia and Prostate Cancer[END_REF] and Li et al in 2021 [START_REF] Li | Tomoelastography Based on Multifrequency MR Elastography for Prostate Cancer Detection: Comparison with Multiparametric MRI[END_REF]. Asbach involved 39 men (25 with benign prostatic hyperplasia + 14 with prostate cancer). He reported that the shear wave velocity in prostate cancer (3.1±0.6 m/s) was higher than in the transition zone (2.8±0.3 m/s) or tended to be higher than in peripheral zone (2.8±0.4 m/s). Li involved 208 men (53 healthy controls + 82 with benign prostatic hyperplasia + 73 with prostate cancer). He reported that the shear wave velocity in prostate cancer (3.4±0.6 m/s) was higher than in benign prostatic disease (2.6±0.3 m/s) and healthy controls (2.2±0.1 m/s). Both of these studies also measured the loss angle which can be related to fluidity ( = 0 corresponding to pure elastic medium, and for a pure fluid). The use of wave velocity and loss angle improved the diagnostic performance of multiparametric MRI (95% specificity vs 77% for mpMRI alone in [START_REF] Li | Tomoelastography Based on Multifrequency MR Elastography for Prostate Cancer Detection: Comparison with Multiparametric MRI[END_REF]).

CONCLUSION

Measurements of prostate stiffness were performed in various ways: in vivo and in vitro, with healthy or cancerous specimens, with standard mechanical assessment or with ultrasound or MR-based elastography, and at different frequencies. So many different factors yield to a high variability in the reported stiffness values. This is demonstrated in Figure 3-3 which reports the elasticity values found in the literature with various techniques at different frequencies. This shows the difficulty to provide a reliable quantitative value. Besides, values are reported with different units and sometimes just the broad term "elasticity" is used. This makes it difficult to compare results between different studies, especially if different mechanical models are used. It is now recommended to report results as shear wave velocity [START_REF] Barr | WFUMB Guidelines and Recommendations on the Clinical Use of Ultrasound Elastography: Part 5. Prostate[END_REF]. However, the literature tends to show that cancers have a different stiffness than healthy tissue.

The clinical use of prostate elastography is aimed at early detection of cancers. The method must be employed in vivo. Ultrasound and MR-based approaches have been investigated. Both were demonstrated to be feasible and provide promising results. International recommendations guide the use of ultrasound elastography. Nonetheless, elastography may not be used alone as a diagnostic tool, but in combination with other parameters.

Recently, international guidelines for prostate cancer screening and treatment recommend to perform mpMRI prior biopsy. This is the perfect occasion to include MRE in the clinical protocol. The patient being already scheduled with an MRI, a short MRE sequence can be added. The radiologist would thus have all of the mpMRI information, but also a 3D stiffness map. The advantages to perform MRE along with mpMRI are a gain in time, no additional device is needed (apart from the vibrator which may remain in the MR room), and elasticity images can already be registered with the mpMRI images.

However, prostate MRE still needs to be perfected. The prostate is heterogeneous, deep in the body, and surrounded by damping tissues. It is thus difficult to excite and to generate displacements of sufficient amplitude and SNR. In the literature, waves can be generated by different devices, invasive or not. However, most in vivo studies used noninvasive actuators. Wave penetration can be improved by using low frequencies (<100 Hz in the literature). However, the quality of the reconstructed velocity depends on the wavelength. Having too large wavelengths decreases the resolution and precision. Also, special attention should be paid to compression waves. MRI having the fantastic ability to acquire motion in all three directions of space, compression is usually handled by using

Tristan DERUELLE 53 the curl operator. However, the derivatives involved may decrease the SNR and therefore decrease the quality of the reconstructed stiffness maps.

The major challenge is thus to have a device that can generate a wavefield in the prostate with sufficient amplitude and frequency, while limiting the generation of compression.

4

The reader is now aware of the need for an efficient wave generation system. Commercially available systems are not suited for the prostate. This chapter explains that legal and security considerations should be thought of even before designing the prototype system. Then, we present the design of a medical device intended to generate waves in the prostate, and we assess its performance on phantom.

Approximately at the same time as Sahebjavaher published his transperineal setup, a first vibration prototype had been developed at our lab by Hoang Dinh during his PhD thesis [START_REF] Dinh | PhD Thesis -Development and evaluation of quantitative parameters of prostate MRI[END_REF]. This device used a piezoelectric actuator and a transperineal approach. It was tested on 2 healthy volunteers, but tended to get dismantled under its own vibrations. The present work is inspired by that prototype.

REGULATIONS

Even before starting to design a medical device (MD), the manufacturer should get familiar with guidelines and regulations of the targeted country. In the present case, the device is not intended to be placed on the market, nor being CE marked. It is purely for investigational purposes. It is intended to be used only by the manufacturers (ourselves) under the supervision of a doctor who is familiar with the concepts of the study. Potential studies would take place in France, and experiments would be performed at Lyon, at the Hôpital Edouard Heriot (HEH) or at the Centre d'Etude et de Recherche Multimodal Et Pluridisciplinaire en imagerie du vivant (CERMEP -imagerie du vivant).

In France, at the time of this writing, there are two main regulations in effect. The first one is the article L1121-1 of the French public health regulations [START_REF]Article L1121-1 -Code de la santé publique[END_REF]. It defines the legal framework related to research involving human people. In particular, researches are divided into three categories: category 1 groups interventional research, category 2, interventional research with minimal risk, and category 3, non-interventional research where all medical acts are performed as per usual. Operations with minimal risk related to research of category 2 are detailed in the annex 1 in the decree of the 12 th of April 2018 [START_REF]Arrêté du 12 avril 2018 fixant la liste des recherches mentionnées au 2° de l'article L. 1121-1 du code de la santé publique[END_REF]. In particular, in includes collecting data using MRI, and allow external mechanical stimulation in compliance with the CE marked device or with guidelines and good practices. MR elastography uses MRI and external vibrations. For the development of the mechanical stimulation device, guidelines and good practices were detailed by Ehman in 2008 [START_REF] Ehman | Vibration safety limits for magnetic resonance elastography[END_REF], based on the European Directive 2002/44/EC on the exposure of workers to the risks arising from physical agents (vibration). MRE studies with an external vibration device would thus be in category 2.

The second main regulation is the European regulation 2017/745 on medical devices [START_REF]and of the Council of 5 April 2017 on medical devices[END_REF].

It defines the normative framework regarding medical devices, and categories the different types of medical devices. In particular, a device aimed at investigating the prostate state in the view of a diagnosis, prevention, monitoring, prediction or prognosis is a medical device. If it requires a power source, it is an active medical device. Annex VIII of the regulation gives the different rules to classify the medical device as class I, IIa, IIb, or III. A non-invasive active device intended to allow direct diagnosis or monitoring of the prostate and intended to supply (mechanical) energy which will be absorbed by the human body is classified as IIa. For MRE, an external vibration device can be placed in this category.

Then, three things are required to start a clinical trial: the agreement from the legal sponsor, the agreement of the Agence Nationale de Sécurité du Médicament et des produits de santé (ANSM), and the agreement from a Comité de Protection des Personnes (CPP). Then, the administrative path to get each agreement will mainly depend on the category of research of the study and the class of the medical device. Those should be carefully justified. Risk analysis and assessment is also requested in compliance with the norm ISO 14971:2019. Risk assessment concerns not only the patient, but also the operator and the environment (e.g. the studied medical device should not affect other devices in the room, and should be immune external disturbances). In particular, electrical safety and electromagnetic compatibility should comply with the IEC60601-1:2005, and be validated by an external laboratory.

With that in mind, the manufacturer can choose the most convenient way to achieve its intended goal, and then build the medical device with the given constraints (electric isolation, biocompatibility, electromagnetic compatibility…).

DESIGN

Specifications

At first, mechanical drivers were intended to induce shear waves by using transverse motion [START_REF] Muthupillai | Magnetic resonance elastography by direct visualization of propagating acoustic strain waves[END_REF], [START_REF] Braun | Electromagnetic actuator for generating variably oriented shear waves in MR elastography[END_REF]. However, it may be uncomfortable for the patient as it is likely to induce friction, skin irritation, and heat. Instead, drivers acting in compression have been developed. It has also been shown that such longitudinal drivers, if their contact surface is small, can induce diffracted shear wave field at an angle around 30° [START_REF] Yin | Diffraction-biased shear wave fields generated with longitudinal magnetic resonance elastography drivers[END_REF]. Our device was thus designed to have a small contact region, and an adjustable angle for the applied vibrations. To accommodate for the different morphology, the height of the contact point should also be adjustable.

In the literature, vibrations are produced by either invasive or non-invasive methods. For our device, we chose a non-invasive setup, that will induce the waves from the exterior of the body. This is more comfortable and acceptable for the patients. It also simplifies the administrative procedures.

Within non-invasive methods, pelvic and transperineal approaches were tested in the literature. From our experience, a pelvic approach with the standard commercial hepatic driver cannot deliver sufficient amplitude of vibration. A pelvic driver dedicated to the prostate would need much higher vibration amplitudes, which may become uncomfortable for the patient. The transperineal approach, proposed by Sahebjavaher et al [START_REF] Sahebjavaher | Transperineal prostate MR elastography: Initial in vivo results[END_REF], is well accepted by the patients. The perineum is well accessible when the patient lies on the MRI bed, and it is relatively close to the prostate.

The device should be able to deliver sufficient displacement within 0-200 Hz. [START_REF] Sahebjavaher | Transperineal prostate MR elastography: Initial in vivo results[END_REF] suggests that displacement of 0.5 μm within the region of interest with 6 dB SNR are usually sufficient for valid velocity reconstruction. Given the values from the literature, we estimate that 100-200 μm of amplitude at the surface of the skin should be enough.

The control of the vibrations must be the most precise and accurate as possible. Control of the frequency is essential for a perfect synchronization with the acquisition. Very fast response time is also required for the control of phase shifts. Piezoelectric elements are well suited for these purposes. However, their small displacement must be amplified for being detected by MRI. The company Cedrat Technologies (Meylan, France) manufactures Amplified Piezo Actuators (APA). They use a stack of piezo elements within a shell whose shape mechanically amplifies the displacements. It is thus possible to reach a few hundreds of micrometers of displacement. Such an APA is presented Figure 4-1.

Finally, a major requirement is that all materials should be non-ferro-magnetic so that they are not attracted by the magnetic field, nor create artifacts in the images. From our experience, APAs are quite fragile. Breakage usually occurs when the actuator is driven at its resonant frequency. A major concern is that the resonant frequency varies with the load applied, i.e. how hard the actuator is pressed against the patient skin. For example, Figure 4-2 (LEFT) presents the frequency response of an APA400 for different forces, corresponding to weights from 0 kg (in blue), to 1 kg (in red). The expected displacement for this model is around 400 μm. At the resonant frequency, the displacement can nearly double. Beyond the resonant frequency, the displacement drops drastically. The resonant frequency is within our frequency range. On the contrary, the APA100 (Figure 4-2 RIGHT) has a resonant frequency beyond the frequency range of use.

Having this behavior in mind, two APAs were chosen. One APA100MML delivering at least 100 μm over the specified frequency range. A custom APA200MML delivering at least 200 μm over the specified frequency range, with a weight of 1 kg. The resonant frequency being slightly above our frequency range, resonance can be used to actually deliver more than the specified displacement. Each line corresponds on a weight, from 0kg (blue) to 1kg (red), every 0.1kg.

Setup presentation

The assembled device is shown on Figure 4-3. The patient lies on the white board, with the device between his legs. The APA vibrator is inside the aluminum box, and powered through the BNC cable. The rod with the ball is rigidly fixed to the APA and transmits longitudinal vibrations. This ball is in contact with the patient (with his underwear), and thus transmits the vibrations to the patient's perineum. A string can adjust the angle of inclination of the device to reach -30 to +30°. The rotation axis can be placed at two different heights depending on the patient's morphology.

To ease patient installation, the bed board can be detached from the rest of the device.

Once the patient is in place on the board, the device is easily clipped and pinned. The outer plastic housing is made of polycarbonate, with a thickness of 6 mm (2 mm for the front and read panel). Walls are assembled using PVC battens, and brass threaded inserts. Non-magnetic stainless screws are used. The rotation axis is made of brass, and held by A4 stainless beta pins (les-inoxydables.com, France).

Figure 4-4 presents in more details the structure of the vibrating part. From the original prototype, the rod, rod guide, and ball were reused. The rod guide, and the ball were 3D printed with common plastic (Stratasys Vero black). In order to limit the friction, low friction guides were added, and the rod was made of PEEK. Regarding the aluminum housing, the back part and front parts also come from the first prototype. They were machined in aluminum 5083 which does not create artefacts on MR images. In between, an aluminum box was added to host the APA ("center part"). Its size was chosen to fit the APA100MML (the largest), or a stack of two smaller APAs. It is interchangeable if future upgrades are required. Here, a Bopla Filotec case (reference 97118080.HMT1) is used. It is made of aluminum AlMgSi0.5 which was found not to create significant artifacts. Other aluminum, copper, or brass alloys were tested (see Appendix B). The general setup of the MRE experiment is explained in Figure 4-5. The computer (PC) control the wave generator (GEN) by sending frequency, amplitude and burst parameters commands. The output signal is amplified by a power amplifier (AMP) which supplies the APA in the device (Vibrator). The power amplifier, LA75C, was also made by Cedrat Technologies to meet the APA power requirements. The APA embeds a strain gauge to provide a feedback of its deformation. This signal also comes through the amplifier and can be visualized with an oscilloscope (OSC) on the PC. For the synchronization, between the acquisition and the vibration, the MR system can provide a sync signal. This signal is sent to the wave generator. An isolation step is added in-between. As of today, our setup is designed to synchronize with GE and/or Siemens MRI scanners. On General Electric scanners, the Spin-Echo (SE) Echo Planar Imaging (EPI) MRE pulse sequence provided by Richard Ehman (Mayo Clinic, USA) outputs a TTL trigger signal on each Repetition Time (TR), i.e. for every RF pulse of the MRI scanner. This J6 trigger signal (Figure 4-6, LEFT) should be connected to the trigger input of the wave generator (through the dedicated isolation box). On Siemens scanners, the SE-EPI MRE pulse sequence provided by Ingolf Sack (Charité Berlin, Germany) uses the 10 MHz sine clock that is built into the MR scanner (Figure 4-6, RIGHT). This clock signal should be connected to the external reference input of the wave generator (through the dedicated isolation box). The use of the clock signals requires that the vibration frequency is set with very high precision on the wave generator, and that the wave generator can use an external clock. In the frequency range 40-140 Hz, the required precision is typically with 6 digits after unity. The Tektronix AFG 3022C was found to achieve the required accuracy.

Isolation considerations

To isolate from the main grid, and make the whole chain unearthed, all electrical parts are powered through an ISO60601-compliant isolation transformer (Tripp Lite IS1000HGDV). The ISO60601 standard describes the requirements for medical devices electrical safety. The additional isolation step (ISO) isolates the sync signal of the scanner, from our setup.

Development of a medical device

Tristan DERUELLE 61 This way, both sides have different grounds, and any issue from one side does not affect the other. Our setup has thus a floating ground and is isolated from all other entries.

The isolation of the 10 MHz sine signal was done by making a simple transformer. Electronic wire with insulation of 6 kV c.c. was chosen (RS PRO UL3239). The ratio of turns of 1:1 was done on a ferrite (TDK HF70T16X8X12), as shown on Figure 4567. A different approach was used for the trigger signal. TTL are short pulses of 0-5 V. It is thus not possible to use the sine isolator. For this TTL signal, a digital isolator was used (Texas Instrument ISO7821DW, 8 kVPK insulation). This component requires to be powered on both sides. For that, two AC/DC converters (MeanWell MPM-05-5, medical grade power supply) were used. By doing so, the input signal and supply power share the same ground, and are isolated from the output signal and supply power. Capacitors were used to clean the power supplies and output signal. A picture of the device is shown on Figure 45678. For visual aid, the position of the AC/DC components were drawn in light blue, and the digital isolator in green. Regarding biocompatibility, common plastics were used. The ball will not be in contact with the patient, but with his underwear. For extra protection, hospital bed paper sheets will be placed on the device.

Tristan DERUELLE 63 For the MRE acquisition, the APA100 was supplied with a sinusoidal vibration of 2 V peakpeak. This represents 1/16 of the maximum power, and 1/4 of the maximum displacement amplitude. This way, the generated displacement field is expected to be similar to what can be achieved in the prostate. A spin-echo EPI MRE sequence was used (MayoClinic, USA). The acquisition parameters were a field of view of 20 cm x 20 cm x 8 cm with an acquisition matrix of 80x80 with 32 slices (2.5 mm isotropic voxels). The motion was acquired in all three directions of space, and with 6 MEG directions (+x, -x, +y, -y, +z,z). Motion sensitivity was improved by using two pairs of MEG, resulting in a motion encoding of 8.78 μm/rad. Motion was captured for 12 different time offsets per cycle. Echo time was 83.7 ms. Repetition time was 3840 ms.

In Figure 4-11, the first two rows represent the raw phase images captured with positive and negative polarities. Waves are already discernable, but there are large phase wraps. The third row presents the phase difference. Waves are already more visible. There are still some phase wraps. The last row is the phase difference, unwrapped with quality guided algorithm. The wraps are correctly corrected.

Tristan DERUELLE 65 The frequency content is analyzed in order to assess the frequency purity of the vibrations. When driven at 100 Hz, the system should only produce 100 Hz vibration, with no or negligible harmonics. The temporal Fourier Transform is taken and the magnitude is summed over the whole volume. Figure 4-12 plots the summed magnitude over frequency for the unwrapped phase differences for the different directions. All of the energy is concentrated at the excitation frequency. No signal distortion is visible in the harmonics, only the fundamental frequency and white noise. Note that on the X component, the spectrum is not perfectly flat on the harmonics. This is still negligible.

From the amplitude of the frequency contents, the phase signal to noise ratio (phase SNR) can be estimated as the ratio of signal power over noise power. On average, the phase SNR in the X, Y and Z displacement fields are respectively 6.6, 7.1 and 9.0 dB.

Reconstruction of velocity maps

To reconstruct the velocity maps, the wave images are first temporally filtered at the excitation frequency (Figure 4-13, top row). This filtered field is referred to as "Full Field".

To remove any compression contribution, the curl operator is applied to this wave field (Figure 4-13, bottom row).

Tristan DERUELLE 67 The propagating wavefield can be viewed on Figure 4-14 for the X component. Up to 15 μm of displacement are measured close to the vibrator, and 2 μm far from the vibrator. This is the typical amplitude of displacement that is generated in the prostate in the literature.

It is to be noted that when looking at the coronal image (Figure 4-15), waves are traveling with an angle around 30°. This is the typical shear diffraction pattern generated by longitudinal driver as described by Yin [START_REF] Yin | Diffraction-biased shear wave fields generated with longitudinal magnetic resonance elastography drivers[END_REF]. Spatial frequency is reconstructed with 3D LFE, and converted to velocity as = . The velocity maps without and with the curl operator are shown on Figure 4-16. The reconstruction from the full field can recover the four inclusions of different stiffnesses. When using the curl field, only the soft ones are visible. Without using the curl operator, a stiffening artefact is present close to the vibrator. This is the well-known compression artefact. Using the curl operator corrects this artefact.

CONCLUSION

A new device is presented to generate waves in the prostate. Attention was put into the ergonomics and safety. Materials were carefully chosen to avoid any image artefacts. The actuator was chosen to be able to deliver enough displacement across the frequency range of interest, while limiting the risk of breakage due to resonance and/or applied load.

The device was tested on a CIRS phantom. The presence of the device did not impair the images. The generated wave field was relatively close to what can be expected in the prostate. All four inclusions were visible in the reconstructed velocity maps. However, the well-known compression artefact (due to the vibrator) was visible. It was removed by the use of the curl operator. The derivatives involved in the curl operator decrease the phase SNR. The quality of reconstructed velocity maps is impacted as a global softening, and the difficulty to recover hard inclusions. This shows the need for better algorithms to handle shear and compression. The previous chapter showed the limits of the curl operator, especially in the presence of noise. This new chapter exposes a new method to split the measured displacement field into the shear and compressional fields. Having the shear displacement fields without additional derivatives can yield to a better and noise-resistant estimation of the shear velocity.

ALGORITHM DERIVATION

In MRE, the displacement can be acquired in all points in a 3D volume. More specifically, the full displacement vector ( , , , ) is known (x,y,z are the three spatial components, and t is time).

The Helmholtz theorem [START_REF] Aris | Vectors, Tensors and the Basic Equations of Fluid Mechanics[END_REF], [START_REF] Gui | A RIGOROUS AND COMPLETED STATEMENT ON HELMHOLTZ THEOREM[END_REF] states that any smooth and rapidly decaying vector field can be decomposed into the sum of an irrotational (curl-free) vector field and a solenoidal (divergence-free) vector field:

= + (5-1)
In the field of mechanics, if the vector field represents a displacement field, then the curl-free component

= (5-2)
represents the underlying dilatational field, i.e. the displacement field caused by the pressure field. The divergence-free component

= (5-3)
represents the underlying shear field, i.e. the displacement field caused by shear stress.

is a scalar potential, and is a vector potential.

The goal is to compute the vector potential and scalar potential . Then, from the knowledge of these potentials, the shear and compression fields can be estimated.

Field separation

Tristan DERUELLE 71 where L is the sparse Laplace matrix, x is the vector unknown ( reshaped as a vector), and b is a vector filled with the right-hand side of the equation. All three components of vector are independent, hence the vector Poisson equation can be solved as 3 independent scalar Poisson's equations.

Scalar potential

In a similar way to , can be expressed as follows.

Applying the divergence operator to equation (5-1) yields:

= ( ) +

(5-12)

= +

(5-13)

The divergence of the shear field is zero. Therefore:

= = ( ) (5-14)
This yields again to the well-known scalar Poisson's equation, solving for and knowing :

= (5-15)
Equation (5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF] can also be solved using any known technique for solving the discrete Poisson equation.

Here also, equation (5)[START_REF] Defossez | Estimations nationales de l'incidence et de la mortalité par cancer en France métropolitaine entre 1990 et 2018[END_REF][START_REF] Stamey | Localized prostate cancer. Relationship of tumor volume to clinical significance for treatment of prostate cancer[END_REF][START_REF] Mottet | EAU -EANM -ESTRO -SUR -ISUP -SIOG Guidelines on Prostate Cancer[END_REF][START_REF] Smeenge | Role of transrectal ultrasonography (TRUS) in focal therapy of prostate cancer: report from a Consensus Panel[END_REF][START_REF] Roberts | Prostate Biopsy-related Infection: A Systematic Review of Risk Factors, Prevention Strategies, and Management Approaches[END_REF][START_REF] Gleason | Classification of prostatic carcinomas[END_REF][START_REF] Salomon | Le score de Gleason pour les nuls[END_REF][START_REF] Wyatt | Heterogeneity in the inter-tumor transcriptome of high risk prostate cancer[END_REF][START_REF] D'amico | Biochemical Outcome After Radical Prostatectomy, External Beam Radiation Therapy, or Interstitial Radiation Therapy for Clinically Localized Prostate Cancer[END_REF][START_REF] Edge | The American Joint Committee on Cancer: the 7th Edition of the AJCC Cancer Staging Manual and the Future of TNM[END_REF] can be solve by discretizing the Laplace operator and solving for .

Deduction of Up and Us

Once and are known, and can be determined using equations (5-2) and (5-3).

Alternatively, one field is computed using one of these equations, and the other is deduced as:

U = U U (5-16)
Or U = U U

(5-17)

Direct solution

Shear displacement field

Alternatively, equation (5-10) can be solved in the frequency domain.

We use the notation FT for the Fourier Transform, kx, ky, kz for the complex spatial wave numbers, and k 2 = kx 2 + ky 2 + kz 2 .

Many authors implemented a Finite Difference Time Difference (FDTD) scheme but either don't implement viscosity or have a compressional velocity in the same order as the shear velocity [START_REF] Catheline | Tomography from diffuse waves: Passive shear wave imaging using low frame rate scanners[END_REF], [START_REF] Jiménez | Ultrasonic evaluation of the hydration degree of the orange peel[END_REF]- [START_REF] Pham-Thi | Two-dimensional complex shear modulus imaging of soft tissues by integration of Algebraic Helmoltz Inversion and LMS filter into dealing with noisy data: a simulation study[END_REF]. Finite Element Method (FEM) based software such as ABAQUS, COMSOL, ANSYS or LS-DYNA [START_REF] Romano | On the noninvasive determination of material parameters from a knowledge of elastic displacements theory and numerical simulation[END_REF]- [START_REF] Mcgrath | Evaluation of wave delivery methodology for brain MRE: Insights from computational simulations[END_REF] implement viscoelastic models. Although those methods are attractive (but require a license), the published studies usually had a ratio of 10 between compressional and shear velocities. Reaching higher ratio is challenging because we usually want to model the wave propagation in free space. This is achieved by increasing the size of the computational grid so that the waves never reach the boundaries, leading to computational load ever bigger. There are studies implementing advanced models for that problem [START_REF] Ye | A selective mass scaling method for shear wave propagation analyses in nearly incompressible materials[END_REF], [START_REF] Ye | Simulation of nonlinear transient elastography: finite element model for the propagation of shear waves in homogeneous soft tissues[END_REF] but their implementation and use become complex.

Our choice turned to the MATLAB toolbox k-wave [START_REF] Firouzi | A first-order k -space model for elastic wave propagation in heterogeneous media[END_REF], [START_REF] Treeby | Modelling elastic wave propagation using the k-Wave MATLAB Toolbox[END_REF]. It implements a Kelvin-Voigt viscoelastic model on a pseudospectral time-domain scheme. Spatial operations are performed in the Fourier domain, while temporal operations are performed using finite difference scheme. K-wave also implements a special Perfectly Matched Layer (PML) using anisotropic absorption. This can efficiently decrease the size of the PML while limiting the amplitude of the reflected waves at the grid edges. The computation can be thus very effective. The toolbox also allows the user to define a custom 2D or 3D geometry with user-defined material parameters. A previous study used this toolbox for elastography simulations [START_REF] Prieur | Simulation of shear wave elastography imaging using the toolbox 'k-Wave[END_REF]. They compared the results with analytical solution and FEM software, and found good accuracy.

Simulations for field separation

Slow compression

To check the capability of the field-separation algorithm, a basic simulation was performed. It consisted in a 3D homogeneous medium with a shear velocity of cS=2 m/s, and a compressional velocity of cP=3 m/s. The source was modelled as a circular membrane (vibrating at 100 Hz) with a diameter of 8cm, similar to that used by Rouviere et al for kidney MRE [START_REF] Rouvière | Magnetic resonance elastography of the kidneys: Feasibility and reproducibility in young healthy adults[END_REF]. Such a vibrator can take the form of a small drum, as shown in Figure 5-1. The membrane is on the bottom side, and can be excited by connecting to an acoustic source.

Field separation

Tristan DERUELLE 75 The behavior of a vibrating membrane is well described [START_REF] Leissa | MEMBRANES[END_REF], [START_REF] Contributors | Vibrations of a circular membrane[END_REF]. For the fundamental mode, the displacement can be written as:

( , ) = sin( )

(5-33)
Where A is the maximum displacement, J0 is the Bessel function of the first kind of order 0, is its first root, and is the radius of the membrane. An example of the resulting membrane is given in Figure 5-2. The field of view of the simulation was 10 cm x 10 cm x 10 cm with 100 pixels in all directions. Perfectly Match Layer (PML) was added on all sides with a thickness of 34 pixels. The particle velocity is simulated and integrated into absolute displacement. The resulting displacement fields were then filtered at the excitation frequency.

The shear field was computed using the finite difference approach, and the compression field as the difference between the full field and the shear field.

We define our frame of reference having the x-direction from left to right, y-direction from top to bottom, and z-direction from front to back.

The resulting displacement fields are shown in Figure 5-3. It can be seen that the full displacement (UX, UZ) presents a mix with some variability in the wavelength (especially visible at the center). In subsequent rows of the same figure, the field is split into the shear field, with short wavelengths, and the compressional field with longer wavelengths.

To verify that these fields are indeed associated with shear and compressional displacement, we check the following constitutive relations:

= (5-34) = (5-35) = 0
(5-36) = 0

(5-37)

The amplitudes of the curl and divergence are also plotted on Figure 5-3. We can see that the curl of the full displacement and of the shear displacement are similar, as expected. Similarly, the divergence of the full field and of the compressional field are similar too. The divergence of the shear field and the curl of the compressional field are very close to zero. The constitutive equations are satisfied. From these results, the algorithm can effectively separate the full field into shear and compression fields.

Fast compression

The next simulation addresses two questions. The first one is to assess the eventual loss in spatial resolution induced by the wave separation algorithm. The second is to demonstrate that the wave separation algorithm is also effective when the compression wavelength is much longer than the field of view.

This simulation consisted in a bilayer medium (see Figure 5-5), with parameters as listed in Table 5-1. In particular, a compressional velocity similar to that of soft tissues was used.

The source was modelled as the same circular membrane as the previous simulation (vibrating at 100 Hz) with a diameter of 8cm.

The field of view was 10 cm x 10 cm x 10 cm with 86 pixels in all directions. Perfectly Match Layer (PML) was added on all sides with a thickness of 34 pixels. Computation time was in the order of 2-3 weeks. The resulting displacement fields are shown in Figure 56. The full field and shear field are relatively similar, but present the main difference close to the source (top of the field of view -FOV). In that region, the wavelength from the full displacement field appears longer than that of the shear. This effect is due to the superposition of the compression field. When looking at the compressional displacement field, we can see the large and rapidly decaying wavelength arising from the membrane. In this scenario, the diameter of the vibrator is small compared to the compressional wavelength (d/ = 0.1/15). Hence, for compression at least, the source behaves like a point source. The compression field should therefore be compared to that of a point source. The compression field (Y component) is indeed similar to the one shown in Figure 5-7 from the article by Catheline et al [START_REF] Catheline | Diffraction field of a low frequency vibrator in soft tissues using transient elastography[END_REF].

The amplitudes of the curl and divergence are also plotted on Figure 56. The constitutive equations are also satisfied. However, the divergence field should reflect the presence of the compression wave (with a wavelength of about 15 m), but small wavelengths are present. They are exactly of the same size as in the curl field, suggesting numerical errors. This shows the limit of the simulation.

Field separation

Tristan DERUELLE 79 Then, the resolution is taken as the distance needed to cover the 5-95% interval of the step and can be computed as = 2 log . .

/ .

The curl, due to the derivatives, is known to increase the weight of higher spatial frequencies, and thus to increase noise and boundaries detection. For comparison of similar degrees of derivation, the curl of the shear field was also computed. The full field and the shear field are thus on the same zero order of derivative, and the curl field (curl of full field) and curl of shear field are on the first order of derivative. The corresponding profiles and fitted sigmoid are shown Figure 56789. A resolution of 1.3cm is measured for the full field, 1.7cm for the shear field, 0.4cm for the curl field, and 0.6cm for the curl of the shear field. It is interesting to note, that, for the full field, the fit also highlights the overestimation of the top layer (2.4m/s). The shear field corrects this bias. Regarding the fit from the curl profiles, they are almost identical. Only the undershoot is less present in the curl of the shear field. We find again that taking a derivative (in the curl) favors the recovery of boundaries. However, this comes at the cost of amplifying noise if present.

Noise immunity

To test the noise immunity, the previous bilayer simulation is used. Different SNR were reached by adding noise to the displacement images: 30, 40, and 60 dB SNRpeak. Peak is referred to the maximum amplitude of the displacement vector max , , .

Figure 5-10 presents an example of noisy data and the associated SNR map. On the left, noise was added to the X component of the displacement vector ( ) to reach 30 dB SNRpeak. On the right, the corresponding SNR map is presented. It is important to note that the SNR is not uniform within the image. As it is referred to the maximum amplitude of the displacement vector, the SNR is at most equal to 30 dB. It is maximum close to the vibration source, and decreases with depth. It is clear that the curl tends to decrease the SNR while the computed shear field has an SNR similar or even slightly better than the full field. even at 40dB, the reconstruction is severely hampered. Regarding the use of the proposed method, the contrast between the two layers is recovered, even at 30 dB. On simulations, the retrieval of the shear field provides a shear velocity estimation as precise and accurate as the curl, while being much more noise resistant. Regarding the resolution, the full field is slightly better by 4mm.

Note that the noise performance using the curl field could be improved (at the expense of spatial resolution) by applying a mean filter or any other processing. An example is provided in Figure 5-14. A spatial moving average of size 7x7x7 pixels was applied prior to the curl operator in order to smooth the noise. No spatial filter was applied to the shear field. Under these conditions, the smooth curl field and the shear field now yield velocity maps having similar SNR and similar spatial resolution (~1.6 cm).

Field separation

Tristan DERUELLE 85

As a conclusion, once spatial resolution has been set equal, we can see that the wave separation algorithm does not improve SNR in the wave velocity images, when compared to the curl algorithm. It can be a little disappointing. Nevertheless, the wave separation algorithm has other benefits. It provides new information that was not available: the shear displacement field and the compression displacement field.

Phantom acquisition

The algorithm was then tested on real MRE experiment. MR experiments were performed on a 1.5 T GE scanner. A CIRS 049 elastography phantom was imaged. T2w image, and the expected shear velocity map are presented For this study, in order to assess the performance of the algorithm alone, a common MRE setup was used. The vibrator consisted in a speaker within a resonant box from which a pipe was connected to transmit acoustic waves to an 8-cm diameter drum-like vibrator.

The vibrator was positioned onto the phantom. The acoustic waves generated by the speaker are amplified by resonance in the box, transmitted in the pipe, and make the drum vibrate on the subject. A schematic representation of the setup is shown Figure 5-16. A spin-echo EPI MRE sequence was used (MayoClinic, USA). The acquisition parameters were a field of view of 24 x 24 x 7.2 cm 3 with an acquisition matrix of 80x80 with 24 slices (3 mm isotropic voxels). The excitation frequency was 100 Hz. The motion was acquired in all three directions of space, and with opposite polarities. The propagation was acquired with 6 MEG directions (+x, -x, +y, -y, +z, -z) and 3 time offsets per cycle. Different levels of phase SNR were achieved by acquiring with 2 pairs of motion encoding gradients, 1 pair of MEG, and one pair with half-amplitude. They are respectively referred to as 0 dB, -3 dB, and -6 dB (here SNR is relative and defined using the 2MEG acquisition as the 0 dB baseline). Echo time was respectively 81.1, 61.1, and 61.1ms. Repetition time was 2880ms.

First, phase difference images were computed for each spatial component. They were then unwrapped using a 3D quality guided algorithm. Due to unwrapping issues, the three first slices were disregarded. In the following, only 21 slices are considered. Data were then filtered at the excitation frequency. Finally, the shear and compression fields were computed. The fields are presented in Figure 5-17. It is interesting to see that the wave fields are similar to those simulated (Figure 56). The compression field presents the same large and decaying wavelength pattern. Surprisingly, short wavelengths are also visible in the compression field, suggesting the existence of two compression waves, a fast one and a slow one.

Field separation

Tristan DERUELLE 87 The wave velocities from the full field, curl field, and shear field were reconstructed with 3D LFE, and are shown in Figure 5-18. Because of the large wavelength in compression, its velocity was not reconstructed. The first thing to notice is that in the reconstructions from the full field, the velocity is systematically overestimated right in front of the vibrator. This is not the case with the curl field and the computed shear field. This overestimation is due to the compression in that region. That contribution is successfully removed by the curl or the computation of the shear field. The next thing to notice is that when decreasing the sensitivity (so decreasing the phase SNR), the reconstructions from the curl field are degraded. This is visible as a global softening and difficulty to recover the hard inclusions. This is also the case for the full field at the sensitivity of -6 dB relative SNR. On the other hand, the reconstructions from the shear field appear much less affected. With high sensitivity, all 4 inclusions are visible, regardless of the type of displacement field. The velocity in the different regions was measured and are reported in Table 5-2. Measurements were also performed using ultrasound with an Aixplorer ® scanner, and the values are also reported here, along with the reference values provided with the phantom certificate. The soft inclusions (Type I & II) are overestimated, and the hard ones (Type III & IV) are underestimated. This is mainly due to reconstructing the velocity inside a small region. In front of the vibrator (see "near field" row), due to compression, the full field presents a higher value than in the rest of the background. This is not the case with the curl or shear field.

Field separation

Tristan DERUELLE 89

ORIGIN OF THE SLOW COMPRESSION WAVE

The origin of the slow waves in the compression field (Figure 5-17) is unknown. Several hypotheses were tested to try to explain this phenomenon.

We first hypothesized that the slow wave pattern was an artefact generated by the wave separation algorithm. Therefore, we took the wave separation algorithm out of the equation, by using only the divergence operator applied on the full field. The divergence of the full field and of the compression field do present the same short waves, as shown in Figure 5-19. Hence, the slow compression waves are not an artefact produced by the wave separation algorithm. We also hypothesized that the slow waves may arise from numerical errors in the computation of the divergence operator. In order to assess this hypothesis, the divergence was computed using two different methods: using a centered finite difference scheme, and with a Fourier based approach (in k-space). With both methods, the same slow waves are present, as shown in Figure 5-20. This suggests that the slow compression waves are not related to computation errors. The next possible reason is an artefact related to the MRE acquisition. To verify that, the experiment was repeated with a 1.5 T Siemens scanner, using a spin-echo EPI MRE sequence (Charité Berlin, Germany). The experiment took place at the CERMEP, using the same excitation setup. Vibration frequency was also 100 Hz. The acquisition matrix was 128x128 with 59 slices, and a FOV of 32x32x14.75 cm 3 (2.5 mm isotropic voxels). Echo time was 80 ms and repetition time was 7950 ms. All three components of the displacement vectors were acquired, for 18 time offsets per cycle. We compare the divergence field from data acquired from two different scanner (and manufacturer), with two different sequences (still both are spin-echo EPI), and with different resolutions. Figure 5-21 presents the two divergence fields. Both fields are similar and present short waves. Even though we cannot rule out an effect of the spin-echo EPI sequence type, these results suggest that the slow compression waves are indeed physically present. Finally, some phase interslice artefacts are present. Those aretefacts can be removed using dejittering algorithms [START_REF] Barnhill | Fast Robust Dejitter and Interslice Discontinuity Removal in MRI Phase Acquisitions: Application to Magnetic Resonance Elastography[END_REF]. Tristan DERUELLE 91 performed using the GE scanner, without and with the use of dejittering. However, as shown in Figure 5-23, jitter is not the cause of the presence of short waves in the divergence. Physical explanations for an additional slow compression wave are sparse. In our situation, the wavelengths are in the same order of magnitude as the object. Consequently, guided waves may be present, as symmetric or antisymmetric modes. The theory predicts that these modes can travel slowly [START_REF] Royer | Elastic waves in solids: Free and Guided Propagation[END_REF], [START_REF] Laurent | In-plane backward and zero group velocity guided modes in rigid and soft strips[END_REF]. However, these modes are a mix of compression and shear, and it is therefore unclear whether these guided modes are associated with slow compression waves or not.

To our knowledge, the only theory that predicts the presence of a second slow compression wave is therefore poroelasticity. It will be developed in the next chapter.

CONCLUSION

A new algorithm of field separation is presented. From a displacement field, the shear and compression fields can be retrieved. Here, the shear field was computed using the finite difference approach. Compression field was computed as the difference between the original field and the shear field. These methods proved to be noise resistant, accurate and precise. Alternatively, the whole process could be performed in the frequency domain. Also, solving the Poisson's equation could be performed with other kind of methods (e.g. iterative algorithms, multigrid solvers...). Further studies are required to analyze and compare the possible performances.

This technique of field separation is beneficial for MRE as it can remove the compressional artefacts that are produced by common excitation system. This results in a more accurate estimation of the shear velocity. The method is noise resistant and allows for recovering structures at low SNR. This can be valuable for reconstructing velocity in deep and not accessible organs where it is difficult to have waves with high amplitude. The downside is that the resolution is slightly decreased. However, thanks to the noise-immunity, it can be possible to use higher frequencies.

Now that separate shear and compression displacement fields can be measured, a new quality metric could be defined as the ratio between shear and compression amplitudes.

A map of shear/compression ratio could be computed. This new quality metric could then be used to remove the compression wave locally, only where needed. Its use in combination with the octahedral shear strain SNR [START_REF] Mcgarry | An octahedral shear strain-based measure of SNR for 3D MR elastography[END_REF] would require further studies.

Finally, the experiments highlighted an unexpected phenomenon. For tissues that are considered purely (visco)-elastic, like the CIRS phantom, the compressional wave is so large that it cannot be recovered. However, in our MRE experiment, the divergence of the acquired field was not zero. In addition, the computed compression displacement field presented slow waves, visible on top of the expected fast compression wave pattern. This result highlights the potential of this algorithm for the analysis of more complex wave phenomena such as poroelasticity. In the next chapter, we will investigate slow compression waves.

Tristan DERUELLE 93

CHAPTER

6 Slow compression wave elastography

6

The previous chapter presented an algorithm to separate the shear and compressional fields. It explained that in elastic tissues, the compressional velocity is too large to be recovered by MRE, and the compressional field is typically unwanted. However, an unexpected slow compression wave was observed. This can be the case in porous materials, where, in addition to the fast compressional wave, a secondary and slow compressional wave propagates. In biological soft tissues, fluid exchange is omnipresent. Whether at the cell or vessel scale, porosity is present. Consequently, a poroelastic model could be well suited for the description of tissues characteristics. This chapter presents those poroelastic waves, and the possible implications using the field separation algorithm.

POROELASTICITY

A porous medium consists in a solid matrix saturated by a fluid, similar to a sponge. This type of media was first studied in the field of geophysics [START_REF] Terzaghi | Theoretical Soil Mechanics[END_REF] with an interest in hydrocarbon fluid in rocks. The main characteristic of porous media is that three different waves propagate: two compressional (one fast, P-I, and one slow, P-II) and one shear wave (S). The first theoretical analysis were performed by Frenkel [START_REF] Frenkel | On the theory of seismic and seismoelectric phenomena in a moist soil[END_REF], [START_REF] Frenkel | On the Theory of Seismic and Seismoelectric Phenomena in a Moist Soil[END_REF], but were completed by Biot [START_REF] Biot | Theory of Propagation of Elastic Waves in a Fluid-Saturated Porous Solid. I. Low-Frequency Range[END_REF], [START_REF] Biot | Theory of Propagation of Elastic Waves in a Fluid-Saturated Porous Solid. II. Higher Frequency Range[END_REF] with viscous effects.

This second compressional wave arises from the fact that the medium is in fact an imbrication of two media that interact between each other. For example, the fluid can navigate between pores, local compression can make the fluid flows further where it can result in local expansion of the fluid space. When one of the two spaces is stressed, it impacts the other space. These microscopic interactions result in the formation of the second compressional macroscopic wave.

The wave velocities depend on many intrinsic parameters as described in Table 6-1. The detailed derivation of the velocities can be found the work by Biot [START_REF] Biot | Theory of Propagation of Elastic Waves in a Fluid-Saturated Porous Solid. I. Low-Frequency Range[END_REF] or Smeulders [START_REF] Smeulders | Experimental Evidence for Slow Compressional Waves[END_REF]. [START_REF] Smeulders | Experimental Evidence for Slow Compressional Waves[END_REF], [START_REF] Plona | Observation of a second bulk compressional wave in a porous medium at ultrasonic frequencies[END_REF]. More recently, Aichele [START_REF] Aichele | Capturing the Shear and Secondary Compression Waves: High-Frame-Rate Ultrasound Imaging in Saturated Foams[END_REF] observed the three waves in a melamine foam phantom using ultrafast ultrasound imaging. The interest in medical imaging is rising [START_REF] Konofagou | Poroelastography: imaging the poroelastic properties of tissues[END_REF]- [START_REF] Perriñez | Contrast detection in fluid-saturated media with magnetic resonance poroelastography[END_REF]. The slow compressional wave has been observed in bone ex vivo [START_REF] Lakes | Slow compressional wave propagation in wet human and bovine cortical bone[END_REF], [START_REF] Hosokawa | Ultrasonic wave propagation in bovine cancellous bone[END_REF] and in vivo [START_REF] Yamamoto | Measurement of human trabecular bone by novel ultrasonic bone densitometry based on fast and slow waves[END_REF]. However, in soft tissues such as brain or liver, their existence remains hypothetical. Hirsch, via MRE, reports that different levels of compression could be measured in the lungs [START_REF] Hirsch | Measurement of vibration-induced volumetric strain in the human lung[END_REF]. Lungs are still quite different from common soft tissues as it contains air. More interestingly, in a different study, he showed cyclic compression in the brain which suggests the presence of poroelastic behavior [START_REF] Hirsch | Compression-sensitive magnetic resonance elastography[END_REF], [START_REF] Hirsch | In vivo measurement of volumetric strain in the human brain induced by arterial pulsation and harmonic waves[END_REF].

When performing an MRE acquisition on a porous media, the MRI will sense the summed displacement of all three waves: the fast compression wave (~1500m/s), the shear wave, and the slow compression wave. The algorithm presented in the previous chapter should be able to split the displacement field into the shear field on one hand, and a field with the contributions from both the fast and slow compression waves on the other hand.

CHARACTERISTIC FEATURES OF THE DIFFERENT WAVES

Before going any further, it is interesting to describe the expected characteristics of each wave.

The shear wave has a velocity of a couple m/s (so a wavelength in the order of the centimeter at 100 Hz). It is expected that the curl of the shear wave is non-zero while its divergence is zero.

The fast compression wave, or P-I wave, is expected to have a velocity in the order of 1500 m/s. This corresponds to a wavelength in the order of 15 m at 100 Hz. The excitation source being much smaller (around 10 cm), it behaves like a point source. The behavior Tristan DERUELLE 95 of a compression wave generated by a point source was described by Catheline [START_REF] Catheline | Diffraction field of a low frequency vibrator in soft tissues using transient elastography[END_REF] as was shown in Figure 5-7, and by the Green's functions [START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF], [START_REF] Catheline | Longitudinal shear wave and transverse dilatational wave in solids[END_REF], [START_REF] Sanchez-Sesma | Retrieval of the Green's Function from Cross Correlation: The Canonical Elastic Problem[END_REF]. It presents a rapidly decaying large wavelength pattern. As a compression wave, its curl is expected to be zero while its divergence is expected to be non-zero. However, because of the very large wavelength, its divergence is expected to be undetectable (almost zero).

The slow compression wave, or P-II wave, is expected to have a velocity similar to, ot only slightly higher than. The expected wavelength at 100 Hz is also in the centimeter scale. Its curl should be zero, while its divergence should be non-zero.

In practice, computing a curl or a divergence involves the computation of derivatives. Derivatives act as high pass filters, removing low frequencies (large wavelengths), and emphasizing higher frequencies (short wavelengths, boundaries, and noise). In the divergence field, the influence of the P-I wave may thus be negligible compared to the influence of the P-II wave.

This chapter aims at separating and observing shear and compression waves, as well as highlighting the presence of the P-I and P-II waves. We investigated four different materials: agar gel, melamine foam, silicone, and finally, human kidney graft in vivo. In the following experiments, the speaker setup was re-used (see 5.3.3).

CASE OF AN AGAR GEL

The first test phantom was a 2% agar gel (98% water). The gel temperature was not monitored, but was somewhere between 4°C and 20°C. Acquisitions were performed at the CERMEP on a 1.5 T Siemens scanner. First, a T2-weigthed image, and a diffusion acquisition were performed. The diffusion acquisition was performed with b-values of 50, 800, and 2000 s/mm². The resulting images are shown Figure 6-1. The Apparent Diffusion Coefficient in the agar gel is 1947±74 μm²/s. This suggests that even though water is constrained inside the agar matrix, molecules can still move relatively easily (as a comparison, pure water at 37°C has an ADC around 3000 μm²/s, brain white matter between 670-800 μm²/s). Regarding elastography, the acquisition was a multi-slice SE-EPI MRE pulse sequence with 2.5 mm isotropic voxels. The field of view was 32 x 32 x 10.25 cm 3 with an acquisition matrix of 128x128x41. Acquisitions were performed for 18 time offsets per cycle, with an excitation frequency ranging between 39 and 132 Hz. Echo time was 80ms and repetition time was 6000 ms. All three components of the displacement vectors were acquired.

Data were unwrapped with a quality guided algorithm, and filtered at the excitation frequency. The resulting field is referred to as 'Full Field'. The full field was separated into the shear and compressional fields. Corresponding curl and divergence fields were also computed. Velocities were reconstructed with 3D LFE.

The different fields (full, shear, and compression) are presented in Figure 6-2 for an excitation frequency of 100 Hz. Several things are to notice. The first one, is that the shear and full displacement fields are relatively similar. The main differences appear in the vicinity of the vibrator, where the wavelengths appear longer in the full field than in the shear field. This discrepancy can be attributed to the presence of the compression wave, whose amplitude is non-negligible in that region. In the compression field, the P-I pattern is visible: large wavelength, rapidly decaying. The next thing to notice is that the curl of the shear field is very similar to the curl of the full field, and the divergence of the shear field is zero. Similarly, the divergence of the compression field is very similar to that of the full field, and the curl of the compression field is almost zero. It is interesting to note that the divergence of the full field is not null, suggesting the presence of a slow compression wave. For a better visualization of the propagation of the compression, Figure 6-3 presents the divergence field for 8 time points over a cycle. Wave velocity was assessed for the different fields. In order to limit noise effect, a spatial mean filter of 3x3x3 pixels was applied to the curl and divergence fields. The full, shear, and compression fields were not spatially filtered. The corresponding velocity maps are shown Figure 6-4 for the different frequencies.

One thing to notice here is that the shear and curl fields both yield the same velocity values. This shows that the shear field is accurately recovered. The shear velocity values are constant across frequencies showing the non-viscous behavior of the agar gel. The full field presents a bias, especially at low frequencies, when the vibration power is maximum.

There are interesting observations to make regarding compression. Velocity measured from the compression field is higher than the velocity measured from the divergence field. We hypothesize that it is due to the presence of the P-I wave in the compression field. In the divergence field, the influence of the P-I is negligible, revealing the P-II wave. The compression and divergence velocities are higher than the shear velocity. This difference decreases when frequency increases. So far, the best explanation for the presence of the P-II wave in this agar gel is poroelasticity. If water can flow inside the gel matrix, it can induce a slow compression wave. Some studies suggest that a P-II wave could theoretically propagate in gels [START_REF] Duda | Fluid flow and interface motion in gels: A finitestrain theory and its application to a channel flow problem[END_REF]. Besides, agar gels are in fact known to be porous, and their porosity is being studied [START_REF] Johnson | Hydraulic permeability of agarose gels[END_REF]- [START_REF] Johnson | A limitation of the Biot-Gardner theory of extensional waves in fluid-saturated porous cylinders[END_REF]. An example of the microscopic structure of the agar gel is given on Figure 6-5. The porous structure is clearly visible and depends on the agar concentration. Tristan DERUELLE 99

FOAM PHANTOM EXPERIMENT

The next experiment is intended to compare our results with those of Aichele [START_REF] Aichele | Capturing the Shear and Secondary Compression Waves: High-Frame-Rate Ultrasound Imaging in Saturated Foams[END_REF]. They observed the poroelastic P-II wave in a porous phantom. To do so, they used a transient vibration and an ultrasound-based imaging technique. In addition, they also measured the dispersion (i.e. velocity vs. frequency) of the shear and the P-II wave in that material.

In our experiment, we used the exact same porous phantom as Aichele [START_REF] Aichele | Capturing the Shear and Secondary Compression Waves: High-Frame-Rate Ultrasound Imaging in Saturated Foams[END_REF] in their publication. It consisted a Basotect ® melamine resin foam immersed and saturated in water. The foam dimensions were 30x19x18 cm 3 .

Acquisitions were performed at the CERMEP, on a 1. Regarding the MRE acquisition, the same speaker setup was used (10 cm diameter vibrator). The foam was placed inside a plastic box and filled and saturated with water.

Repetitive manual compression was performed to drain bubbles out of the foam as much as possible. The T2-weighted image (Figure 6-6A) shows that some bubbles remained. Pictures of the setup are provided in Figure 6-7. The MRE acquisition was a multi-slice SE-EPI with 2.5 mm isotropic voxels. The field of view was 32 x 32 x 14.75 cm 3 with an acquisition matrix of 128x128x59. Acquisitions were performed for 18 time offsets per cycle, at excitation frequencies between 39 and 132 Hz. Echo time was 80ms at all frequencies. Repetition time was 7970ms at all frequencies, except at 39 Hz (7980ms) and 100 Hz (7950ms). All three spatial components of the displacement field were acquired. Data were first unwrapped with 3D quality guided algorithm, and filtered at the excitation frequency. The resulting field is referred to as 'Full Field'. The full field was separated into the shear and compressional fields. Corresponding curl and divergence fields were also computed. Velocities were reconstructed with 3D LFE.

Figure 6-8 presents the full, shear, and compression displacement fields (100 Hz), along with the amplitude of their respective curl and divergence. It is interesting to see how different the full and shear fields are. The small wavelengths stand out in the shear field.

Here also, there is a good visual match between the curl of the full and shear fields, and between the divergence of the full and compression fields. However, when looking at the compression field, it seems that the P-I wave pattern encompasses the totality of the field of view (especially visible on the Y component of the Compression field). Because of the predominance of the P-I pattern, it is difficult to visualize an eventual slow compression wave. In order to investigate the eventual presence of a P-II compression wave, the divergence is computed, and then averaged over a 7x7x7 neighborhood. This kernel is equivalent to a cube of side 1.75 cm, which is still small compared to the wavelengths. The resulting divergence field is shown in Figure 6-9 for a few time points. The zero crossing of the divergence is propagating across the field of view, and is highlighted in the figure (dashed line). The divergence of the P-I wave is expected to be negligible, and not visible. Moreover, in a P-I wave, all voxels within the field of view are expected to be more or less in phase with each other, and the zero crossing should therefore be stationary. Consequently, the presence and the propagation of the zero crossing may be representative of a slow compression wave front. Note also that the associated wavelength appears larger than that in the shear field. Overall, these observations are consistent with the presence of a slow compression wave Wave velocity was assessed for the different fields. In order to limit noise effect, a mean filter of 7x7x7 pixels was applied to the curl and divergence fields. The full, shear, and compression fields were not filtered. The corresponding velocity maps are shown Figure 6-10 for the different frequencies. Some stripes artifacts are visible, which may arise from small errors in the unwrap procedure, or from the multicoil recombination. These stripes disappeared when switching the multicoil reconstruction mode from "adaptative combination" to "sum of squares". But then the new MRE acquisition suffered from an even more severe problem. Phase vortices appeared in the images, inducing large void zones in the velocity reconstruction (see . Therefore, we went back to the original "adaptative combination" multicoil reconstruction method.

Velocities from the shear and curl field are identical. So are the velocities from the compression and divergence fields (apart from very low frequencies). Velocity from the full field is in-between the shear and compression velocities. It is important to notice that compression velocity is systematically greater than the shear velocity. It is also visible that the wave penetration depth decreases when frequency increases. At 132 Hz, penetration depth is roughly halved compared to lower frequency. At high frequency, e.g. 132 Hz, the region over which the wave may be visible decreases. Only a portion of the wavelength may thus be detectable, leading to a biased measurement of the velocity. The region where the velocity could not be reconstructed due to the phase vortex is circled by a dashed line.

Velocity values measured in the experiment are reported in Table 6-2 and Figure 5-12 along with the ones measured by Aichele by ultrasound [START_REF] Aichele | Capturing the Shear and Secondary Compression Waves: High-Frame-Rate Ultrasound Imaging in Saturated Foams[END_REF]. Regarding the shear velocity, up to 100 Hz, there is a good agreement between the two method. Beyond 100 Hz, wave penetration is not sufficient to correctly recover the velocity. Regarding compression velocity, the same behavior is observed for frequency beyond 100 Hz. For lower frequencies, there is also a good agreement, even though our values are a bit higher than those from Aichele. However, on average, the ratio of compression velocity over shear velocity is very close to that found by Aichele (1.4 in our experiment vs roughly 1.3 in their article). In this experiment, a porous media phantom, with known values of shear and slow compression velocities was used. Shear and compression fields were computed, and velocity maps reconstructed. The measured shear and compression velocities are in good agreement with the previous measurements performed with a different technique. This confirms that the slow poroelastic compression wave can be detected using MRE, and corresponding velocity can be assessed.

IN VIVO OBSERVATION

Human body is composed of 50-60 % of water distributed in intracellular fluid (2/3) and extracellular fluid (1/3) [START_REF] Henderson | Organization and composition of body fluids[END_REF]. Water content varies between different organs and tissues, and water transport can occur at the membrane scale, or at the capillaries/vessels scale.

While organs are thus known to be porous, the dynamic of propagation of poroelastic wave in soft tissues is still not well studied.

In this paragraph, we look at MRE images of kidney graft. The raw wave images stem from a previous study conducted in 2012 on renal allograft elastography (ClinicalTrials.gov Identifier NCT01234363). The study was approved by the French National Authorities (AFFSAPS identifier 2010-A00130-39) and the Ethics Committee, and all patients gave informed consent. Among all participants, three patients were randomly selected for this analysis. The same excitation setup as for the agar or foam was used. The study was conducted on a 1.5T GE MRI scanner. The vibrator was positioned on the skin right above the graft. The propagation was acquired with 6 MEG directions (+x, -x, +y, -y, +z, -z). 3 with an acquisition matrix of 64x64 and a slice spacing of 11mm (voxels 5.63x5.63x11mm 3 ). However, the output images were 256x256, with a subsequent voxel size of 1.4x1.4x11mm 3 . To avoid issues from heavily anisotropic voxels, images were resized back to 5.6mm isotropic voxels (64x64x22 matrices). A picture of the experimental setup is provided Figure 6-13A. During surgery, the renal grafts are located in the abdomen, as can be seen in the transverse view shown in in Figure 6-13B. The drum-like vibrator (arrows) was positioned anteriorly, facing the graft.

Tristan DERUELLE 105 Phase differences between MEG of opposite polarities were first computed, and then unwrapped with a 2D quality guided algorithm. Data were then filtered at the corresponding excitation frequency. Shear and compression displacement fields were computed, and velocity maps were reconstructed with 3D LFE.

Figure 6-14 presents the full, shear, and compression displacement fields (72 Hz), along with the amplitude of their respective curl and divergence. Here also, there is a good visual match between the curl of the full and shear fields, and between the divergence of the full and compression fields. However, the curl of the compression fields is not zero very close to the vibrator. We hypothesize that it is an error due to large vibration amplitude that could not be unwrapped properly in that specific region. Here, due to the complex geometry of the human body, it is much harder to visualize clear propagating wave fronts, especially after a few centimeters when the amplitude has been attenuated. It is however still possible to reconstruct velocity images, as shown in Figure 6-15 and Figure 6-16.

The position of the graft is given in Figure 6-13A. Shear and compression velocities increase with frequency (Figure 6-15). Compression velocity is again higher than the shear velocity. It is to be noted that not the whole graft is reconstructed. The region very close to the vibrator could not be recovered due to the unwrapping issues. It is also interesting to see that other organs could also be reconstructed. In those organs, further from the vibrator, the compression velocity is also higher than the shear velocity.

The same behavior is observed in the velocity maps reconstructed from the curl and divergence fields (Figure 678910111213141516). However, the reconstructed values are lower than with the shear and compression fields. This may arise from noise, but also from errors in the derivatives, especially at high frequency, when the number of points per wavelength decreases.

Tristan DERUELLE 107 Values are reported for the three patients, at the different frequencies. There is a good agreement between the different subjects. Only subject #1 presents lower shear velocity at high frequency. Nonetheless, compression velocity is higher than shear velocity. The shear and compressional fields were successfully split. Compression velocity could be measured in the kidney graft and is indeed higher than the shear velocity. It might however be overestimated because of the proximity of the vibrator, and thus by P-I wave.

Other organs, further from the vibrator are also visible. In these organs, the compression velocity also higher than the shear velocity, while the influence of the P-I wave is negligible.

This has several clinical implications. The first one is that a P-II wave can be measured.

The second is that measuring the P-II (and S) wave velocity can bring new information to the radiologist. Indeed, the P-II velocity is related to many porosity parameters. A change in the velocity can indicate a change in some of these parameters, and consequently a change in the state of the porosity of the organ. Some of these parameters are the ratio between solid and fluid volumes, tortuosity (if fluid ducts are rather straight lines or complex paths), permeability, densities…

SPECIAL CASE OF THE SILICONE PHANTOM

The last phantom used consisted in a silicone gel. Silicones are known to be very elastic, and stable with time and temperature. They consist in silicon-oxygen polymers of varying length and with varying organic group attached on the Si. The backbone chain can either be linear or cyclic. A very wide variety of silicones can thus be thought up. To avoid creep under traction, the different chains are cross-linked, generating a 3D molecular mesh. This can be done by the addition of a catalyst. The advantage of using such silicones is that it is elastic, visible under MRI, and free of water or any fluid flow. It is thus a good candidate for a purely elastic (or viscoelastic) phantom. No poroelastic waves wave should be visible.

The silicone used was a mixture of Bluesil ESA 6018 A&B (Elkem, Bluestar) with 450g of each component, cured at room temperature. The resulting gel is transparent, homogeneous, with no apparent bubbles as can be seen in Figure 6-18. The same MRE speaker setup was used (10cm vibrator). The phantom was imaged with a SE-EPI MRE sequence. The acquisition parameters were a field of view of 32x32x11.25 cm 3 with an acquisition matrix of 128x128x45 (2.5 mm isotropic voxels). Displacements were recorded for excitation frequencies between 39 Hz and 132 Hz. Motion was acquired in all three directions of space and captured for 18 different time offsets per cycle. Repetition time was 7970ms (7980 for the excitation at 39 Hz, and 7950 at 100 Hz), and echo time was 80ms.

The displacement being encoded in the phase, the first step was to unwrap the data. This was done with the use of the 3D quality guided algorithm. Then, the data were filtered at the excitation frequency. The shear field and the compressional field were calculated, along with their divergence and curl. The different velocity maps were reconstructed with 3D LFE. However, after smoothing the divergence field with a mean filter, a small propagating wave is brought out of the noise. Figure 6-21 presents the divergence field for a couple of time points. The origin of this wave is unknown yet. A mean filter (3x3x3 pixels) was applied to the curl and divergence fields. The resulting velocity maps were computed and are shown Figure 6-22. The shear and curl velocities are similar. The curl velocity is slightly lower, due to the noise still present in the image. Shear velocity increases with frequency, indicating that the silicone is viscoelastic. Velocity from the compression field is completely dominated by the P-I pattern and is thus not reliable. When estimating the velocity from the divergence field, the behavior is different from the agar. Here, both the shear and the compression velocities seem to increase with frequency. It should be noted that noise still impacts the quality of the reconstructions. This experiment used a silicone phantom. With the intrinsic structure of the silicone, and its very low ADC (1000x lower than in the agar or foam), it appears reasonable to consider it as non-porous. Compared to the agar gel, no slow waves are visible in the compression displacement field. However, the divergence presents a slow wave. Its origin is yet unknown and will require further studies.

CONCLUSION

In this chapter, we investigated the presence of a slow compression wave in various media. The compression field was computed using the field separation algorithm. Such compression field contains influence of both fast and slow compression waves.

Measuring the slow compression wave speed might thus be biased by the presence of the P-I wave. As paradoxical as it can sound, using the divergence field can help minimizing the influence of the P-I wave. The slow P-II wave was observed in an agar gel, which might be considered porous. It was also observed in a melamine foam. Measured velocities are in agreement with previous measurements performed with ultrasound. Difference between shear and compression velocities were measured in vivo in the kidney graft. Finally, a slow divergence wave was observed in a silicone phantom. From our point of view this phantom is not porous, suggesting other unknown phenomena might be at play.

Conclusions

Tristan DERUELLE 113 CHAPTER 7 Conclusions

7

The prostate is a gland that is part of the male reproductive system. When growing old, chances are that a cancer will appear. This is one of the most prevalent cancer in men.

Early detection is preferable to follow its evolution and pick the appropriate treatment option. Given its deep location in the body, and the surrounding protective tissues, access to the prostate is limited. One of the first screening tests is digital rectal examination. During this examination, only the superficial part is accessible, and if it appears relatively stiff, suspicion is raised and additional tests may be required. In the cancer detection pathway, multiparametric MRI is recommended prior biopsies. MRI can easily image the whole prostate, and spot suspicious areas. However, the analysis of MR images is quite difficult and can bring many false positives. Biopsies are required to reveal the ground truth. Nonetheless, in addition or in place of blind biopsies, guided biopsies can be performed based on the MR results. Improving prostate imaging could help detect and characterize prostate cancers while limiting the number of biopsies.

Elastography is a technique that can virtually palpate an organ and draw a stiffness map (shear modulus or shear velocity). In line with digital rectal examination, prostate elastography could provide a quantitative stiffness map of the whole prostate. In the literature, it has been shown that prostate elastography can provide significant piece of information for cancer detection. While it could be performed using ultrasound and an endorectal transducer, the operation is quite invasive and limited in imaging depth. Given that the patient undergoes MRI as part of the standard care, an MR elastography sequence could easily be added to the protocol.

MR elastography is based on the generation of waves in the body, imaging these waves by MRI, and analyzing the wave field to evaluate mechanical properties. One major challenge is to be able to generate waves in the prostate that can be detected by MRI. The second challenge is that the wavefield presents contributions from shear and compression phenomena. While compression wavelength in elastic media is too large to be measured, it can still bias the estimation of the shear modulus. Compression contribution should thus be removed without degrading the shear modulus estimation.

In this work, we designed a new prototype of medical device to induce mechanical waves for prostate Magnetic Resonance Elastography. It was designed to overcome the difficulties to vibrate the prostate (stiff, deep organ protected by damping tissues) while keeping patient's comfort and safety first. To do so, a non-invasive transperineal approach was chosen. It was designed while having in mind the European regulation. Attention was put into choosing the appropriate materials for MR compatibility. Vibrations are produced by a stack of piezoelectric elements in order to have a precise control of the vibrations and a very fast response time. Additional constraints in vibration amplitude, operating frequency, and applied load were considered. This prototype was tested on a phantom. Resulting images did not present artifacts, and the phantom stiffness geometry could be reconstructed. We are currently in the process of seeking official regulatory approval to test the prototype device in healthy volunteers and in patients.

While developing the device, the issue of compression has been addressed. A new algorithm was developed to separate the measured displacement field into the shear displacement field on one hand, and the compressional displacement field on the other hand. This algorithm was tested on simulations and experimentally on a phantom. Results show that with the new the shear field, the shear velocity can be accurately reconstructed, without artifacts from the compression. The algorithm is resistant to noise, which can be valuable for prostate imaging. Indeed, due to damping, the amplitude of the shear waves might be decreased, resulting in a loss of SNR. In elastic media, the compression field cannot provide further information. The compression wavelength is too large to be recovered. However, surprisingly, looking at the compression field in the phantom, small waves were present. This is all the more interesting that these waves are also present in the divergence field, suggesting the presence of a slow compression wave.

To our knowledge, there are not many theoretical explanations for the presence of a second and slow compression wave. The best explanation so far stems from the theory of poroelasticity. In poroelastic media, three waves can travel: the common shear and fast compression waves, but also a slow compression wave resulting from interactions between the fluid and solid matrix of the media. The velocity of this slow compression wave is higher than that of shear wave, but in same order of magnitude. Using the previous algorithm and the divergence, slow compression waves could be observed, and compression velocity estimated in different phantoms and in vivo. This shows that magnetic resonance compression elastography is feasible.

The presence of slow compression waves gives a new impetus to compression elastography. We cannot rule out other phenomena, but in the case of poroelasticity, the wave velocity is related to several parameters such as porosity, tortuosity or permeability just to name a few. These are as many new biomarkers that are available. The proposed field separation algorithm led us to properly analyze the shear wave without the influence of the fast compression wave. Thanks to the compression field, we realized that a slow compression wave could also be present. However, analysis of its velocity suffers from the presence of the fast compression wave. The use of the divergence could help mitigate its effect. The present work used MRI to acquire displacement fields. Other sequences and other modalities could be used to derive more intrinsic porosity parameters. Future studies should show if change in porosity parameters are associated with biological

Once the N directions are determined, the radial filters can be built using a similar expression as in equation (A-6). is the angle relative to the i th direction and can be computed as acos . can be defined as the angle between the i th direction and one of its closest neighbors. It should be noted that because of the imperfect distribution, and because an ellipsoidal gaussian is used, there may be some overlap between some of the gaussians. The sum of all filters would thus not be perfectly flat. Further work should focus on defining custom shaped gaussian based on the group of nearest neighbors.

Nonetheless, the method can still be used as shown on Figure A-6. Compared to using the common orthogonal directions, using 20 directional filters distributed with the Fibonacci lattice produces correct values. Artefacts present with the orthogonal directions are reduced. Impact on resolution, and true benefit on real data remains to be investigated. Because more directions are used, computed time is also increased.

A Extension of multi-directional filtering in 3D

Tristan DERUELLE 133

B

MR images are very sensitive to materials present at proximity. Juste by their presence, some objects can introduce subsequent artifacts in the images. While developing our device, we were led to test the compatibility of different materials. The present chapter reports the pieces that were tried and their impact on image quality. To do that, a CIRS 049 phantom was imaged with a T2w sequence at the CERMEP on a Siemens 1.5 T scanner.

Results are reported in Table B-1. All images are presented with the same dynamic range.

The results are interesting. Copper clearly artefatcs the image. So is aluminum 6063. This alloy is widely used. Depending on the required piece, it can be difficult to find a different alloy. Other aluminum alloys tested produce much less artifacts. Brass is also a good candidate. Here also the algorithm works well, and the results are very close to the bilayer simulation (see paragraph 5.2.2.2).

Figure 1 - 1 :

 11 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 1 - 2 :

 12 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 2 - 2 :

 22 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 2 - 7 :

 27 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 2 - 8 :

 28 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 2 - 9 :

 29 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 3 - 1 :

 31 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 3 - 2 :

 32 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 4 - 5 :

 45 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 4 - 9 :

 49 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 4 - 11 :

 411 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 4 - 12 :

 412 Figure 1-1: Anatomy plate (figure from [3]) .............................................................................Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide [5] ........................................................................................................................................................ Figure 1-4: Gleason Grade (figure from [11], [12]) .................................................................. Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress .................................................... Figure 2-2: Scheme of uniaxial stress .......................................................................................Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34]) ...................................................................................................Figure 2-4: Compressional and Shear wave behavior (figure from [35]) ............................ Figure 2-5: Steps of Magnetic Resonance Elastography ........................................................ Figure 2-6: Scheme of phase accumulation with MEG ...........................................................Figure 2-7: Pixel example for quality guided unwrap .............................................................Figure 2-8: Example of unwrap ..................................................................................................Figure 2-9: Example of unwrap #2 ............................................................................................Figure 2-10: Scheme of LFE algorithm processing .................................................................. Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76]) .......................................................................................................................... Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). ........ Figure 3-3: Elasticity in the human prostate (figure from [86]) ............................................. Figure 4-1: Picture of an APA ...................................................................................................... Figure 4-2: Frequency response for the APA for different weights ...................................... Figure 4-3: Prostate MRE device ................................................................................................ Figure 4-4: Scheme of the vibrating part .................................................................................. Figure 4-5: MRE setup .................................................................................................................Figure 4-6: MR output signals ....................................................................................................Figure 4-7: Picture of sine signal insulator ...............................................................................Figure 4-8: Picture of the TTL insulator along with electronic diagram ...............................Figure 4-9: CIRS representation .................................................................................................Figure 4-10: Setup of the vibration device with a phantom ..................................................Figure 4-11: Phase fields in the CIRS phantom .......................................................................Figure 4-12: Frequency content of the displacements in the CIRS phantom ..................... Figure 4-13: Displacement and curl fields in the CIRS phantom ..........................................Figure 4-14: Wave propagation in the CIRS phantom ............................................................Figure 4-15: Displacement field in the CIRS phantom, Coronal view ...................................Figure 4-16: Reconstructed velocity maps ...............................................................................

Figure 5 - 1 :

 51 Figure 5-1: Example of drum-like vibrator ............................................................................... Figure 5-2: Example of vibrating Membrane ........................................................................... Figure 5-3: Homogeneous Simulation -Full displacement, shear and compressional fields, with associated norms of curl and divergence ....................................................................... Figure 5-4: Homogeneous Simulation -Velocity maps .......................................................... Figure 5-5: Bilayer Simulation Scheme ..................................................................................... Figure 5-6: Bilayer Simulation -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) ......Figure 5-7: Displacement amplitudes of a compressional wave generated by a point source in an agar-gelatin phantom (figure from [129]) .........................................................Figure 5-8: Bilayer Simulation -Velocity maps ........................................................................ Figure 5-9: Median bilayer profile and fitted sigmoid ............................................................ Figure 5-10: Example of an SNR map ....................................................................................... Figure 5-11: Full field, Curl Field, and Shear Field at the different levels of noise ............. Figure 5-12: Reconstructed velocities [m/s] at the different levels of noise ....................... Figure 5-13: Boxplots of the reconstructed velocities at the different levels of noise ...... Figure 5-14: Comparison of reconstructed velocity from the smooth curl field and shear field ................................................................................................................................................ Figure 5-15: T2w, elasticity, and ADC maps in the CIRS phantom ........................................Figure 5-16: Schematic representation of the speaker setup ...............................................Figure 5-17: Full field, shear field, and compression fields in the CIRS phantom ..............Figure 5-18: Velocity reconstruction in CIRS phantom [m/s] ................................................Figure 5-19: Divergence field in the CIRS phantom ................................................................ Figure 5-20: Divergence field computed with two different algorithms .............................. Figure 5-21: Divergence field computed from data from two different scanners ............. Figure 5-22: Sagittal View -Dejittering ..................................................................................... Figure 5-23: Divergence field computed with and without jitter .......................................... Figure 6-1: T2w and ADC in the agar gel .................................................................................. Figure 6-2: Displacement fields in the agar gel ....................................................................... Figure 6-3: Divergence field in the agar gel ............................................................................. Figure 6-4: Velocity maps in the agar gel for different frequencies .....................................Figure 6-5: Cryo-SEM micrographs of different agarose hydrogels (figure from [158]) ....Figure 6-6: T2w and ADC in the foam phantom ......................................................................Figure 6-7: MRE acquisition setup for the foam phantom ..................................................Figure 6-8: Displacement fields in the foam phantom .........................................................Figure 6-9: Divergence field in the foam phantom ...............................................................Figure 6-10: Velocity maps in the foam phantom for different frequencies ....................Figure 6-11: Example of void in the velocity reconstruction induced by phase vortex ...Figure 6-12: Shear and Compression velocities measured in the foam ............................ Figure 6-13: In vivo experiment ...............................................................................................Figure 6-14: Displacement fields in vivo .................................................................................Figure 6-15: Shear and compression Velocity maps in vivo ................................................Figure 6-16: Curl and divergence Velocity maps in vivo .......................................................

Figure 6 - 9 :

 69 Figure 5-1: Example of drum-like vibrator ............................................................................... Figure 5-2: Example of vibrating Membrane ........................................................................... Figure 5-3: Homogeneous Simulation -Full displacement, shear and compressional fields, with associated norms of curl and divergence ....................................................................... Figure 5-4: Homogeneous Simulation -Velocity maps .......................................................... Figure 5-5: Bilayer Simulation Scheme ..................................................................................... Figure 5-6: Bilayer Simulation -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) ......Figure 5-7: Displacement amplitudes of a compressional wave generated by a point source in an agar-gelatin phantom (figure from [129]) .........................................................Figure 5-8: Bilayer Simulation -Velocity maps ........................................................................ Figure 5-9: Median bilayer profile and fitted sigmoid ............................................................ Figure 5-10: Example of an SNR map ....................................................................................... Figure 5-11: Full field, Curl Field, and Shear Field at the different levels of noise ............. Figure 5-12: Reconstructed velocities [m/s] at the different levels of noise ....................... Figure 5-13: Boxplots of the reconstructed velocities at the different levels of noise ...... Figure 5-14: Comparison of reconstructed velocity from the smooth curl field and shear field ................................................................................................................................................ Figure 5-15: T2w, elasticity, and ADC maps in the CIRS phantom ........................................Figure 5-16: Schematic representation of the speaker setup ...............................................Figure 5-17: Full field, shear field, and compression fields in the CIRS phantom ..............Figure 5-18: Velocity reconstruction in CIRS phantom [m/s] ................................................Figure 5-19: Divergence field in the CIRS phantom ................................................................ Figure 5-20: Divergence field computed with two different algorithms .............................. Figure 5-21: Divergence field computed from data from two different scanners ............. Figure 5-22: Sagittal View -Dejittering ..................................................................................... Figure 5-23: Divergence field computed with and without jitter .......................................... Figure 6-1: T2w and ADC in the agar gel .................................................................................. Figure 6-2: Displacement fields in the agar gel ....................................................................... Figure 6-3: Divergence field in the agar gel ............................................................................. Figure 6-4: Velocity maps in the agar gel for different frequencies .....................................Figure 6-5: Cryo-SEM micrographs of different agarose hydrogels (figure from [158]) ....Figure 6-6: T2w and ADC in the foam phantom ......................................................................Figure 6-7: MRE acquisition setup for the foam phantom ..................................................Figure 6-8: Displacement fields in the foam phantom .........................................................Figure 6-9: Divergence field in the foam phantom ...............................................................Figure 6-10: Velocity maps in the foam phantom for different frequencies ....................Figure 6-11: Example of void in the velocity reconstruction induced by phase vortex ...Figure 6-12: Shear and Compression velocities measured in the foam ............................ Figure 6-13: In vivo experiment ...............................................................................................Figure 6-14: Displacement fields in vivo .................................................................................Figure 6-15: Shear and compression Velocity maps in vivo ................................................Figure 6-16: Curl and divergence Velocity maps in vivo .......................................................

Figure 6 - 11 :

 611 Figure 5-1: Example of drum-like vibrator ............................................................................... Figure 5-2: Example of vibrating Membrane ........................................................................... Figure 5-3: Homogeneous Simulation -Full displacement, shear and compressional fields, with associated norms of curl and divergence ....................................................................... Figure 5-4: Homogeneous Simulation -Velocity maps .......................................................... Figure 5-5: Bilayer Simulation Scheme ..................................................................................... Figure 5-6: Bilayer Simulation -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) ......Figure 5-7: Displacement amplitudes of a compressional wave generated by a point source in an agar-gelatin phantom (figure from [129]) .........................................................Figure 5-8: Bilayer Simulation -Velocity maps ........................................................................ Figure 5-9: Median bilayer profile and fitted sigmoid ............................................................ Figure 5-10: Example of an SNR map ....................................................................................... Figure 5-11: Full field, Curl Field, and Shear Field at the different levels of noise ............. Figure 5-12: Reconstructed velocities [m/s] at the different levels of noise ....................... Figure 5-13: Boxplots of the reconstructed velocities at the different levels of noise ...... Figure 5-14: Comparison of reconstructed velocity from the smooth curl field and shear field ................................................................................................................................................ Figure 5-15: T2w, elasticity, and ADC maps in the CIRS phantom ........................................Figure 5-16: Schematic representation of the speaker setup ...............................................Figure 5-17: Full field, shear field, and compression fields in the CIRS phantom ..............Figure 5-18: Velocity reconstruction in CIRS phantom [m/s] ................................................Figure 5-19: Divergence field in the CIRS phantom ................................................................ Figure 5-20: Divergence field computed with two different algorithms .............................. Figure 5-21: Divergence field computed from data from two different scanners ............. Figure 5-22: Sagittal View -Dejittering ..................................................................................... Figure 5-23: Divergence field computed with and without jitter .......................................... Figure 6-1: T2w and ADC in the agar gel .................................................................................. Figure 6-2: Displacement fields in the agar gel ....................................................................... Figure 6-3: Divergence field in the agar gel ............................................................................. Figure 6-4: Velocity maps in the agar gel for different frequencies .....................................Figure 6-5: Cryo-SEM micrographs of different agarose hydrogels (figure from [158]) ....Figure 6-6: T2w and ADC in the foam phantom ......................................................................Figure 6-7: MRE acquisition setup for the foam phantom ..................................................Figure 6-8: Displacement fields in the foam phantom .........................................................Figure 6-9: Divergence field in the foam phantom ...............................................................Figure 6-10: Velocity maps in the foam phantom for different frequencies ....................Figure 6-11: Example of void in the velocity reconstruction induced by phase vortex ...Figure 6-12: Shear and Compression velocities measured in the foam ............................ Figure 6-13: In vivo experiment ...............................................................................................Figure 6-14: Displacement fields in vivo .................................................................................Figure 6-15: Shear and compression Velocity maps in vivo ................................................Figure 6-16: Curl and divergence Velocity maps in vivo .......................................................

Figure 6 - 12 :

 612 Figure 5-1: Example of drum-like vibrator ............................................................................... Figure 5-2: Example of vibrating Membrane ........................................................................... Figure 5-3: Homogeneous Simulation -Full displacement, shear and compressional fields, with associated norms of curl and divergence ....................................................................... Figure 5-4: Homogeneous Simulation -Velocity maps .......................................................... Figure 5-5: Bilayer Simulation Scheme ..................................................................................... Figure 5-6: Bilayer Simulation -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) ......Figure 5-7: Displacement amplitudes of a compressional wave generated by a point source in an agar-gelatin phantom (figure from [129]) .........................................................Figure 5-8: Bilayer Simulation -Velocity maps ........................................................................ Figure 5-9: Median bilayer profile and fitted sigmoid ............................................................ Figure 5-10: Example of an SNR map ....................................................................................... Figure 5-11: Full field, Curl Field, and Shear Field at the different levels of noise ............. Figure 5-12: Reconstructed velocities [m/s] at the different levels of noise ....................... Figure 5-13: Boxplots of the reconstructed velocities at the different levels of noise ...... Figure 5-14: Comparison of reconstructed velocity from the smooth curl field and shear field ................................................................................................................................................ Figure 5-15: T2w, elasticity, and ADC maps in the CIRS phantom ........................................Figure 5-16: Schematic representation of the speaker setup ...............................................Figure 5-17: Full field, shear field, and compression fields in the CIRS phantom ..............Figure 5-18: Velocity reconstruction in CIRS phantom [m/s] ................................................Figure 5-19: Divergence field in the CIRS phantom ................................................................ Figure 5-20: Divergence field computed with two different algorithms .............................. Figure 5-21: Divergence field computed from data from two different scanners ............. Figure 5-22: Sagittal View -Dejittering ..................................................................................... Figure 5-23: Divergence field computed with and without jitter .......................................... Figure 6-1: T2w and ADC in the agar gel .................................................................................. Figure 6-2: Displacement fields in the agar gel ....................................................................... Figure 6-3: Divergence field in the agar gel ............................................................................. Figure 6-4: Velocity maps in the agar gel for different frequencies .....................................Figure 6-5: Cryo-SEM micrographs of different agarose hydrogels (figure from [158]) ....Figure 6-6: T2w and ADC in the foam phantom ......................................................................Figure 6-7: MRE acquisition setup for the foam phantom ..................................................Figure 6-8: Displacement fields in the foam phantom .........................................................Figure 6-9: Divergence field in the foam phantom ...............................................................Figure 6-10: Velocity maps in the foam phantom for different frequencies ....................Figure 6-11: Example of void in the velocity reconstruction induced by phase vortex ...Figure 6-12: Shear and Compression velocities measured in the foam ............................ Figure 6-13: In vivo experiment ...............................................................................................Figure 6-14: Displacement fields in vivo .................................................................................Figure 6-15: Shear and compression Velocity maps in vivo ................................................Figure 6-16: Curl and divergence Velocity maps in vivo .......................................................

Figure 6 - 17 :

 617 Figure 6-17: Wave velocities in the kidney graft for the different patients and frequencies ...................................................................................................................................................... Figure 6-18: Picture of silicon phantom ................................................................................. Figure 6-19: T2w and ADC in the silicone phantom .............................................................. Figure 6-20: Displacement fields in the silicone phantom .................................................. Figure 6-21: Divergence along time in the silicone phantom ..............................................Figure 6-22: Velocity maps in the silicone phantom .............................................................

Figure A- 1 :

 1 Figure A-1: 2D Orthogonal filters ............................................................................................. Figure A-2: Representation of spherical coordinates ........................................................... Figure A-3: 12 directional filters in 2D .................................................................................... Figure A-4: Distribution of points over a sphere ................................................................... Figure A-5: Example of 3D directional filter ........................................................................... Figure A-6: Reconstruction with more filters from synthetic data ..................................... Figure D-1: Homogeneous Simulation, cP = 1500 m/s -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) .................................................................................................................... Figure D-2: Homogeneous Simulation, cP = 1500 m/s -Velocity maps .............................

Figure 1 - 1 :

 11 Figure 1-1: Anatomy plate (figure from [3])

Figure 1 - 2 :

 12 Figure 1-2: Scheme of a normal prostate with different zones (figure adapted from [4]) PZ = Peripheral Zone; CZ = Central Zone; TZ = Transition Zone; AFS = Anterior Fibromuscular Stroma.

Figure 1 - 3 :

 13 Figure 1-3: Number of new cases and deaths of different cancers in men, worldwide[5]

Figure 1 - 4 :

 14 Figure 1-4: Gleason Grade (figure from[START_REF] Gleason | Classification of prostatic carcinomas[END_REF],[START_REF] Salomon | Le score de Gleason pour les nuls[END_REF])Usually, multiple grades are present in a sample. The Gleason's score has thus been defined as the sum of the two most present grades. For example, if the two most present grades are 3 and 4, a score of 7 is attributed (if the grade is homogeneous e.g. 3, then a score or 2*3=6 is attributed). The International Society of Urological Pathology (ISUP) score. The ISUP score was introduced in 2014 to simplify reporting. It ranges from 1 to 5, where ISUP 1 corresponds to Gleason 3+3, and ISUP 5 is the most aggressive (corresponding to Gleason 4+5 or higher).

Figure 2 - 1 :

 21 Figure 2-1: LEFT: Isotropic pressure; RIGHT: Shear stress

Figure 2 - 2 :

 22 Figure 2-2: Scheme of uniaxial stress

Figure 2 - 3 :

 23 Figure 2-3: Variation of Bulk, Shear, and Young's moduli for various materials. Units are in Pa. (Figure adapted from [34])

Figure 2 - 4 :

 24 Figure 2-4: Compressional and Shear wave behavior (figure from [35])

Figure 2 - 5 :

 25 Figure 2-5: Steps of Magnetic Resonance Elastography

Figure 2 - 6 :

 26 Figure 2-6: Scheme of phase accumulation with MEG

Figure 2 - 7 :

 27 Figure 2-7: Pixel example for quality guided unwrap H1 and H2 are the two orthogonal differentiators, and D1 and D2 the diagonal differentiators. This can be easily extended to 3D (3 orthogonals, 10 diagonals).

Figure 2 -

 2 Figure 2-8 presents the results of different algorithms on synthetic spherical wave. Figure 2-8A shows results in a noise-free scenario, and Figure 2-8B in a noisy scenario. The full dynamic phase data (blue line with * markers) is wrapped to the [line). The wrapped image is unwrapped using the quality guided algorithm (yellow solid line). If the seed pixels belong to a wrapped region, the unwrapped image can present a general offset (Figure2-8B). This offset may be removed with the Fourier transform, or by subtracting the average value (yellow dashed line). The image unwrapped with the Laplacian algorithm (green line) is also very close to the input data. Some deviation can still be present if the data is insufficiently sampled. If noise is added (Figure2-8B), both algorithms successfully unwrap the data. However, the Laplacian method may not recover the full amplitude.

Figure 2 - 8 :

 28 Figure 2-8: Example of unwrap Images of (un)wrapped synthetic data with corresponding profile. TOP: From left to right: Original full dynamic phase, wrapped phase to [-; ], unwrapped phase with quality guided (QG) algorithm, unwrapped phase with quality guided algorithm minus average of the resulting unwrapped image (QG -AVG), unwrapped phase with the Laplacian algorithm using Discrete Cosine Transform (DCT). BOTTOM: Profile along the middle column of each image. A without noise, B with uniform random noise.

Figure 2 - 9 :

 29 Figure 2-9: Example of unwrap #2 Images of (un)wrapped synthetic data. From left to right: Original full dynamic phase, wrapped phase to [-π; π], unwrapped phase with quality guided algorithm, unwrapped phase with the Laplacian algorithm using Discrete Cosine Transform.

  -10.

Figure 2 -

 2 Figure 2-10: Scheme of LFE algorithm processing

Figure 3 - 1 :

 31 Figure 3-1: SWE prostate measurements in 57-year-old man with a PSA level of 6.62 ng/mL (figure from [76])The top image shows an orange-coded stiff area (solid-line ROI; mean elasticity, 54 kPa) incidentally discovered in the left base paramedian peripheral zone with no abnormal features at B-mode US (bottom image). Another ROI (dashed line; mean elasticity, 18.4 kPa) was located in the lateral peripheral zone of the same sextant. The biopsy performed in the stiffest area revealed an 11-mm Gleason 7 adenocarcinoma.

Figure 3 - 2 :

 32 Figure 3-2: Apical imaging plane of subject with prostate cancer (figure from [89]). PI-RADS 5, Gleason: 4+4=8. ADC map (a); averaged MRE magnitude image (b); and elasticity map depicting shear wave speed c (c). The red arrows point to the area of restricted diffusion in the ADC map (a) and to the area of increased stiffness in the c-map (c).

Figure 3 - 3 :CHAPTER 4

 334 Figure 3-3: Elasticity in the human prostate (figure from [86])

Figure 4 - 1 :

 41 Figure 4-1: Picture of an APA

Figure 4 - 2 :

 42 Figure 4-2: Frequency response for the APA for different weights LEFT: Frequency response for the APA400MML. RIGHT: Frequency response for the APA100MML. Each line corresponds on a weight, from 0kg (blue) to 1kg (red), every 0.1kg.

Figure 4 - 3 :

 43 Figure 4-3: Prostate MRE device

Figure 4 - 4 :

 44 Figure 4-4: Scheme of the vibrating part

Figure 4 Figure 4 - 6 :

 446 Figure 4-5: MRE setup

Figure 4 - 7 :

 47 Figure 4-7: Picture of sine signal insulator

Figure 4 - 8 :

 48 Figure 4-8: Picture of the TTL insulator along with electronic diagram

4. 3 PERFORMANCE 4 . 3 . 1

 3431 Wave acquisitionThe performances of the setup were assessed on an elastic CIRS ® 049 phantom. The experiment was conducted at HEH, on a 1.5 T GE scanner. This phantom consists in a tissue-mimicking medium with 4 spherical inclusions of different stiffnesses. An image and three orthogonal views of the T2w MRI image are presented on Figure4-9. The 4 inclusions, of different elasticity are visible. The two inclusions on the left are stiffer than the background, and the two on the right are softer.

Figure 4 - 9 :

 49 Figure 4-9: CIRS representation A: Phantom image; B: Orthogonal views of the T2w image.

Figure 4 -

 4 Figure 4-10 shows a picture of the experiment. The phantom is placed on the bed board, in contact with the device's ball. Contact is maintained using an elastic belt.A T2w image was acquired (Figure4-9B). While the device is in the vicinity of the phantom, there is no visual image degradation due to its presence.

Figure 4 - 10 :

 410 Figure 4-10: Setup of the vibration device with a phantom

Figure 4 - 11 :

 411 Figure 4-11: Phase fields in the CIRS phantom X, Y, and Z phase fields measured in the CIRS. X is from left to right, Y is top to bottom, Z is front to back. 1 st row represents the raw phase fields captured with positive polarity of the MEG. 2 nd row represents the raw phase fields captured with negative polarity of the MEG. 3 rd row represents the phase difference between positive and negative polarities. 4 th row represents the unwrapped phase difference.

Figure 4 - 12 :

 412 Figure 4-12: Frequency content of the displacements in the CIRS phantom

Figure 4 - 13 :

 413 Figure 4-13: Displacement and curl fields in the CIRS phantom TOP: Full displacement field, in micrometers. BOTTOM: Curl field, dimensionless.

Figure 4 - 14 :

 414 Figure 4-14: Wave propagation in the CIRS phantom Wave propagation in the X direction is represented for 8 times points over a vibration cycle. Displacements are in micrometers.

Figure 4 - 15 :

 415 Figure 4-15: Displacement field in the CIRS phantom, Coronal view X Displacement field from Full Field, in micrometers. Arrows indicate ±30°.

Figure 4 - 16 :

 416 Figure 4-16: Reconstructed velocity maps Velocity maps reconstructed from the Full Field (LEFT) and from the Curl Field (RIGHT). Units are in m/s. TOP row shows the axial cut as captured at the acquisition. BOTTOM row shows the corresponding coronal cut. The overestimation due to the vibrator is circled by a dashed line.

Figure 5 - 1 :

 51 Figure 5-1: Example of drum-like vibrator

Figure 5 - 2 :

 52 Figure 5-2: Example of vibrating Membrane

Figure 5 - 3 :

 53 Figure 5-3: Homogeneous Simulation -Full displacement, shear and compressional fields, with associated norms of curl and divergence The columns represent the x-component, y-component, z-component of displacement, amplitude of the curl, and of the divergence. The rows represent TOP: full displacement, CENTRAL: shear displacement, and BOTTOM: compressional displacement. The slice 33/100 is shown.

Finally

  , the velocities are reconstructed with 3D LFE. The velocity maps are presented in Figure 5-4. Velocity from the full field (A) is not uniform with an average value around Field separation Tristan DERUELLE 772.8 m/s, in-between the shear and compression velocities, and a dispersion (standard deviation) of 0.5 m/s. On the other hand, the shear velocity (B) is unbiased, centered around its expected 2.0 m/s, and more precise, with a standard deviation of 0.2 m/s. The compressional velocity (C) is also more homogeneous and centered around the 3.0 m/s. There are artefacts in the top corners mainly because there is almost no displacement in these regions.

Figure 5 - 4 :

 54 Figure 5-4: Homogeneous Simulation -Velocity maps Reconstructed velocities at slice 33/100 from LEFT: Full Field; MIDDLE: Shear Field; RIGHT: Compressional Field. Values are reported as median ± standard deviation. Units are in [m/s].

Figure 5 - 5 :

 55 Figure 5-5: Bilayer Simulation Scheme

Figure 5 - 6 :

 56 Figure 5-6: Bilayer Simulation -Full displacement, shear and compressional fields (units are in [m]), with associated norms of curl and divergence (units are dimensionless) The slice 33/86 is shown. The columns represent the x-component, y-component, z-component of displacement, amplitude of the curl, and of the divergence. The rows represent TOP: full displacement, CENTRAL: shear displacement, and BOTTOM: compressional displacement.

Figure 5 - 7 :

 57 Figure 5-7: Displacement amplitudes of a compressional wave generated by a point source in an agar-gelatin phantom (figure from [129]) LEFT: Isoamplitude RIGHT: gray level image

Figure 5 - 9 :

 59 Figure 5-9: Median bilayer profile and fitted sigmoid LEFT: median vertical profile of the reconstructed velocity. RIGHT: fitted sigmoid.

Figure 5 - 10 :

 510 Figure 5-10: Example of an SNR map LEFT: X component of the displacement vector with 30dB of noise. RIGHT: corresponding SNR map, in dB.

Figure 5 -

 5 Figure 5-11A presents the x-component (left column), the y-component (middle column), and the z-component (right column) of the full displacement field. The four rows nd 30 dB SNRpeak.Figure 5-11B and Figure 5-11C show the curl field and the computed shear field.

Figure 5 -

 5 Figure 5-11A presents the x-component (left column), the y-component (middle column), and the z-component (right column) of the full displacement field. The four rows nd 30 dB SNRpeak.Figure 5-11B and Figure 5-11C show the curl field and the computed shear field.

Figure 5 - 11 :

 511 Figure 5-11: Full field, Curl Field, and Shear Field at the different levels of noise Slice 29/86. A: Full Field; B: Curl Field; C: Shear Field. For each subplot, the columns refer to the x-, y-, and z-components, and the rows refer to the ∞, 60, 40, and 30dB (top to bottom). Full and shear fields are in µm. Curl field is dimensionless. The corresponding velocity maps are shown Figure 5-12. Without preprocessing, the two media are well reconstructed when the SNRpeak is at least 40 dB. With the curl operator,

Figure 5 - 12 :Figure 5 -

 5125 Figure 5-12: Reconstructed velocities [m/s] at the different levels of noise

Figure 5 - 13 :

 513 Figure 5-13: Boxplots of the reconstructed velocities at the different levels of noise

Figure 5 - 14 :

 514 Figure 5-14: Comparison of reconstructed velocity from the smooth curl field and shear field Prior applying the curl operator, a 7x7x7 spatial averaging was performed. Averaging was not performed for the shear field. Reconstructed velocities are in [m/s]. Slice 33/86 is presented Corresponding median profiles are plotted on the right (orange: mean filter + Curl; green: shear field).

Figure 5 -

 5 15A-B. Diffusion images with b-values of 50, 400, and 800 s/mm² were also acquired. The resulting Apparent Diffusion Coefficient (ADC) map is shown Figure 5-15C. ADC is 1450±40 μm²/s in the background, and 1400±10 μm²/s in the hard inclusion (for comparison, the ADC is around 1300 μm²/s in benign regions in the prostate peripheral zone, and often drops below 1000 μm²/s in cancerous lesions).

Figure 5 - 15 :

 515 Figure 5-15: T2w, elasticity, and ADC maps in the CIRS phantom A: T2w image. B: expected shear velocity map. C: ADC map.

Figure 5 - 16 :

 516 Figure 5-16: Schematic representation of the speaker setup

Figure 5 - 17 :

 517 Figure 5-17: Full field, shear field, and compression fields in the CIRS phantom Slice 15/21. Displacement fields are in micrometers [µm].

Figure 5 - 18 :

 518 Figure 5-18: Velocity reconstruction in CIRS phantom [m/s] Slice 15/21. Units are in [m/s].

Figure 5 - 19 :

 519 Figure 5-19: Divergence field in the CIRS phantom Divergence of the full field (LEFT column) and of the compression field (RIGHT). Both show the slow compression waves

Figure 5 - 20 :

 520 Figure 5-20: Divergence field computed with two different algorithms LEFT: Divergence computed using centered Finite Difference scheme. RIGHT: Divergence computed using Fourier-based method.

Figure 5 - 21 :

 521 Figure 5-21: Divergence field computed from data from two different scanners LEFT: Divergence computed from data acquired on a GE scanner. RIGHT: Divergence computed from data acquired from a Siemens scanner. Centered Finite Difference scheme is used for the computation of the divergence.

Figure 5 -

 5 22 presents the Sagittal view of the acquisition Field separation

Figure 5 - 22 :

 522 Figure 5-22: Sagittal View -Dejittering X component of the Full Field in sagittal view. Jitter is present (see arrow). Dejittering successfully removes the stripes artifacts.

Figure 5 - 23 :

 523 Figure 5-23: Divergence field computed with and without jitter LEFT: Divergence computed without dejittering. RIGHT: Divergence computed with dejittering. In both cases, the divergence is computed with centered Finite Difference.

Figure 6 - 1 :

 61 Figure 6-1: T2w and ADC in the agar gel A: T2w image. Arrows indicate the extremities of the vibrator for of the MRE experiment. B: ADC map. C: Histogram of the ADC values taken in the blue rectangle.

Figure 6 - 2 :

 62 Figure 6-2: Displacement fields in the agar gel TOP ROW: Full Field. MIDDLE ROW: Shear Field. BOTTOM ROW: Compression Field. The from LEFT to RIGHT: X-, Y-, and Z-components of the displacement (in µm), followed by the amplitude of the curl and of the divergence (dimensionless). Slice 24/41 is represented. Excitation frequency was 100 Hz. The in-plane phantom dimensions were 11x13 cm.

Figure 6 - 3 :

 63 Figure 6-3: Divergence field in the agar gel Divergence of the full wave field at 8 points over a vibration cycle. A mean filter is applied (3x3x3 pixels).

Figure 6 - 4 :

 64 Figure 6-4: Velocity maps in the agar gel for different frequencies Reconstructed velocity maps (in m/s) from the full field (1 st row), shear field, (2 nd row), curl field with mean filter (3 rd row), compression field (4 th row), and divergence field with mean filter (5 th row). Each column corresponds to a vibration frequency. Slice 24/41 is represented.

Figure 6 - 5 :

 65 Figure 6-5: Cryo-SEM micrographs of different agarose hydrogels (figure from [158]) (a) 0.3 wt% ; (b) 1 wt% ; (c) 2.9 wt%. Samples are imaged using a secondary electron detector (SEI mode) at an acceleration voltage of 5 kV

5 T

 5 Siemens scanner. First, a T2weigthed image, and a diffusion acquisition performed. The diffusion acquisition was performed with b-values of 50, 800, and 2000 s/mm². The resulting images are shown Figure 6-6. The Apparent Diffusion Coefficient in the foam phantom is 1976±69 μm²/s. It is similar to that in the agar gel.

Figure 6 - 6 :

 66 Figure 6-6: T2w and ADC in the foam phantom A: T2w image. Arrows indicate the extremities of the vibrator in the MRE experiment B: ADC map. C: Histogram of the ADC values taken in the blue rectangle.

Figure 6 - 7 :

 67 Figure 6-7: MRE acquisition setup for the foam phantom A: The foam is placed in a plastic case an immersed with water. The vibrator is positioned on top of the foam. B: Contact is maintained with an abdominal belt and abdominal coil.

Figure 6 - 8 :

 68 Figure 6-8: Displacement fields in the foam phantom TOP ROW: Full Field. MIDDLE ROW: Shear Field. BOTTOM ROW: Compression Field. The from LEFT to RIGHT: X-, Y-, and Z-components of the displacement (in µm), followed by the amplitude of the curl and of the divergence (dimensionless). Slice 30/59 is represented. Excitation frequency was 100 Hz. The in-plane phantom dimensions were 30x19 cm.

Figure 6 - 9 :

 69 Figure 6-9: Divergence field in the foam phantom Divergence of the full wave field. 5 out 12 points over a vibration cycle. A mean filter is applied (7x7x7 pixels). For better visualization, a part of a zero-crossing front is drawn as a dashed line. Slice 30/59 is represented. Excitation frequency was 100 Hz.

Figure 6 - 10 :

 610 Figure 6-10: Velocity maps in the foam phantom for different frequencies Reconstructed velocity maps (in m/s) from the full field (1 st row), shear field, (2 nd row), curl field (3 rd row) compression field (4 th row), and divergence field (5 th row). Each column corresponds to a vibration frequency. Slice 30/59 is represented.

Figure 6 - 11 :

 611 Figure 6-11: Example of void in the velocity reconstruction induced by phase vortex

Figure 6 - 12 :

 612 Figure 6-12: Shear and Compression velocities measured in the foam

Figure 6 - 13 :

 613 Figure 6-13: In vivo experiment A: The vibrator is positioned right above the renal graft. Contact is kept with abdominal belt. B: MR magnitude image of an axial plane. The renal graft is drawn in red. Arrows indicate the extremity of the vibrator.

Figure 6 - 14 :

 614 Figure 6-14: Displacement fields in vivo TOP ROW: Full Field. MIDDLE ROW: Shear Field. BOTTOM ROW: Compression Field. The from LEFT to RIGHT: X-, Y-, and Z-components of the displacement (in µm), followed by the amplitude of the curl and of the divergence (dimensionless). Slice 13/22 is represented. Excitation frequency was 72 Hz.

Figure 6 - 15 :

 615 Figure 6-15: Shear and compression Velocity maps in vivo Reconstructed velocity maps (in m/s) from the shear field (TOP row), and compression field (BOTTOM row). Each column corresponds to a vibration frequency. Slice 13/22 is represented.

Figure 6 - 16 :

 616 Figure 6-16: Curl and divergence Velocity maps in vivo Reconstructed velocity maps (in m/s) from the curl field (TOP row), and divergence field (BOTTOM row). Each column corresponds to a vibration frequency. Slice 13/22 is represented.

Figure 6 -

 6 Figure 6-17 presents boxplots of the velocity measurements in the whole kidney graft.Values are reported for the three patients, at the different frequencies. There is a good agreement between the different subjects. Only subject #1 presents lower shear velocity at high frequency. Nonetheless, compression velocity is higher than shear velocity.

Figure 6 - 17 :

 617 Figure 6-17: Wave velocities in the kidney graft for the different patients and frequencies Wave velocity (m/s) is plotted against frequency (Hz) for the different patients. Values are measured in the whole graft. The dot represents the mean value. LEFT: values reconstructed from the shear field; RIGHT: values reconstructed from the compressional field.

Figure 6 - 18 :

 618 Figure 6-18: Picture of silicon phantom

Figure 6 - 19 :

 619 Figure 6-19: T2w and ADC in the silicone phantom A: T2w image. Arrows indicate the extremities of the vibrator in the MRE experiment B: ADC map. C: Histogram of the ADC values taken in the blue rectangle.

Figure 6 -

 6 Figure 6-20 presents the full, shear, and compression displacement fields (100 Hz), along with the amplitude of their respective curl and divergence. Short wavelengths are visible on the shear fields, while only large and rapidly decaying patterns are visible on the compressional fields. The curl of the compressional field is almost zero indicating there is very little shear contribution left. The divergence of the compression field is very close to that of the full field. The divergence field is very low, suggesting there is very little compression.

Figure 6 - 20 :

 620 Figure 6-20: Displacement fields in the silicone phantom TOP ROW: Full Field. MIDDLE ROW: Shear Field. BOTTOM ROW: Compression Field. The from LEFT to RIGHT: X-, Y-, and Z-components of the displacement (in µm), followed by the amplitude of the curl and of the divergence (dimensionless). Slice 24/45 is represented. Excitation frequency was 100 Hz. The in-plane phantom dimensions were 11x10 cm.

Figure 6 - 21 :

 621 Figure 6-21: Divergence along time in the silicone phantom Divergence of the full wave field. 4 out 12 points over a vibration cycle. A mean filter is applied (3x3x3 pixels). Slice 24/45 is represented. Excitation frequency was 100 Hz.

Figure 6 - 22 :

 622 Figure 6-22: Velocity maps in the silicone phantom Reconstructed velocity maps (in m/s) from the full field (1 st row), shear field, (2 nd row), curl field (3 rd row) compression field (4 th row), and divergence field (5 th row). Each column corresponds to a vibration frequency. Slice 24/45 is represented.

  Figure A-5 shows an example of directional filter from a series of 20 filters.

Figure A- 5 :

 5 Figure A-5: Example of 3D directional filter

Figure A- 6 :

 6 Figure A-6: Reconstruction with more filters from synthetic data

Table B- 1 :

 1 Tested Material along with corresponding phantom image Material tested MR images No material is tested Reference image of the phantom alone

Figure D- 2 :

 2 Figure D-2: Homogeneous Simulation, cP = 1500 m/s -Velocity maps Reconstructed velocities at slice 33/86 from A: Full Field; B: Shear Field; C: Compressional Field. Values are reported as median ± standard deviation. Units are in [m/s].

Table of figures

 of Tristan DERUELLE

Table of tables

 of Tristan DERUELLE

Table of tables

 of

Table 5 -

 5 1: Bilayer simulation properties .................................................................................. Table 5-2: Measured velocities (m/s) in the CIRS phantom ................................................... Table 6-1: Poroelastic parameters ............................................................................................ Table 6-2: Measured velocity values in the foam phantom [m/s] ...................................... Table B-1: Tested Material along with corresponding phantom image ............................

	Abbreviations
	ADC	Apparent Diffusion Coefficient
	ANSM	Agence Nationale de Sécurité du Médicament et des produits de santé
	APA	Amplified Piezo Actuator
	ARF	Acoustic Radiation Force
	CAD	Computer Aided Diagnostic
	CERMEP Centre d'Etude et de Recherche Multimodal Et Pluridisciplinaire en imagerie du
		vivant
	CPP	Comité de Protection des Personnes
	DCE	Dynamic Contrast Enhanced
	DRE	Digital Rectal Examination
	EPI	Echo Planar Imaging
	FDTD	Finite Difference Time Difference
	FEM	Finite Element Method
	FOV	Field of View
	GE	General Electric
	HEH	Hôpital Edouard Heriot
	HIFU	High Intensity Focused Ultrasound
	ISUP	International Society of Urological Pathology
	LFE	Local Frequency Estimation
	LNE	Laboratoire National de métrologie et d'Essais
	MD	Medical Device
	MEG	Motion Encoding Gradient
	mpMRI	multiparametric MRI
	MR	Magnetic Resonance
	MRE	Magnetic Resonance Elastography
	MRI	Magnetic Resonance Imaging
	PIRADS	Prostate Imaging Reporting and Data System
	PML	Perfectly Match Layer
	PSA	Prostate Specific Antigen
	ROI	Region of Interest
	SE	Spin Echo
	SNR	Signal to Noise Ratio
	TE	Echo Time
	TR	Repetition Time
	US	Ultrasound
	ZC	Central Zone
	ZP	Peripheral Zone
	ZT	Transitional Zone

Table 5 -

 5 1: Bilayer simulation properties

		Medium1: Top layer	Medium2: Bottom layer
		1000	1000
		4	1 6
		0	0
	c	1500	1500
	cS	2	4

Table 5 -

 5 2: Measured velocities (m/s) in the CIRS phantom

	Values are reported as mean ± standard deviation
	INCLUSION MANUFACTURER		FULL FIELD	CURL FIELD	SHEAR FIELD
	Near Field	2.45	2.44±0.32	2.34±0.22 2.58±0.26	2.15±0.22 2.07±0.22	2.32±0.17 2.25±0.17
	Type I	1.41	1.54±0.66	1.64±0.18	1.49±0.13	1.67±0.14
	Type II	1.87	1.87±0.37	1.94±0.09	1.77±0.07	1.90±0.07
	Type III	3.56	3.08±0.58	2.64±0.10	2.41±0.18	2.46±0.07
	Type IV	4.89	3.71±1	2.96±0.21	2.59±0.29	2.67±0.15

Table 6 -

 6

	1: Poroelastic parameters
	Symbol Property
	k 0 Permeability
	K b Matrix Bulk modulus
	K f Fluid Bulk modulus
	K s Solid Bulk modulus
	N Shear modulus
	Tortuosity
	Dynamic fluid viscosity
	f Fluid density
	s Solid density
	Porosity
	Experimental observations of this wave are scarce. Most of the observations were
	performed in rock-like media

Table 6 -

 6

			Hz	55 Hz	68 Hz	84 Hz	100 Hz	116 Hz	Hz
	cS	MRE 5.2±0.7 7.2±1.2 8.4±1.4 US 5.6 6.5 7.3	9.2±1.3 8.2	9.2±1.7 9.1	8.6±2.4 7.1±2.8 10.1 11
									9
		US	6.5	7.9	9.0	10.4	11.8	13.2	14.6

2: Measured velocity values in the foam phantom [m/s] Shear (cS) and compression (cP) velocities measured in this MRE experiment, and by Aichele [141] using ultrafast ultrasound (US). Values are reported as mean ± standard deviation, in m/s. c MRE 7.7±1.3 9.4±1.7 11.2±1.9 13.5±2.8 14.1±3.2 14.2±4.5 8.3±4.

  8 time points were recorded over a cycle. The whole acquisition was repeated at different frequencies: 40, 72, 105, and 138 Hz. Echo Time (TE) and Repetition Time (TR) were 75 and 1641 ms at 40 Hz, 81 and 1521 ms at 72 Hz, 64 and 1579 ms at 105 Hz, and 84 and 1518 ms at 138 Hz. The Field Of View was 36 x 36 x 12cm

Remerciements

Tristan DERUELLE

Computation of the potentials

Integral solutions

From there, and can be computed by their 3D integral forms [START_REF] Gui | A RIGOROUS AND COMPLETED STATEMENT ON HELMHOLTZ THEOREM[END_REF]:

(5-4)

( ) = 1 4

(5-5)

Where is the nabla operator with respect to r' the integration variable, is the outward surface normal, V is the bounded domain, and S the surface that encloses V.

For each space point, the whole integrals must be re-computed making this method quite slow and computer-demanding.

Alternatively, the equations can be driven further.

Vector potential

Applying the curl operator to equation (5-1) yields:

= ( ) +

(5-6)

= +

(5-7)

The curl of the dilatational field is zero. Therefore:

= =

(5-8)

Using the vector identity for the curl of curl, and using the fact that is also divergencefree [START_REF] Gui | A RIGOROUS AND COMPLETED STATEMENT ON HELMHOLTZ THEOREM[END_REF]:

This yields to the well-known vector Poisson's equation, solving for and knowing :

=

(5-10)

Equation (5-10) can be solved using any known technique for solving the discrete Poisson equation. Several methods are described in a book by Iserles [START_REF] Iserles | A First Course in the Numerical Analysis of Differential Equations[END_REF].

In the present method, we discretize the Laplace operator (using centered finite differences) and solve the system for each component:

Solving for :

(5-20)

Solving for :

And finally:

Finally, is obtained by an inverse Fourier Transform.

Equation results in a very simple and fast solution.

Here, the derivation was done with the Fourier Transform, but it could be performed with others (Discrete Cosine Transform (DCT), Discrete Sine Transform (DST)…).

Tristan DERUELLE 73

Compression displacement field

In the spectral domain, can be computed as follow.

Solving for :

(5-27)

(5-28)

Solving for :

And Finally:

This equation can also be solved using other transforms.

We note that summing equations (5-25) and (5-32) returns indeed = + .

APPLICATION ON SIMULATED DATA

Choice of the simulation tool

The algorithm was first tested on simulated data with known and controlled velocities.

Effort was put into choosing the right simulation tool. Both shear and compression should be handled and coupled. In real tissues, the speed of sound is greater than the shear wave speed by three orders of magnitude. The software should be able to approach such extreme velocity difference while keeping a reasonable computing time. Finally, for future work, it should be possible to add viscosity.

The ideal solution would be to use analytic solutions from Green's function [START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF]- [START_REF] Kausel | Fundamental Solutions in Elastodynamics[END_REF], but it is limited to homogeneous media. compressional velocity (C, 5.8 and 6.6 m/s) are completely erroneous compared to the true value of 1500 m/s. This is because the LFE algorithm is based on wavelength estimation, and it is not possible to recover the extremely large compressional wavelength (15m for a FOV of 10cm). The field separation algorithm works well. It is efficient at retrieving the shear field and compressional field. Close to the source, the shear velocity can be more accurately recovered but the compressional velocity is not expected to be recovered in common soft tissues. Note that similar results are obtained with a homogeneous medium (see Appendix D).

PERFORMANCES

The performances in computing the shear field are compared to the full field and to the curl field (standard to remove the compression contribution). The previous part showed the precision and accuracy of the new technique. This part will focus on resolution and noise immunity.

Resolution

The resolution between the different methods is assessed. This is done using the bilayer simulation. A sigmoid function is fitted to the median vertical profile:

Where A is the amplitude of the step, C is a constant offset, k controls the slope of the step, and x0 the position of the step.

Conclusions

Tristan DERUELLE 115 tissues evolution such as cancer development, liver fibrosis, Alzheimer's disease progression…

APPENDIX

A Extension of multidirectional filtering in 3D

A Reconstruction of tissues properties was shown to be improved if data were first filtered along multiple directions, processed, and then combined to make the final map. Manduca et al [START_REF] Manduca | Spatio-temporal directional filtering for improved inversion of MR elastography images[END_REF] used orthogonal filters with 2D LFE, and Tzschätzsch et al. [START_REF] Tzschätzsch | Tomoelastography by multifrequency wave number recovery from time-harmonic propagating shear waves[END_REF] used 12 directions with 2D kMDEV. However, MRE data are usually acquired in 3D, and it is preferable to perform the reconstruction with 3D algorithm to account for the actual shape of the wave field. Orthogonal filters are easily transposed to 3D. Using more than the 6 orthogonal directions is more challenging as it requires to evenly divide the 3D space into N directions. This is the problem that is addressed in this appendix.

A.1 ORTHOGONAL DIRECTIONS

Making orthogonal directional filters is relatively straightforward:

Where is the unit position vector and is the direction vector (e.g. along ).

For example, for the filter in the +x (1;0;0) direction, the resulting expression is:

and + .

It can be demonstrated by common trigonometry that the sum of all directional filters has no angular dependence (and in this case, equals 1). This ensures that no specific direction is preferred. This can be easily extended to 3D by introducing the cartesian expression of the spherical coordinates:

= cos( ) cos( ) cos( )sin( ) sin( )

Where is the azimuth in the x-y plane, and is the elevation away from the x-y plane: The dot product from equation (A-1) now becomes:

= cos( ) cos( ) ( ) + sin( ) sin( )

It can also be demonstrated that the sum of the directional filters does not depend on nor .

A.2 INCREASING THE NUMBER OF DIRECTIONS A.2.1 The 2D case

In 2D, instead of the four previous orthogonal directions, one may want to use more be easily achieved, and the i th direction is given by:

Tristan DERUELLE Tzschätzsch et al. [START_REF] Tzschätzsch | Tomoelastography by multifrequency wave number recovery from time-harmonic propagating shear waves[END_REF] used a radial gaussian type distribution centered on the i th direction.

= ( )

(A-6) Because the variance equals the angle step, the sum of all gaussians will be flat, without any preferred direction.

A.2.2 The 3D case

In 3D, it is more complex to define the different directions. In 2D, it is rather easy to evenly distribute N points on a circle. In 3D, distributing N points over a unit sphere has been the subject of numerous studies in various domains. It joins the Thomson and Tammes problems. Except for a few particular cases, there is so far no exact solution, and depends on the assessment of uniformity.

The Thomson problem is to evenly distributes N electrons on a sphere by minimizing the electrostatic potential energy [START_REF] Thomson | XXIV. On the structure of the atom: an investigation of the stability and periods of oscillation of a number of corpuscles arranged at equal intervals around the circumference of a circle; with application of the results to the theory of atomic structure[END_REF]. Common algorithms to estimate the solutions are based on local optimization. The Tammes problem is a particular case of the Thomson problem by using spherical arrangement (packing circles) [START_REF] Tammes | On the origin of number and arrangement of the places of exit on the surface of pollen-grains[END_REF]. It has also been studied while analyzing the structure of the sunflower [START_REF] Vogel | A better way to construct the sunflower head[END_REF]. The latter showed that the Fibonacci sequence can be used to distribute points over a disc. Martin Roberts [START_REF] Roberts | How to evenly distribute points on a sphere more effectively than the canonical Fibonacci Lattice[END_REF] explains how A Extension of multi-directional filtering in 3D

Tristan DERUELLE 131 the Fibonacci lattice can be used to distribute points over a sphere by minimizing the minimum nearest neighbor distance or the average nearest-neighbor distance. Both techniques yield similar results. The particle repulsion algorithm can only provide results for N>13 and requires computer resources to iteratively optimize the solution. The Fibonacci algorithm accept any number of points, is very easy and fast to compute, but in some cases, singular points can appear at the poles of the sphere.

D

This simulation consisted in a 3D homogeneous medium, with a shear velocity of cS = 2 m/s and a compression velocity cP = 1500 m/s. The source was modelled as the same circular membrane as the previously described, with a diameter of 8 cm and vibrating at 100 Hz.

The field of view was 10 cm x 10 cm x 10 cm with 86 pixels in all directions. Perfectly Match Layer (PML) was added on all sides with a thickness of 34 pixels. Computation time was in the order of 2 weeks.

The resulting displacement fields are shown in Figure D-1. The patterns are similar to the bilayer simulation (Figure 56). The full field and shear field are relatively similar, but present the main difference close to the source (top of the field of view -FOV). In that region, the wavelength from the full displacement field appears longer than that of the shear. This effect is due to the superposition of the compression field. When looking at the compressional displacement field, we can see the large and rapidly decaying wavelength arising from the membrane.

The constitutive equations are also satisfied. However, similar to the bilayer simulation (Figure 5-6), small wavelengths are visible in the compression field. They are of the same size as in the curl field, suggesting numerical errors. This shows the limit of the simulation. The same behavior is visible on the velocity maps (Figure 5678). The velocity from the full field (A) is slightly overestimated (2.2 m/s) and presents some artefacts. The reconstruction from the shear field (B) is homogeneous and centered around the expected 2.0 m/s. The presence of the compression biases the estimation of the shear velocity. The values of the compression velocity (C, 5.8) are completely erroneous compared to the true value of 1500 m/s. This is because the LFE algorithm is based on wavelength estimation, and it is not possible to recover the extremely large compressional wavelength (15m for a FOV of 10cm).