Keywords: Small animal SPECT, Dynamic images, Respiratory motion, Respiratory gating, 99m Tc-pertechnetate biodistribution, Compartmental analysis . . . . . . . . imaging system . . . . . . .

The aim of this thesis is to investigate temporally evolving phenomena with the use of single photon emission computed tomography (SPECT).

The rst part of this thesis addresses the problem of respiratory motion in SPECT images of mice. The presented method permits us to detect the respiratory motion in 4D SPECT images, to extract an intrinsic respiratory signal and to determine the widest possible phase of the respiratory cycle without movement. The data recorded during these motionless phases are then used to reconstruct a single 3D SPECT image without motion artefacts per acquisition. Reconstructed motionless SPECT images present a good compromise in terms of statistics and accuracy of the measurements with respect to basic 3D SPECT and 4D SPECT images.

In the second part, we study the iodide uptake kinetics in the stomach of mice through the study of 99m Tc-pertechnetate biodistribution with the use of 4D SPECT images. To understand the biological role of the iodide accumulation in the stomach, we modelled the phenomenon with a compartmental analysis approach using a simplied two-compartment (stomach wall and cavity) model with one input (blood). Time activity curves (TAC) of each compartment are deduced from observations and a rst estimation of the parameters was obtained.

The context of this thesis is the study of longitudinal phenomena. Our aim is to use imaging in order to understand the evolution of biological processes.

The objective is to provide tools to perform longitudinal studies of these processes.

In order to follow a phenomenon over time, we need to gather enough information, but also we need this information to be as accurate as possible. That is why, our objective in this thesis is to provide to biologists tools to extract as many information as possible from the imaging data at their disposal. Image processing techniques already provide solutions to dicult problems. Our approach consists here in developing simple methods based on some of the image processing techniques. These methods are dedicated to the resolution of precise biological problems.

The context of this thesis is also a collaborative work between two research teams: the ASCLEPIOS team at INRIA that focuses on biomedical image analysis and physiological systems modelling among others, and the TIRO team at Commissariat à l'Energie Atomique et aux énergies alternatives (CEA) / Université de Nice Sophia Antipolis (UNS) / Centre de lutte contre le cancer Antoine Lacassagne (CAL) that focuses on iodide metabolism: thyroid function, radiotoxicology, reporter gene imaging and treatments in oncology. This collaboration was materialised by the purchase of an innovative small animal camera as part of the contrat de projet État-Région (CPER), application le Telius proposed by INRIA and its collaborators.

Consequently, the biological problems we address in this thesis are related to the in vivo tracking of iodide uptake in cells with small animal single photon emission computed tomography (SPECT) imaging, and in particular mice imaging. It is of particular interest to remind that longitudinal imaging is really important for longitudinal studies comparing to the animal sacrices that would be required otherwise.

To summarise, the direct impact of this thesis is to provide tools to precisely characterise the studied longitudinal phenomena, observed on mice, through imaging. Most problems we encountered are related to the movement of the anaesthetized small animal. Getting accurate measures from images in the presence of such movement is among the main challenges addressed in this thesis.

In this thesis we illustrate how to deal with evolution problem with two dierent contributions:

A method that address respiratory motion for imaging of tumours. We propose a procedure to reconstruct SPECT motionless images with a good signal-to-noise ratio.

A representation of the stomach tolerant to its deformations and to animal motion in order to understand its iodide uptake. We propose a way to measure this uptake in the stomach with a layer-based approach that is much more tolerant to stomach contraction and relaxation than a classical voxel-based approach.

Material and experiments

The pre-clinical camera used in this thesis allows for 4D SPECT images acquisition. The four dimensions refer to the three spatial dimensions plus a fourth temporal dimension. This temporal dimension is characterized by the fact that several acquisitions are processed in a row, at regular temporal intervals (minutes, days or weeks). SPECT images allow the observation of the evolution of biological phenomena with respect to the biodistribution of radiopharmaceutical injection.

In the case of temporal interval on a daily or weekly scale, SPECT is dedicated to the study of processes such as tumour growth. For example, the tumour uptake sites can be localized and quantied, and by extension, their growth can be estimated.

In the case of short temporal interval on a minute or hour scale, SPECT allows the observation of radiopharmaceutical biodistribution variations after injection.

Respiratory motion

The rst challenge addressed in this thesis is the reconstruction of high contrasted and quantitatively reliable SPECT images of peritoneal tumours of anaesthetized mice. Among the dierent physical artefacts, respiratory motion is the most important one [START_REF] Hwang | Assessment of the sources of error aecting the quantitative accuracy of SPECT imaging in small animals[END_REF].

Respiratory motion may aect the tomographic reconstruction of positron emission tomography (PET) or SPECT images, which subsequently impair quantitative measurements, e.g. in the upper abdomen area. Respiratory signal phase-based gated reconstruction addresses this problem, but deteriorates the signal-to-noise ratio (SNR) ratio and other intensity-based quality measures. In this thesis, we propose a method for the reconstruction of 3D images that are not impaired by the respiratory motion and keep a high SNR.

The expected impact of this work is to provide a method that guarantees motionless image reconstruction and that allows accurate measurement and understanding tumour growth. The biological objective here is to be able to track, to model and ideally to predict the evolution of tumours in order to treat them more eciently.

Dynamic distribution of iodide in stomach

The second study is incorporated in the framework of the study of iodide metabolism in the murine stomach. The idea is to understand the gastric cycle of iodide. The rst objective is to automatically segment an organ in layers in a set of successive three (spatial) dimensions (3D) SPECT images, where the organ is animated by a random movement. The second objective is to model the dynamic distribution of iodide with the use of compartmental analysis, where each compartment is a functional compartment at the cellular scale.

In order to understand the physiological role of iodide, we rst study the process of iodide uptake. Then we elaborate a model of the stomach using compartmental analysis. Our contribution concerning this phenomenon is more a study of feasibility and a methodology than the exact characterisation of the process. Our results can be summarised as follows.

We realise a regional decomposition of the (3D) stomach in layers in order to be tolerant to its contraction and relaxation, and to the animal movement.

We elaborate a compartmental model for the iodide transport in the stomach.

We dene mixing coecients that allow us to compute the activity in each compartment from the measured activity in each layer.

We design a combination of analytical and numerical resolution steps allowing the full characterisation of the iodide uptake.

Finally, we compute the transfer coecients for our compartmental model; those coecients are the parameters of the compartmental model. The results are reliable enough for a rst study.

Manuscript organisation

This manuscript revolves around two parts demonstrating the two main contributions:

1. The rst part deals with the respiratory motion and motionless reconstruction of small animal SPECT images.

2. The second part describes the modelling of the iodide biodistribution in the murine stomach.

Chapter 2 gives a background on the biological and technical context of this thesis. It also describes the inherent material and animal model constraints.

Part I focuses on respiratory motion and consists of two chapters. Chapter 3 presents a state of the art about respiratory motion correction methods.

Chapter 4 presents in its entirety the journal paper [BMG + 13], we where we propose a dedicated method for optimal motionless reconstruction of 3D SPECT images.

Part II focuses on the iodide biodistribution in the murine stomach and consists of two chapters. In Chapter 5, we describe qualitatively the biodistribution of the iodide analog 99m Tc-pertechnetate in the stomach based on the observations of dynamic SPECT images. In Chapter 6, we characterise the iodide uptake phenomenon with a simplied compartmental model.

Lastly, Chapter 7 concludes this thesis with the list of contributions and directs us towards the feasible perspectives to this work. A list of publications is provided page 133.

The most useful appendices of this thesis are the following. Appendix A gives additional results for the motionless reconstruction method. Appendix B describes the principles of compartmental analysis.

Introduction

During this thesis, a great deal of eort has been done in order to work on common subjects with a team of biologists. From the denition of the subject to the interpretation of the results, a permanent dialog has been maintained between the researchers of the two domains. On the one side, biologists shared their experience, knowledge, and working constraints related to the study of the small animal. The advantage of this dialog and experience sharing was to dene and to set up a custom-made experimental protocol that would allow us to answer the questions that are of real interest for the biologists. On the other side, we also investigated the material at disposal -a SPECT/CT camera -in order to determine the material constraints.

In this chapter, the idea is to describe at best the tools, material, and techniques that were used in the presented works, and to dene their limits.

First, the small animal imaging camera will be described in order to determine the material constraints that have to be considered in the protocol elaboration Chapter 2. Technical context for the experimental protocol elaboration (Section 2.2). These constraints are to inject the smallest dose (radioactivity or X-rays) to the animal, to use the shortest acquisition sequence, and to acquire image without artefacts that could lead to wrong analysis and interpretation. Second, the creation of an appropriate animal model is also a complex task for biologists (Section 2.3). The constraints are numerous and mainly concern how to successfully elaborate an animal model with tumour growth and how to maintain this animal alive over several weeks with repetitive acquisitions so that it is possible to image the tumours along the time.

The knowledge of both material and animal constraints were taken into account for the denition of the nal experimental protocol.

2.2 Material constraints: Pre-clinical imaging system 1

The imaging system used for this work is a coupled SPECT/CT camera dedicated to small animal (eXplore speCZT-CT 120, GE Healthcare/Gamma Medica (London, ON, Canada)). This camera consists of three distinctive parts: a small animal table, the SPECT module, and the CT scanner (see Figure 2.1(a)). The rst module, the animal table, is used for the animal preparation. The table is the place where the anaesthetised animals are prepared and installed for acquisitions. The animal preparation and monitoring steps are described in Section 2.3 below.

The two imaging modalities are presented separately in this section in order to understand the functioning of each modality, and to identify the constraints of their use in an imaging protocol. A SPECT/CT camera permits to image two dierent types of information from the same subject. The idea for using such a multimodal camera is to combine both modalities in order to have more detailed information and images of the studied phenomenon.

However, in this camera, the SPECT module is located in front of the other module CT scanner, and therefore they are not localised on the same gantry.

The acquisition of each modality can not be processed simultaneously and images of each modality can only be acquired one after another. Nevertheless, SPECT and CT images corresponding to the same acquisition sequence can be mapped together afterwards, after registration using the global rigid transform between the two gantries. This rigid transform is established from preliminary test acquisitions. 1 A big thank you to Philippe Franken for is help at all time for camera understanding and handling, and support during acquisitions.

Single photon emission computed tomography (SPECT)

Camera description The SPECT module allows for imaging the biodistribution of radiopharmaceutical. The resulting SPECT images exhibit physiological and functional information about biological processes. This information is related to radiopharmaceutical distribution and accumulation through the subject body.

In order to acquire biodistribution images, the camera is equipped with a rotating seven-pinhole collimator (Figure 2.1(b)) and a stationary gantry of gamma photon detectors. The ring of detectors (Figure 2.1(d)), made of ten panels of four semi-conductor CZT modules, will record all events with a one-ms sampling.

The CZT modules serve for the detection of the gamma photons, also called scintillation. Usually, the scintillation phenomenon describes the gamma photon absorption by the NaI(T1) crystals in the camera [START_REF] Buvat | Single Photon Emission Computed Tomography[END_REF]. The crystal scintillates in response to incident gamma radiation. When a gamma photon is emitted by the animal and hits the crystal, it induces an electron loose from an iodine atom in the crystal, and a ash of light is produced when the electron again nds a minimal energy state. Then, the ash of light (called event) is detected by photomultiplier tubes (PMTs) behind the crystal and recorded. In our case, the scintillation process with the CZT modules is different from the one with scintillating classical NaI(T1) crystals. The energy of the detected photons is directly converted into electrical energy. The detection information is recorded with a one-ms sampling under the form of a list, called listmode, that contains the characteristics (time, detection location, energy, gating trigger) of each detected photons.

Then, in order to reconstruct an image from the listmode, we need to know the direction of each detected photons and the location of their emission in the animal. The direction can be deduced with the use of a collimator. The collimator is a lead tube with slits or pinholes, inside which one the animal is positioned. In our case, the collimator is a seven-pinhole collimator. The principle of collimation with a pinhole (or cone-beam) collimator is simple.

Gamma photons will be emitted in each directions from the animal but only the ones that will go through the collimators pinholes will hit the detectors.

The direction of the detected photon is thus determined from the position of the detection and the pinhole. The collimator selects the gamma photons according to their direction. The other photons will be stopped by the lead collimator in order to specically detect photons emitted according to a conebeam geometry. Since the collimator is going to select photons according to a certain direction with respect to each pinhole, the collimator is rotating in order to acquire data from each angle around the animal projections. In the case of the seven-pinhole collimator, the system acquires 51 projections for a single eld of view (FOV). One FOV measures 32 mm with mice collimator.

In addition to the photon direction selection, a selection of the detected photons according to their energy (spectrometry) is processed in order to discard the photons that were Compton scattered. This selection is done by choosing an energy window around the energy of the emitted photons.

For this work, the radiopharmaceutical that has been chosen is the 99m Tcpertechnetate, and the gamma photons are emitted with a 140 keV energy.

The width of the energy window of 99m Tc ranges from 126 to 154 keV with this camera.

Consequently, SPECT acquisition with the eXplore speCZT camera consists in recording information about gamma photons that are emitted during an user-dened time over the 360°around the animal and detected on the CZT detectors. The SPECT images are reconstructed with a block-iterative process, the ordered subset expectation maximization algorithm [START_REF] Malcolm | Ordered Subsets of Projection Data[END_REF], that has been proved to be fast. From seven to nine subsets and from six to nine iterations were used for the algorithm, the number of EM equivalent updates ranged from 54 to 56.

However, a certain number of physical phenomena aect the SPECT imaging [START_REF] Buvat | Single Photon Emission Computed Tomography[END_REF], and thus degrade the quality of the obtained image, such as attenuation, scatter, partial volume eect, as well as by motion artefacts.

Attenuation is characterised by a diminution of the energy of the photons.

Scatter is characterised by the deviation of the emitted photons. The two steps of cone-beam collimation and spectrometry take care of the attenuation and scatter eects. Therefore, respiratory motion compensation is one of the main obstacles of accurate quantication in micro-SPECT (µ-SPECT).

Radiopharmaceuticals Dierent factors have to be taken into account in order to choose the radiopharmaceutical:

the cells or molecules that are targeted; the half life of the radioisotope;

the ease of supply.

For this work, the radiopharmaceutical that has been chosen is the 99m Tcpertechnetate. Considering its chemical properties, 99m Tc-pertechnetate is considered as an iodine analog. That makes of it the perfect radiopharmaceutical for our purpose, which is tracking tumourous cells that express sodium/iodide symporter (NIS) protein (see Section 2.3). Indeed, as an iodine analog, 99m Tc-pertechnetate is transported into the body through the same paths, it is mediated by NIS protein, and it accumulates in these tissue.

99m Tc has a short half life (around 6 hours), which is perfect for frequent animal imaging and radioactive waste management. As the radioactivity will decay quickly, it can be considered as null after only a few days. Injected animals can then be safely considered not having radiopharmaceutical anymore within a week. This short half life make the storage of radioactive contents easier: the radioactive waste (syringe, paper towel or even dead animal) have to be put in quarantine for only a short period of time (typically a week).

Gating or not? During SPECT acquisition, data are recorded continuously.

If the region imaged is animated by a motion, one can desire to reconstruct images without motion artefact in order to to study the biodistribution of the radiopharmaceutical. When the motion is cyclic, motionless images can be reconstructed by using only a fraction of data at a specic motionless period of time of gated acquisition.

Retrospective gated reconstruction is possible when the acquisition has been triggered by an external monitoring system. External monitoring systems which are used to record a signal can be associated with a physiological motion (see Section 2.3.2 for further details on physiological signal extraction).

Gated SPECT images are then reconstructed within a specic time window of the periodic motion. The shorter the window of reconstruction, the less the motion artefacts in the image, and the smaller the signal to noise ratio. In order to prevent from a too high noise level, a compromise has to be done between the window of reconstruction and the motion artefacts in the image.

For further details on motions that can aect SPECT imaging and the method for reconstructing free of motion image, see Chapter 3.

Computed Tomography (CT)

Camera description The CT scanner allows for anatomical imaging of small animals. The CT imaging part of this scanner consists of a rotating system made of one X-ray tube and one at panel of detectors as presented in Figure 2.1(c). 3D images are reconstructed from a set of X-ray projections that are acquired at dierent angles around the subject. Image units are expressed in Hounseld units (HU) as for clinical CT scanners. This means that each kind of tissue has a dierent radiodensity with respect to water according to the expression: HU = 1000 * µ X -µ water µ water where µ X is the average linear attenuation coecient in the image voxel X.

Thus water has a zero HU. Air, that is not attenuating the X-rays has a radiodensity of -1000 HU. While bones which is attenuating X-rays, has a high radiodensity from 700 HU to 3000 HU. The distinction between soft tissue and organs is harder since they all contains water in similar proportion.

They have radiodensity ranged from 15 to 300 HU. In our case, we are not interested in bones, so we would like to favor soft tissue contrast with: fast acquisition; low dose.

Acquisition settings There are several predened acquisition settings for the CT scanner. Each setting allows for either faster acquisition, better image resolution, lower noise (HU), or lower dose (mGy). Thus a compromise between information, image resolution, acquisition time, and dose (X-rays radiation) has to be done (see the Appendix D for more details on each available setting).

The fastest sequences are the ones done only over 180°, with 220 projections, while other sequences acquire between 360 and 1200 projections over the 360°. These are numbers of projections required for a single FOV. If the region of interest is larger than one FOV, the number of projections required for the reconstruction is increased. In our case, the CT FOV is around 80 mm-length, it is large enough to image both lungs and abdomen.

Besides the number of projections required for CT image reconstruction, the properties of the X-rays tube (voltage and current) can also be adjusted.

The tube current controls the dose that is received by the subject, while the voltage only controls the energy of the emitted X-rays and does not aect the dose. A lower voltage means less penetrating photons and a better contrast in soft tissues. A higher voltage means more penetrating photons and the resulting image with a better contrast in bone. The dose is proportionate to the tube current by the acquisition time. The acquisition time is a multiple of the number of X-rays shoot by the duration of one shoot. Thus, for lower dose and better contrast in soft tissue, the setting with the lowest tube current and tube voltage would be the best. For the three constraints mentioned previously, the best setting is the Fast Scan setting (less than 2 min, tube current: 32 mA, tube voltage: 70 kV, and low dose: 36 mGy).

Gating or not? However, Fast Scan setting does not take into account motion information during the acquisition and reconstructed lungs images will be impaired with a blurred peri-diaphragmatic area. It is possible to obtain gated CT images with another CT acquisition setting, called Cardiac setting. This setting presents many advantages but in the scope of our study, the duration of a 4D is too long and increases incidence of animal mortality.

This setting is not further detailed here since we did not investigate motion correction method based on dynamic CT.

For the protocol, we had the choice of acquiring either a high resolution CT image in order to extract some information about the lesion growth, or a lower resolution and faster image that will be used only for rough anatomical information. In our context, the benet from CT image is low since the lesions that are studied are located in the abdomen, among soft tissues, and are not visible with a basic CT image. As a consequence, only a Fast Scan setting has been planned in the protocol for longitudinal studies.

Animal model constraints 2

For this interdisciplinary work, the presence of biologists is required and needed at all steps of the process, from the protocol elaboration to image acquisitions via the animal preparation. The protocol elaboration from the biologist side consists in dening everything that is related to the animal model preparation and the special care for it, before, during, and after imaging. This is also required for regulatory reasons.

Animal preparation

Animal model elaboration Elaborate the animal model consists in associating a specic mice lineage with a tumourous cell lineage that has been introduced in its body in order to generate life like tumours. The tumourous cell lineage is the result of in vitro development and studies beforehand. Cell lineage that are considered here have been modied in order to express a specic protein, the NIS. This protein is naturally expressed in thyroid cells and it is responsible for the uptake of iodine into these cells. The idea with this modication is to use radiopharmaceuticals based on radioisotopes (such as 99m Tc or 123 I) in order to localise these cells wherever it goes, xes and grows in the body. This is called reporter gene imaging.

Then, the lineage of the mice that will receive these cells has to be carefully chosen in order to allow the development of the tumour. Preliminary biological investigation on dierent mice lineage showed that syngenic tumourous cells expressing NIS protein caused a violent immune response at early stage of the development of the tumour. This reaction seems to be due to the Natural Killer (NK) lymphocytes that are able to distinguish altered cells from normal by recognising changes of a surface molecule. These NK lymphocytes are known to prevent proliferation of tumours and infected cells by destroying both of them. Thus, if we want to study the growth of tumours, we need to select a lineage with impaired immunodecient NK lymphocytes. Table 2.1 presents a list of mice lineage with its immune system characteristics. According to the constraints on immune characteristics, the most eligible mice lineage are Fox-chase SCID (®) Beige and NOD SCID. For the following, the mice lineage that has been chosen by biologist is NOD SCID mice.

Next, the tumourous developing site depends on the type of injection of the tumourous cells into the animals. Lesions will develop in the peritoneal cavity when cells are injected intraperitoneally. If injected intravenously, tumourous cells are distributed through all the body, and tumour might be developing in the lungs. When subcutaneous or orthotopic injection are done, tumourous 2 A big thank you to Julien Guglielmi and Audrey Lamit for their support and help at all time before and during acquisitions. Animal daily care Once the animal model has been set up, biologists have to take care of the animals in order to prevent from death and to ensure that the animals can be imaged according to the protocol until the end. Longitudinal study of animal requires the animal to be under daily surveillance. The daily surveillance can be completed by additional palpation in order to detect any disorder in the animal that could be linked with the tumourous cells injection and that could aect the future acquisitions. Acquisition time and frequency are determined according to the knowledge of biologists based on preliminary studies. Acquisition protocols are planned in order to minimise the side eects of an early imaging or frequent imaging. Biologist surveillance also aims to determine when the subject has to be sacriced in order to prevent the animal pain.

Animal preparation for imaging When an imaging session is planned, biologists are here to deal with the animal preparation before the acquisition.

The preparation is divided into three steps: radiopharmaceutical injection;

animal anaesthesia; animal positioning on the bed with motion and control sensors (anaesthesia, temperature).

Animal monitoring during acquisition

After the animal preparation section, this section is about the animal monitoring during acquisition. The more dicult part is to be able to control the anaesthesia and the breathing rhythm at all times of the acquisition, with the less side eects possible. The respiratory rhythm is particularly important since it indicates the state of the mice: deeply asleep, lightly asleep or almost awake. The respiratory motion is also important since the studied specimens point of interest are located in the area of the diaphragm (tumorous lesions in the lungs or abdomen). Therefore, this physiological periodic motion has to be studied and controlled, and eventually corrected. An heaviest anaesthesia can be used but it aects the animal ventilation.

Free breathing anaesthesia vs ventilation

These kind of anaesthesia consists of intra-peritoneal or intra-muscular injection of ketamine mixture. The mixture is injected only once before surgery or acquisition, sends animals to deep sleep, but the duration of sleep is not controlled and recovery is long. However, the animal can be placed in any position which allows better imaging or surgery without additional anaesthesic control [DCMW + 04]. In order to keep animal alive with a regular respira- tory rate, mechanical ventilation can also be required. Ventilation is mainly used for a fully controlled respiratory rhythm but it is heavy and invasive (deterioration of airways with tubes). This does not allow for repeated acquisitions. Thus, heavy anaesthesia (especially when coupled with ventilation) is not recommended on a regular basis as planned for longitudinal studies.

With supercial anaesthesia (gas inhalation such as isourane), the sleep is light and free breathing acquisition can be done with continuous anaesthetic gas supply. In that case, animal will recover fast, and repeated and longer acquisitions can be planned.

For the animal care, supercial anaesthesia is preferred.

Anaesthetic gas eect on breathing rhythm Although less invasive, to use anaesthesic gas also induce side eects on the respiratory behaviour. The respiratory-like signal can be deduced from various physiological information that can be observed with a range of tools: the deformation of the ribs and thorax with a video-based system, the alteration of the heart electrical activity with ECG, the variation of pressure applied to a pressure sensor attached to the chest or the abdomen, the variation of the air volume in the lungs measured with a spirometer, or the temperature variation of the air that is going through the lungs. Once extracted this respiratory signals are used for two things. First, it is used as an information about the breathing rhythm of the animal, according to which the anaesthetic gas ratio is adjusted in the air that the animal breathes. It is also used for triggering image acquisition.

In the coming paragraphs, dierent respiratory motion sensors that can be used to monitor the mice are described.

Opto-electronic or video-based system ing the movements of several markers positioned on the chest or on the abdomen of the subject with camera. In the case of a system called real-time position management (RPM), they use infrared light to track the infrared reective markers set up on the subject abdomen during an acquisition. In that case, the respiratory like signal is extracted from the spatial position variation of the markers.

ECG electrocardiogram (ECG) consists in recording and interpreting

the electrical activity of the heart over a specied period of time with the use of electrodes attached to the surface of the skin across the thorax and chest.

It has been shown that the electrical activity signal can also be interpreted in order to extract a respiratory like signal [START_REF] Moody | Derivation of respiratory signals from multi-lead ECGs[END_REF][START_REF] Travaglini | Respiratory signal derived from eight-lead ECG[END_REF]. Indeed, the heart is also submitted to respiratory motion and the heart electrical activity is impaired by respiratory disturbances. The cushion contains a sensor that detects and measures this pressure variation. A higher pressure applied to the cushion is associated to an increased volume in the lung and thus to the inspiration phase while a lower pressure is associated to the expiration phase.

Spirometer The spirometer [CJP + 04, GJQ + 11, JPK + 08] is another tool that records the volume of air exchanged by the lungs (inspiration and expiration) over a specied period. This information can be exploited in order to obtain a signal representative of the respiration motion. The signal is resulting from the measurement of the ow rates through a sensor in the tool that is xed to the head. Positive ow rates are then associated to expiration while negative ones are associated to inspiration.

Temperature This is another method that can be linked with the spirometer. Instead of measuring the ow rate, it consists in recording the temperature of the air incoming and outgoing from the nose during breathing [NE08, BPD + 10]. The temperature is compared to the air temperature in the room. A lower temperature is associated with the inspiration (air going from the room into the lungs) while higher temperature is associated with expiration. These method is used both for human and for animal imaging.

For all methods, the obtained respiratory-like signal can be used online (while recording) in order to control gated acquisition. Triggers can be generated at specic time of the respiratory cycle thanks to an external signal monitoring software. Those triggers are used for prospective gating acquisition of dynamic images in the case of CT. In the case of emission tomography (ET), the triggers are used for retrospective gated reconstruction.

Among all these systems, a few have been adapted for small animal, but all systems require specic tools. For the following, the respiration motion sensor that has been chosen is the pressure sensor, although this is probably not the most accurate. It has been chosen since it is a part of a mice monitoring set that is already on disposal, and that the biologist already know how to use for the control of the well-being of the animal during the acquisition (temperature, breathing rhythm). The signal that is obtained from this pressure sensor associated with the monitoring system can also be used for the detection of the respiratory phases. Indeed, the respiratory signal can be split into dierent temporal phases according to the amplitude variation associated with the motion phases. An ideal phase splitting would be that these phases were dened in such a way that the organ can be considered as motionless during each phase.

Additionally, we will see in Chapter 4 that the precision of the sensor is not necessarily as crucial here as in other works. Indeed, in our work, we only use the signal for SPECT imaging for a preliminary listmode gating step:

the delimitation of each consecutive respiratory cycle based on an amplitude motion thresholding. The suppression of the motion artefacts themselves will be done by a specic method that we will detail in Chapter 4.

Conclusion

One of the strong points of this work is its interdisciplinary nature, between image analysis and biology. The material presented in this chapter is the link between those two research domains. In this chapter, we focused on the material and biological aspects, and presented many related constraints. We summarise below the constraints, the decisions we took for this study and their consequences.

Animal model First, our objective is to study tumour growth, and in order to do this, we chose to observe cells expressing the protein NIS. For this, we have seen that the animal lineage had to be highly immuno-depressed (NK lymphocytes); we chose the lineage NOD-SCID.

We chose to inject the cells intra-peritoneally inducing tumour growth in the peritoneal area, subject to signicant movements.

Anaesthesia We chose isourane as anaesthetic gas which is the most exible anaesthesic technique since we want to maximise the duration of the experimental period and to optimise the quality of acquisitions during this period. The chosen anaesthesia is the less invasive and allows fast recovery for the health of the mouse.

Frequency of acquisition Here a compromise has to be made: frequent acquisitions gather more longitudinal information, but may interfere with the tumour growth. We chose to acquire images once a week.

Motion Due to the fact that imaged tumours are in the peritoneal area, the acquisition will be impaired by respiratory induced movements of the imaged area. Also the chosen anaesthesia, isourane, induces a change in the respiratory pattern. Consequently, the movement will have to be studied, monitored, and taken into account to improve the accuracy of our acquisitions.

Motion signal detection To monitor the animal motion, a pressure cushion has been chosen. Due to the particular respiratory motion pattern and the importance of the induced movement, we only use this sensor for detecting breathing cycles.

The following constraints depend on the imaging modality chosen.

For SPECT imaging, the problem here is to choose the adequate radiopharmaceutical product. It has to be adapted to the camera and to the cells we want to observe. Among the available products, we chose 99m Tcpertechnetate. Since 99m Tc half life is relatively short (6 hours), a compromise has to be done between the dose of injected radioactivity and the duration of the acquisition. If the dose is too low, the acquisition will have to be longer in order to obtain an image with a good signal-to-noise ratio. If the dose is too high, this could interfere with the tumour growth.

For CT imaging, the only constraints concern the camera settings that must be chosen carefully to ensure a good contrast and accuracy in the region of interest. The setting that has been chosen is Fast Scan, corresponding to a dose of 36 mGy for one FOV, a tube current of 32 mA, and a tube voltage of 70 kV. It has a good compromise between acquisition duration and irradiation of the animal.

Some works consider the joint use of CT and SPECT imaging, the main advantage is that CT image could be used for anatomical registration while SPECT provides biodistribution information [SDP + 03]. However, the anatomical information included in the CT is not sucient. If more detail on abdominal area is needed, the use of CT contrast agent could be considered, but this will also increase the diculties in preparing the animal (multiple injection with unknown side eects). Thus, we directed our investigations toward a method dedicated to SPECT, neglecting CT image information.

All the constraints and resulting choices presented in this chapter will have an impact on the methods that will be considered for the motion correction as we will see in the following chapter.

Introduction

To acquire static images from a moving target is a real challenge. It depends on the relative speed between the observed motion and the imaging modality.

The diculty lies in the determination of a type of acquisition that is adapted to both the motion of the organ of interest and the modality.

In clinical and pre-clinical imaging, dierent types of motion have been identied [START_REF] Menke | Compensation methods for head motion detected during PET imaging[END_REF]: involuntary subject motion with only one occurrence, physiological cyclic motion such as cardiac and respiratory motions, and physiological non cyclic motion induced by digestion for example. Involuntary subject motion cannot be predicted in advance but one can prevent the subject from moving by using anaesthesia and mold especially for infants, elderly, or animals. If the acquisition duration is longer than the motion duration, as it is the case in small animal single photon emission computed tomography (SPECT) imaging, specic methods have to be developed. Cardiac and respiratory motion have been frequently studied since there are cyclic and the duration of a cycle is shorter than some acquisition sequence.

Two types of acquisition approaches have been developed. The motion can be considered as soon as the image acquisition with gating approaches, or a posteriori during the image reconstruction. The choice is done according to the modality characteristics and the type of observed physical phenomenon:

X-rays for computed tomography (CT), gamma photons for SPECT, radiofrequency waves for magnetic resonance imaging (MRI).

Generally gated acquisitions are frequently used for CT since images are reconstructed from a set of X-rays projections. Each projection is obtained from a short time X-ray shoot at dierent angular views around the subject, and the X-ray shoots can be gated according to the signal (gated acquisition).

Concerning emission tomography (ET) images, the acquisition process is usually longer than the respiratory motion and motion is taking into account more frequently at the reconstruction step. In SPECT imaging, the observed phenomenon is one single gamma photon emitted by radioactive decay. In positron emission tomography (PET), the observed phenomenon are two 511 keV gamma photons being emitted at almost 180 degrees to each other by annihilation of one positron emitted by decay. The reconstructed tomography images will then give the biodistribution of the radioisotope, and as a consequence the distribution of the targeted cells. Radioactive decay is a continuous process and not controllable, so there are emission of gamma rays at any time of movement if one. Thus, the acquisition and reconstruction processes have to take into account the motion of the organ of interest at the time of acquisition.

As in this thesis we work with a pre-clinical SPECT/CT camera, we will focus this chapter on motion correction methods dedicated to these two modalities. First, we detail the motion eects in ET that justify the necessity to nd an appropriate motion correction method. Then, we go through dierent solutions for ET that have been developed and proposed in the literature. Finally, we discuss the advantages and drawbacks with respect to our subject in order to introduce the method that we present in Chapter 4 Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT.

Motion eects in emission tomography (ET) images

ET give images of radiopharmaceutical biodistribution in the subject. This modality permits to localise group of cells by their ability to accumulate radiopharmaceuticals.

As previously mentioned, the acquisition is done during a continuous period of time dependent on the activity of the injected radiopharmaceutical.

Thus, when there is a radiopharmaceutical uptake in a region submitted to motion, the source is moving and the emitted photons coming from the same source will hit the detector at dierent locations. The reconstruction will lead to a blurred image and quantication will underestimate the real target uptake.

In the context of our work, we wanted to be able to detect small intra peritoneal lesions as soon and as small as possible, and to follow quantitatively their growth along the time. Chapter 4 will describe our approach for motion detection in this context. If the lesion is submitted to motion, the emitted gamma photons will be spread around in a wider region of the image and the average activity in that region will then be reduced. As a consequence, the lesion can either be not detected or wrongly quantied. The rst case is frequently met when the average activity is too low and the hotspot cannot be distinguished from the noise.

When the average activity is still above the noise (second case), it is still possible to detect the lesion but the location is less accurate, the denition of the size will lead to a larger lesion volume and quantication will be biased.

Since a longitudinal study based on ET will consist in comparing lesion size, the quantication measure have to be reproducible especially for a group study. If the motion was identical for all subjects and at all times, the impairment will be the same and the measures would be comparable. However, respiratory motion is cyclic but not very reproducible, motion correction methods have to be developed in order to avoid the reduce the bias on the size measurement.

Methods for image acquisition and reconstruction of moving subjects

This section provides a more detailed state of the art to complement the Section 4.1.

How to get rid of motion eect in emission tomography imaging?

We identied two main techniques to prevent motion from aecting image acquisition and reconstruction: motion gating, which consists in acquiring or reconstructing images when the subject is not moving or when there is momentarily no motion tering each three (spatial) dimensions (3D) image of a sequence to the rst image [START_REF] Bai | Regularized B-spline deformable registration for respiratory motion correction in PET images[END_REF], and adding up all registered images to a single 3D image [START_REF] Bai | Motion correction and attenuation correction for respiratory gated PET images[END_REF].

In this section, we present motion gating approaches rst, then motion compensation methods, and nally one major variant of the motion correction methods.

For human, a trivial solution is to stop the considered motion at the time of acquisition, this is the purpose of breath-holding and this assumes that the duration of acquisition is short enough. The drawbacks of their approach is that it only works for short-time acquisition, basically CT, and cannot obviously be applied to animals.

For small animals, an other cheap method is to constrain the respiratory motion by banding the chest [DCMW + 04]. However chest banding changes the geometry of thorax, and aects the respiratory motion. So the motion might be reduced enough for the lesion detection from the operator point of view (no dorso-ventral or ribs enlargement). Indeed, the motion is just constrained and it might enlarge the inside motion in the craniocaudal axis.

That is why more sophisticated methods such as gating appear to be necessary for accurate quantitative study. In theory, any combination of gating methods is possible but according to the imaging modalities, some are preferred to others. is performed, assuming one-ms is a short time with respect to the respiratory cycle. Then post-reconstruction analysis on this dynamic PET permits the tracking the source motion. The position variation is assumed to be synchronous to respiratory motion. This data is then assimilated to respiratory signal. Gates or bins are then determined as previously described and data corresponding to specic position of this source during the respiratory cycle were used for retrospective reconstruction.

The advantage of this approach is indeed that it works for imaging without anatomical reference. It is also the best solution when one extracts the information from the imaging data, since it is harder to analyse in real time the imaging data in order to determine the gating. It is compatible with imaging without external monitoring system whether it be a choice or a constraint. The drawback for most of retrospective gating method is then to acquire enough data in order to make sure that the amount of data that will be kept for gated reconstruction will be sucient to not aect the quality of the resulting image. A third alternative motion correction approach consists of motion blur reduction. Xu et al. [START_REF] Xu | Respiratory motion blur identication and reduction in ungated thoracic PET imaging[END_REF] consider that the motion blur can be modelled as a point spread function (PSF). The 2D PSF is estimated directly from the reconstructed PET image. PSF parameters are estimated from image intensity derivation along dierent angular directions. Then the motion blur is corrected with a deconvolution algorithm using the resulting 2D PSF and wavelets.

Motion compensation

As a conclusion, motion correction approaches advantages are that they aim at using all data available for reconstruction, even though the motion itself is not well known, thanks to a wide range of registration methods (both for raw data or images of a sequence). However, for ET, these approaches often require additional correction steps for making sure that nal reconstructed image is correct in term of estimated activity in lesions.

Discussion

Remembering that our study is a preliminary work for lesion study, we want to nd the best approach that is also as simple and fast as possible with respect to the tools and constraints that we have. For ET, the ideal method would be a motion compensation method since this kind of approaches uses The complexity for elaborating the motion model encouraged us to look for a simpler method.

all
Motion correction methods also tend to use all the data. First registration step is usually done on anatomical modality images (e.g. CT): either on the data that will be used for the following, or on synchronised data with the other modality acquired simultaneously. This is then perfect when both coupled modality can be acquired simultaneously. However, in our case, the dynamic CT cannot be acquired at the same time as SPECT. Then, the registration and following sum up steps have to be processed directly on the dynamic SPECT, without the anatomical reference. All together this might degrade the main information that needs to be preserved: the activity. As a consequence, an additional interpolation step has to be considered and carefully chosen [HSO + 95]. The approach that is developed for gated images could be considered in our case at the condition that the gating information is reliable, and that all data used for the reconstruction of each phase correspond exactly to the same motion phase.

In our case the gating information is not as reliable as expected: a time shift between respiratory-like signal recorded thanks to the pressure sensor and the motion measured on data itself was observed. This suggests that pressure variation and the intrinsic motion information are not synchronous.

The question is then which motion information to believe. As the motion is impairing data, one decided to nd an approach that takes into account the intrinsic motion information.

This led us then to consider the motion gated techniques. Although the advantage is that gating approaches are rather simple, the diculty lies in the determination of the gates. If not well dened, the gated acquisition might not correspond exactly to the same motion phase for each cycle. Then, the reconstructed image will be impaired by motion. This problem can occur when the physiological signal is not perfectly cyclic (both with respect to amplitude or phase) and when acquisition are done in a prospective gating mode. This might be often the case in CT. In the case of retrospective gating reconstruction methods (for ET), only a small fraction of the available data is used for each gate. The reconstructed images will have a low SNR that is inappropriate for activity quantication. However regarding our technical constraints (see previous Chapter 2) and the nature of the animal respiratory signal (gasps), a method using gating approaches is the best choice. The main drawback is that gating in SPECT means using a fraction of available data.

What we decided is then to elaborate a gating method that will maximise the fraction of available data used with respect to the physiological signal and without a complex registration step.

In the next chapter, we will describe our method, based on well dened gating parameters. Among others, our method achieves in using around 65% of available data for the motionless reconstruction of small animal 3D SPECT image.

Chapter 4

Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT Respiratory motion can blur the tomographic reconstruction of PET or SPECT images, which subsequently impair quantitative measurements, e.g. in the upper abdomen area. Respiratory signal phase-based gated reconstruction addresses this problem, but deteriorates the signal-to-noise ratio and other intensity-based quality measures. This article proposes a 3D reconstruction method dedicated to micro-SPECT imaging of mice. From a 4D acquisition, the phase images exhibiting motion are identied and the associated list-mode data are discarded, which enables the reconstruction of a 3D image without respiratory artefacts. The proposed method allows a motion-free reconstruction exhibiting both satisfactory count statistics and accuracy of measures. With respect to standard 3D reconstruction (NG3D) without breathing motion correction, an increase of 14.6% of the mean SUV has been observed, while, with respect to a gated 4D reconstruction (G4D), up to 60% less noise and an increase of up to 124% of the SNR have been demonstrated.

Introduction

Molecular pre-clinical imaging is a major research tool which provides noninvasive in vivo information on cellular processes and allows longitudinal studies [PGW10, FAMH08, KC08]. In oncology, obtaining measurements of tumour characteristics is mandatory. Emission tomography (ET) quantication is hampered by physically induced biases such as attenuation, scatter and partial volume eect, as well as by motion artefacts. In clinical imaging, the physical biases can be corrected in Positron Emission Tomography (PET) as well as in Single-Photon Emission Computed Tomography (SPECT) [START_REF] Ritt | Absolute quantication in SPECT[END_REF]. However, there are few studies addressing specically the quantication issues in pre-clinical imaging. Concerning physical biases, Hwang et al. [START_REF] Hwang | Assessment of the sources of error aecting the quantitative accuracy of SPECT imaging in small animals[END_REF] showed that for technetium-99m ( 99m Tc), attenuation and scatter errors are reduced in small animal SPECT compared to clinical SPECT. Therefore, respiratory motion compensation is one of the main obstacles of accurate quantication in micro-SPECT (µ-SPECT).

Computed Tomography (CT), PET or SPECT images are produced by a tomographic reconstruction from the projections of the object of interest, with the implicit assumption that the imaged object remains still during the projection acquisition. Motivated by the respiratory and cardiac motion observed in clinical practice, a number of motion handling methods have been developed, with the goal of reconstructing images free of motion artefacts. A trivial approach consists in controlling the considered motion during the acquisition, by breath-holding, for example. This is only valid for human imaging and short-time acquisition techniques such as CT or Magnetic Resonance (MR).

Motion artefacts can also be handled by gating protocols. An additional signal, which is considered representative of the motion of interest, is recorded and synchronised with the data acquisition. When a motion model is available, some authors proposed to incorporate it into the tomographic reconstruction procedure in order to take into account all the acquired data to reconstruct a 3D still image. In human imaging, cardiac [START_REF] Blondel | 3D tomographic reconstruction of coronary arteries using a precomputed 4D motion eld[END_REF] and respiratory [LLC + 07, RMK + 07] motion models have been developed. In addition to the diculty in dening an accurate motion model, such methods also imply a huge computational cost.

In human imaging, motion compensation has been studied in CT or PET [START_REF] Sa Nehmeh | Respiratory motion in positron emission tomography/computed tomography: a review[END_REF][START_REF] Lucignani | Respiratory and cardiac motion correction with 4D PET imaging: shooting at moving targets[END_REF]. In pre-clinical imaging, the methods are similar. For In this article, we propose an amplitude-based gating reconstruction method in mice µ-SPECT. The signal obtained from an external pressure sensor is ltered and allows the reconstruction of an initial 4D image, which is then analysed to detect the respiratory induced motion. We take advantage reconstruction in micro-SPECT of the particular breathing pattern of anaesthetized mice to detect motioncorrupted phases, and subsequently use all raw data corresponding to the motionless phases to reconstruct a single motionless 3D image.

Materials and methods

Animal model

The dierent methods were tested on data obtained on 3 female mice of a mouse model. The mouse model consisted of NOD-SCID mice with intraperitoneally injection of 1.5 million cells of rat colonic adenocarcinoma expressing the Natrium Iodide Symporter (PROb-mNIS) in 50µL of Phosphate Buered Saline. SPECT/CT acquisitions were performed after 1, 2, and 3 weeks of growth, and PROb-mNIS peritoneal nodules of NIS-transfected cells were observed. However, we considered only 8 acquisitions out of the 9 for our study, since the Biovet pressure signal failed for one acquisition.

Animal housing and procedures were conducted according to the guidelines of the French Agriculture Ministry and were approved by the local ethics committee.

Imaging protocol

We acquired experimental animal images using a dedicated small animal SPECT/CT scanner (eXplore speCZT CT 120, GE Healthcare Bioscience, London ON, CA).

The SPECT imaging part of this camera consists of a xed full-ring of detectors coupled to a rotating 7-pinhole collimator for mice. The ring of detectors is composed of 10 panels of 4 CZT detectors. SPECT acquisitions lasted around 18 min. The output of the camera is a list-mode record of all CZT detected events : each event consists of the detection time (ms), the gating indices, and the properties of the detected photon. SPECT reconstructions used data from the list-mode in the energy window from 125 to 150 keV, according to the energy of 99mTc.

Mice were injected intra-peritoneally 180 MBq in 400 µL of 99m Tcpertechnetate. Five minutes later, the animals were anaesthetized with isourane (1.3% v/v) (Baxter, France). The animals were placed prone and freely breathing (without mechanical ventilation) on the bed and kept sedated during the whole imaging protocol, using inhaled anaesthetic.

A monitoring system (BioVet, m2m Imaging Corporation, Newark, USA)

was used with a pneumatic pressure sensor. The pressure sensor was placed under the animal abdomen. The pressure signal was recorded (1 ms-sampled rate) and used to monitor the animal anaesthesia.

Body temperature was maintained at 37

• Celsius and anaesthetized gas rate was regularly controlled to keep the respiratory rate between 60 and 100 breaths per minute. Once a stabilized respiratory rate was reached, the SPECT acquisition started.

At the beginning of each acquisition, an amplitude threshold was arbitrarily xed in order to detect the falling edge of the pressure signal variation.

Each time the signal passed this threshold, a trigger was recorded in the listmode SPECT acquisition data to be used for reconstruction. We assumed that this trigger corresponded to the same time in the cycle for each respiratory cycle. A respiratory cycle is then dened by the interval between two consecutive triggers.

Reconstruction algorithm

SPECT reconstruction was done using the Ordered Subset Expectation Maximization (OSEM) algorithm [START_REF] Malcolm | Ordered Subsets of Projection Data[END_REF]. We used from 7 to 9 subsets and from 6 to 9 iterations, the number of EM equivalent updates ranged from 54 to 56. For a given set of data, the number of subsets and iterations were kept constant for all the reconstruction methods.

"Standard" reconstruction methods

In this section, we describe the dierent kinds of SPECT reconstruction methods for the list-mode, resulting in either 3D or 4D images. We start with the reconstruction that are provided by the camera and then we introduce the improvements which lead to our nal reconstruction scheme: a breath-hold like 3D image.

Non-gated 3D reconstruction

Non-gated 3D reconstruction (NG3D) consists in reconstructing a single 3D image using all data recorded in the list-mode within the range of energy corresponding to the radiopharmaceutical used. The raw data (list-mode) used is the same as for NG3D reconstruction.

Proposed Method

The proposed method is a combination of three main steps that are described here: a cycle ltering which consists in cycle pre-selection according to duration, a cycle resampling, and an image-based motionless phase detection.

This method denoted by breath-hold like 3D reconstruction (BH3D) gives 3D

images without motion artefact.

Cycle selection

Cycle resampling

Obviously, even after the cycle selection step, the remaining respiratory cycles still exhibit some length variations.

About half of the cycles are longer than C s and their trigger time t j satises t j+1 -t j > C s . As a consequence, the counts after the last gating window of this cycle and before the next trigger signal (i.e. in the interval [t j + C s , t j+1 [) will not be used at all in the reconstruction process. This is illustrated by Conversely, the other half of the cycles are shorter than the average cycle duration C s and their trigger time t k satises t k+1 -t k < C s . As a consequence, the counts of the last gating windows dened for this cycle (in the interval [t k+1 , t k + C s [) also belong to the rst gating windows of the next cycle, meaning that those counts will be used twice in the reconstruction process for dierent phase images. This is illustrated by gure 4.2(a).

These observations motivate a temporal resampling of the counts so that the interval between two trigger signals is equal to the cycle average duration C s .This ensures that each event will be used once and only once in the reconstruction process. Therefore the detection time t of an event located between the trigger signals t i and t i+1 will be changed to This method is based on the observation that, around the borders of high The workow for detecting motionless phases can be described as followed:

t = (i -1)C s + (t -t i ) C s t i+1 -t i
1. Intensity normalisation of the 3D images I n , n ∈ [1, N ] resulting from the 4D reconstruction according to the proportion of available counts that has been used for reconstruction.

∀x : În (x) = I n (x)
total number of counts number of counts used to reconstruct I n 2. Spatial regularization with a 3D Gaussian lter (σ = 1.00) of the normalised images În yielding Ĩn .

3. Volume of interest denition. Since we are interested in abdominal lesions, we restrict the above analysis to a sub-volume containing only the abdominal area.

Generation of an amplitude image. A 3D image of intensity variation is

generated from the images Ĩn . For a given voxel x, the uptake amplitude reconstruction in micro-SPECT can be simply computed by subtracting the minimal uptake value from the maximal one, resulting in an amplitude image A.

∀x :

A(x) = max{ Ĩn (x)} n∈[1,N ] -min{ Ĩn (x)} n∈[1,N ]
5. Generation of a mask M of regions of potential motion. This mask aims at isolating voxels with large variations of time-course uptake values.

This mask is obtained by thresholding the amplitude image A at 15% of its maximum value, yielding segmented regions with potential spatial displacement.

6. Detection of the gasping breaths. We observe that around 2/3 of the respiratory cycle do not exhibit signicant motion. Since N phases are reconstructed, with N = 15, this correspond to around 10 phases without motion while the motion-induced large uptake variations are then mostly due to the other 5 phases. Therefore, for each voxel x, we compute the amplitude after discarding every possible sequence of 5 consecutive phases. We dene S m as a sequence of 5 consecutive phases:

S m = {k ∈ [1, N ]|∃j ∈ [0, 4] : k -m ≡ j[N ]}
The discarded sequence yielding the smallest amplitude is associated to a possible motion under the considered voxel. This sequence is dened by its rst phase m(x) computed by:

m(x) = arg min m∈[1,N ] max{ Ĩn (x)} n∈[1,N ]\Sm -min{ Ĩn (x)} n∈[1,N ]\Sm We cumulate the obtained sequences S m(x) into a histogram h. ∀n ∈ [1, N ], h(n) = card{x ∈ M/n ∈ S m(x) }
The resulting histogram is automatically thresholded using the Otsu method [START_REF] Otsu | A Threshold Selection Method from Gray-Level Histograms[END_REF], yielding two groups of phases. The smallest group (typ. SUV peak , SUV mean , and volume) were not directly comparable since they vary from lesion to lesion. Therefore they were compared using ratio paired t-tests that correspond to paired t-tests on their logarithms (base 10). illustrates both biological and imaging challenges as there are numerous small peritoneal carcinoma lesions, distributed throughout the abdomen.

Results

Lesions distribution

A total of 76 foci of 99m Tc-uptake were detected on the 8 acquisitions (3 to 14 lesions per acquisition). They were located in the abdomen, showing NISexpressing tumour tissue, and corresponded to peritoneal carcinoma nodules.

These nodules were small, their size (equivalent spherical diameter) ranged from 1.09 to 2.34 mm (average 1.55 ± 0.28 mm). Twenty-seven lesions were located in the upper abdominal area (above the stomach), 31 were located in the middle abdominal area (at stomach level), and 18 were located in the lower abdominal area (below the stomach). Changing the trigger parameters on the pressure signal (rising/falling edge, with/without delay) for the consecutive SPECT acquisitions of the same animal did not change the location of the peak (data not presented here). The peak of the motion detected in the 4D CT appeared at the same location as in SPECT, with a temporal shift with respect to the pressure signal (gure 4.5).

Respiratory signal analysis

Last, we computed the lung volume from the CT images along the phases with an ad-hoc and simple method (Gaussian smoothing, threshold between -650 HU and -250 HU, and morphological closing). Again, the peak of volume change appeared at the same location. There is a discontinuity in the lung volume variations (at phase #10): the corresponding reconstructed CT image appears to exhibit dierent characteristics which explained the lung segmentation dierences.

Sensitivity of the image-based motion detection with respect to noise

We investigated the ability of the proposed method to detect the motion in noisy images. To that end, we extracted from one 4D SPECT sequence a series of 3D images containing a moving lesion. This lesion has been thresholded, and dierent (white and gaussian) noise levels have been added to the binary images. It appears from our experiments that the image-based motion detection method is always successful in detecting the motion phases for SNRs above 2.5.

Image statistics

The number of counts used for each reconstruction method were compared to the total number of acquired counts. Table 4.2 presents the gures concerning the average percentages on the set of 8 acquisitions. By denition, the 3D reconstruction method NG3D used 100% of the counts.

15 3D phase images were reconstructed by the G4D method. The average percentage of counts used for the reconstruction of each 3D phase image is then around 1/15th of the counts (6.65%). However, because of the variation reconstruction in micro-SPECT of respiratory cycle duration, 3.64% of the counts are not taken in account while 3.51% of them are used twice (i.e. for the reconstruction of two dierent phases). If the total number of used counts is dened as the sum of the number of counts used for each phase reconstruction, it appears that this total number may exceed 100% if the number of counts used twice is larger than the number of omitted counts. As a consequence, the average total of counts used for the reconstruction of the 15 phases is 99.76% with a standard deviation of 2.88%.

The twice used counts were not evenly distributed among phases, and are more likely to occur in the rst phases: an average of 31.44% (resp. 9.67%) of counts used for phase 1 (resp. phase 2) were already used in a preceding phase. Twice used counts appeared up to the 9th phase.

The cycle selection step discarded an average of 3.25% of the counts. After the update of the average respiratory cycle duration, there were still 3.74% of omitted counts and 2.73% of twice used counts in a 4D reconstruction without cycle resampling. The trends in the repartition of the twice used counts among the phases were similar to those observed for G4D: e.g. an average of 26.86%

(resp. 7.61%) of counts used for phase 1 (resp. phase 2) were already used in a preceding phase.

After cycle resampling, all the retained counts are used once and only once in the 4D reconstruction G4DSR. Each reconstructed phase of G4DSR used around 6.45% of counts per reconstructed phase, which corresponded to a total of 96.75% counts over the 15 phases.

For the BH3D reconstruction, the counts from the equivalent of 9 or 10 phases have been used. On average 64% of the total number of counts have been used for the reconstruction of one 3D image.

Qualitative comparison

Each acquisition was reconstructed using the three methods NG3D, G4D and BH3D, that have been described in section 4.2.

The rst comparison was a visual comparison of the lesions in images reconstructed with the dierent methods. the outlines of the lesions are less blurred and the volumes are smaller than in NG3D image and closer to G4D end-of-exhalation image, while the lesion uptake is higher than in NG3D image with an equivalent noise level. -0.64% and -1.2%) (p = 0.9873 and p = 0.5608).

Volumes measured with BH3D method (see gure 4.7 (b)) were signicantly lower (-14.41%) than with NG3D (p < 0.001) and signicantly higher (+12.04%) than with G4D (p < 0.001), but were not signicantly dierent (-1.3%) to 6 phases-G4D (p = 0.655).

Lesion SNR on BH3D images was signicantly higher (+124% and +55%) (see gure 4.7 (c)) than on both 15 phases-and 6 phases-G4D images (p < 0.001) and closer to NG3D lesion SNR, while it was still signicantly smaller than on NG3D (-5.40%) (p = .0035). in the case of G4D with 15 phases (respectively 6 phases).

b Noise measure refers to the SD in a homogeneous liver area.

Discussion

The motivation of this work is tumour detection and quantication in oncology small animal SPECT imaging. Respiratory motion is a major source of impairment not only for lung lesions but also for abdominal lesions as in our peritoneal carcinosis model. The gating reconstruction method was able to reconstruct motionless images but at the cost of a degraded SNR. Therefore we proposed a retrospective amplitude-based data selection that aimed to maximise the data available for reconstruction, hence keeping a high SNR, while suppressing respiratory motion artefacts.

The proposed method was possible thanks to the breathing rhythm particular to anaesthetised mice, which is characterised by gasping breaths separated by long intervals without respiratory motion. The average duration of the motionless baseline is 9.75 phases (for a 15phases reconstruction), which corresponds to 65% of the respiratory cycle, which indicates that the same percentage of acquired counts may be used for the reconstruction of a motion free image.

We observed that the true motion, as observed in the CT or SPECT respiratory gated reconstruction (3D+t images), is not synchronised with the pressure signal, but appeared temporally shifted. This shift represents around a fourth of the cycle length (around 187 ms). It has been identied as a reconstruction in micro-SPECT The proposed image-based motion detection method is simpler than those proposed in e.g. [BMME + 08] or [START_REF] Paul | Retrospective data-driven respiratory gating for PET/CT[END_REF]. Indeed, both approaches aimed at reconstructing a respiratory representative signal from the images while we only focused on the detection of the phases exhibiting a motion with respect to baseline. More precisely, [BMME + 08] followed the centre of mass of a (manu- ally specied) tumour in a 4D series and then used the measured displacement as a trigger signal. This approach required some manual interaction and is dedicated to the reconstruction of one single tumour. The proposed method is closer to [START_REF] Paul | Retrospective data-driven respiratory gating for PET/CT[END_REF] since we study the variation of counts induced by an object motion. They proposed to integrate this variation in the sinogram to get a respiratory-like signal while we simply identied the largest variations at the voxel level. As with theirs, our method is versatile and can deal either with emission tomography or CT data. The proposed method appears to be quite robust with respect to noise, and is still able to detect the motion in noisy challenging images (SNR of 2.5), noisier than the pre-clinical images (see the SNRs in Table 4.3).

The average time shift between the image-based detected motion and the amplitude peak in the pressure signal was -3.875 phases (for 15 phases-G4D), which corresponds to 25% of the respiratory cycle. Rather than using a 4D CT image for the motion detection, that would imply a large acquisition time and a subsequent radiation dose, we rely on a 4D SPECT image. Although the operator tried to stabilize the respiratory cycle of the anaesthetised mice, some variation still occurred (Table 4.1). As a consequence, some counts may be skipped from the reconstruction of gated images and others will be used twice, mostly in the rst phases to be reconstructed. This motivated the equalization (done by resampling) of all respiratory cycles before reconstruction.

The proposed data selection allowed to retain 64% of all detected counts to reconstruct one single 3D image (BH3D). These counts correspond to the motionless baseline, hence the proposed reconstruction method simulated a breath-hold acquisition.

Using all acquired counts in the reconstruction would have required to either co-register all 3D phase image and then average them (e.g. [START_REF] Dawood | Lung motion correction on respiratory gated 3-D PET/CT images[END_REF]),

or to incorporate a motion model into the reconstruction process (e.g.

[RMK + 07]). This would then either rely on some non-linear registration method, whose validation for such noisy data is challenging, or on the calculation of a personalized respiratory motion model. 4D CT may oer a means for such computation, but will still depends on the used registration algorithm and at the cost of a high radiation dose of the mice. Therefore, the proposed method oers a means toward the reconstruction of a still image, while not requiring the validation of any registration method.

Measures of image noise (gure 4.7d) and lesion signal-to-noise ratio (gure 4.7c) demonstrated that BH3D is of better quality than a 3D phase image extracted from a 4D reconstruction (G4D, 6 and 15 phases reconstruction have been tested), and even suggested that it is of similar quality to the 3D image reconstructed with all the available counts (NG3D). Measures of SUV mean are signicantly larger in BH3D than in NG3D, and close to the ones from the phase image (G4D). Similarly, volume measures are smaller in BH3D than in NG3D, and close to G4D. These gures suggested that the blurring eect due to the respiratory motion has been compensated for. As expected, using less phases in G4D reconstruction yielded better quality images. The advantage of the proposed method is to discard prospectively the motion impaired data before reconstruction, while one has to identify retrospectively (possibly with a monitoring system) the motionless phases in G4D reconstruction. Last, the optimal G4D will reconstruct 2 phases, each of them from 50% of the data, reconstruction in micro-SPECT assuming only one of them is impaired with motion, while our method used 64% of the data.

More thorough validation could be conducted with simulations. First, one would have to animate a numerical mouse phantom with the same particular respiratory pattern that has been observed (gasps followed by a rest period).

Such a model has not been identied. The MOBY phantom [STF + 04] imple- mented a respiratory motion similar to humans, and would not be adequate.

Moreover, we have no indication whether such a model may produce the same time shift as we observe between the pressure signal and the abdominal motion. Second the acquisition can be simulated (e.g. with the GATE software

[SRK + 06]
), but this will also require to model our camera. Such a consid- erable amount of work is certainly worthwhile for a ner characterization of the benets of the proposed method. However, the presented gures already suggest there is an the overall benet.

Conclusion

We developed a reconstruction method dedicated to anaesthetised freebreathing mice. It enabled the reconstruction of a breath-hold like acquisition that is comparable to a non-gated reconstruction in terms of noise measure and signal-to-noise ratio, and intermediary between gated and non-gated reconstruction for lesion-based measurements (SUV peak , SUV mean and lesion volume).

Overall, the proposed method improves the quality of pre-clinical images, and the precision of the quantitative measurements they provide. This method is promising for more challenging studies concerning organs or tumour affected by motion using the emission tomography modalities in the pre-clinical research.

In addition to the results presented in this chapter, we provide complementary results in Appendix A. This appendix illustrates dierent reasons that led us to the approach presented in this chapter, and illustrates with additional tests the eciency of our approach.

Chapter 5

Dynamic SPECT analysis of 99m Tc-pertechnetate uptake in stomach: biodistribution study By studying 99m Tc-pertechnetate uptake kinetic, we mean to analyse spatially and temporally the 4D SPECT images. This requires the acquisition of dynamic SPECT, then extraction of information from the sequence. The rst objective is achieved in this chapter while the modelling aspect is addressed in Chapter 6.

Dynamic SPECT acquisition and analysis

Dynamic SPECT or 4D SPECT consists of a sequence of 3D images. The four dimensions refer here to the three spatial dimensions plus the time of image acquisition with respect to the time of the iodide analog 99m Tc-pertechnetate injection. During the acquisition sequence, the animals are maintained under anaesthesia at all time. Indeed, the length of acquisition could not exceed 140 min in order to prevent animal from premature death.

The temporal and spatial study of iodide uptake is a biodistribution study.

Biodistribution refers to both static and dynamic distribution of compounds within a biological system or within an organism. In this work, it is not the iodide biodistribution that is studied but the 99m Tc-pertechnetate biodistribution. As 99m Tc-pertechnetate is a substitute of iodide with respect to the NIS protein, we assume that 99m Tc-pertechnetate biodistribution reects the iodide transport mechanism.

99m Tc-pertechnetate is transported like iodide but it is not incorporated like iodide.

99m Tc-pertechnetate is thus privileged for the study of iodide dynamic distribution. Moreover, we assume that the transformation of iodide into organic compound does not occur in this gastric biodistribution. At rst sight, it consists of a voxel-based analysis: the activity measure in the dierent regions of the stomach along the time. This task could be done easily if the stomach was not deforming and if there was an anatomical reference at all times, given by computed tomography (CT) imaging. However, as previously mentioned in the rst part of this thesis, the abdominal organs are not distinguishable in non contrasted CT images.

Contrast agent could be used in order to see the stomach but this would require to alternate SPECT and CT acquisitions in order to have a 4D image for each modality. Such 4D CT images would have allowed the longitudinal registration of the 4D SPECT on the anatomical information. However, this kind of protocol makes more complex the acquisition and image analysis process. That is why, in a rst phase, we investigate only 4D SPECT images to study the biodistribution, without any anatomical reference for the position or shape of the stomach.

Challenges

The challenges here are to identify the two functional compartments of the stomach in the images, and to understand what is going on between these two compartments. In others words, the rst challenge of this study is to nd an appropriate segmentation method for organs in movement. The second challenge is to obtain qualitative data from this segmentation that allow the compartmental analysis.

The diculties that arise are that the stomach is moving according to the animal movement, respiration, and deformation induced by digestion. Although the animal is under anaesthesia, it can still slightly move. Thus, an appropriate automatic segmentation method should be able to integrate all of that.

Region-based analysis of the 99m Tc-pertechnetate biodistribution

The task consists in developing an appropriate method for stomach segmentation and activity tracking dedicated to 4D SPECT. The ideal approach would be a voxel-based tracking approach as the animals are not moving under anaesthesia during the whole acquisition. Such a voxel-based approach could be done without any additional registration for static object. However, as the stomach is subject to deformation along the time, images should be registered in order to compensate for the organ deformation. The problem with registration is that it is usually based on intensity. As no anatomical CT image for each frame of the dynamic SPECT was provided, the registration algorithm should be proceed directly on the SPECT images. In the case of a activity progressive uptake in the organ of interest, the voxel intensities in SPECT images change over the time and this might bias the outcome of the registration. Therefore, voxel-based tracking approach might not be appropriate. For example, a particular voxel that corresponds to the stomach wall in the rst image of the sequence might not correspond to the stomach wall for the next images. From this point of view, we considered a region-based approach for the extraction of the time activity curves. A set of time activity curve will be obtained with each curve associated to one region. In a rst time, the activity will be the average activity in the each region.

Biological context

The biologists study the iodide uptake and elimination processes in the body.

Iodide accumulates in dierent organs of the body, not only the thyroid but also salivary glands, or stomach. NIS protein mediates all this accumulation in these tissues. In this context, biologists from TIRO team focus their study on the role of iodide uptake in the extrathyroidal tissue.

Here, we focus on the murine stomach. Biologists conducted ex vivo investigations with immunohistochemistry (IHC) image. From these investigations, they identied NIS protein in glandular cells of the gastric wall. Since, NIS protein mediates the iodide uptake, biologists expect that iodide has a specic function in the stomach. However this hypothesis has to be conrmed by kinetic studies. These studies have to be conducted in vivo and SPECT imaging was proved to be the appropriate tool.

Biological objectives

The motivations for this study are biological. The main biological goal is the understanding of the iodide gastric cycle.

99m Tc-pertechnetate is a substitute of iodide with respect to NIS protein.

In the sense that NIS mediates The wall of the gastric stomach, the gastric mucosa, is the mucous membrane layer of the stomach which contains the glands and the gastric pits. The gastric epithelium is responsible for most digestive, absorptive and secretory processes. In human, gastric glands are organized in order to secrete a mucus, enzyme or acids that all contribute to the digestive process. The secreted mucus contains among other iodide secreted by NIS expressing cells in its epithelium. In comparison, the forestomach does not present NIS expressing cells in its epithelium and is used for food storage and digestion process. Biological investigations with murine IHC imaging highlight the fact that NIS proteins are expressed by cells from the innermost layer of glandular stomach called epithelium as we can see on Figure 5.2. This can be characterized by an important iodide uptake observed with the SPECT imaging.

Stomach SPECT imaging with 99m Tc-pertechnetate

The iodide biodistribution can be observed and studied in vivo with SPECT imaging of mice using the 99m Tc-pertechnetate that is an iodide substitute. In this sense, it is transported by same routes than iodide. First, it is mediated by NIS proteins and thus is accumulated and diused in the whole cellular space.

Then, in the case of stomach, the glandular cells that accumulate 99m Tcpertechnetate also secrete it in the gastric mucus. This 99m Tc-pertechnetate also diuses in the whole stomach cavity.

High resolution SPECT imaging with 99m Tc-pertechnetate allows now the observation of 99m Tc-pertechnetate distribution and allow the identication of both tissue and mucus as well as the quantication of the iodide uptake 99m Tc-pertechnetate uptake in stomach: biodistribution study capacity in tissue.

Material

Animal model

No specic mouse lineages were required for this work. Two common laboratory inbred strains (i.e. with pure genotype), C57Bl/6 and Balb/c, have been used. A set of fourteen 4D SPECT images were acquired with six to 24 frames. The description of these acquisitions is detailed in Table 5.1. Five animals were injected intra-peritoneally between 110 MBq and 180 MBq in 400 µL of 99m Tc-pertechnetate. The other nine animals were injected subcutaneously between 90 MBq and 110 MBq in 400 µL of 99m Tc-pertechnetate.

The radioactivity was obtained from a freshly eluted 99Mo/99mTc generator.

Then the animals anaesthetized with isourane (2-chloro-2-(diuoromethoxy)-1,1,1-triuoro-ethane)(1.3% v/v) (Baxter, France) and maintained asleep during the whole imaging protocol, using inhaled anaesthetic.

Animal housing and procedures were conducted according to the guidelines of the French Agriculture Ministry and were approved by the local ethics committee.

Dynamic SPECT imaging protocol

The acquisition process of dynamic SPECT images is very similar to the acquisition process described in Section 4.2. However, the temporal dimension is not at the respiratory cycle scale but at the scale of several dozen of minutes.

Thus, we set up a sequence of six to 24 acquisitions distributed along 80 min to 145 min. As we are not initially interested in respiratory gated acquisitions, images are reconstructed using the basic non-gated 3D reconstruction scheme (see Section 4.2).

During the whole protocol, animals were maintained asleep and their breathing rhythm was controlled with the respiratory monitoring system.

Besides the temporal setting of the protocol, the eld of view were carefully chosen in order to see the whole stomach. In most acquisitions (11 out of 14), the eld of view also contains the heart. However, the larger the eld of view, the longer the acquisition. The three dynamic acquisition with 24 frames were acquired without the heart. All protocols also include a CT image for anatomical reference, although it is not helping much with the stomach localization without contrast agent. It is principally used to localize the left ventricle. Indeed the heart and in particular the left ventricle is acquired to measure the blood activity function. This function is essential for the Three sequences include an additional SPECT image and a CT images acquired after contrast agent absorption. Food mixed with barium sulphate was given to mice after the dynamic SPECT imaging. This required to remove the mice from the bed, ll them with the barium food, to re-anaesthetized and to replace in the cradle. This step was done only after the radioactivity dynamic observation since the barium sulphate can aect the iodide uptake kinetic in stomach. These additional images were acquired to investigate the feasibility of stomach imaging combining both modalities SPECT and CT.

However, the impact of this barium sulphate ingestion is not further studied in this thesis since the images acquired afterwards exhibited large changes.

These changes would have required additional registration that has not been performed.

Method and qualitative results

In this section, we describe the method that we used to analyse the dynamic SPECT images in order to observe qualitatively the iodide biodistribution in the stomach. There are two challenges for this method that we had to consider:

1. the progressive uptake in the stomach wall coupled with a progressive diusion of the activity in the stomach cavity, and 2. the progressive deformation of the stomach due to the combination of digestion, air bubble displacement, and respiration.

Despite these diculties, the method is kept as simple as possible with several assumptions. The assumptions that have been done are:

The stomach has a bean shape. It can be decomposed into layers.

In a sequence of images, layers can be associated with respect to the distance to the border.

These assumptions are debatable and will be discussed in Section 5.5. Stomach segmentation An intermediary dynamic image is obtained after applying a Gaussian ltering with σ = 1.00 to the original image. This step is done in order to smooth irregularities that might aect the intensity-based segmentation. This ltered image will be used only for the stomach segmentation. Then the stomach is segmented by thresholding the intensity at the level determined by the Otsu method [START_REF] Otsu | A Threshold Selection Method from Gray-Level Histograms[END_REF]. The level is the threshold that minimises the standard deviation in each group. This level is computed independently for each image of the sequence. For each, a mask of the higher intensity is obtained. The largest connected component is kept. The mask corresponds to the stomach wall at early times. For later times, because of the diusion phenomenon in the cavity, the mask includes also the stomach cavity. In order to segment both wall and cavity, we apply a 3D morphological closing with the 3D extension of the 5x5 structuring element in Figure 5.4.

Layer-based decomposition of the stomach

For the temporal observation of the stomach activity a region-based analysis has been chosen. This approach can be also designated as a layer-based approach. The stomach is divided in to layers modelled on onion layers. Each 99m Tc-pertechnetate uptake in stomach: biodistribution study stomach region or layer consists of the set of voxels that are equidistant to the outside of the stomach (background). Then, the idea is to associate each layer along the time, starting from the outer layer, and then to observe the temporal activity variation for each layer.

0 1 1 1 1 √ 2 √ 2 √ 2 √ 2 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 1 √ 2 √ 2 √ 2 √ 2 √ 3 √ 3 √ 3 √ 3 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 √ 5 (
The distance transformations of the segmented stomach was computed with a 3D neighbourhood mask. According to [START_REF] Borgefors | Distance transformations in digital images[END_REF], the optimal distances in 2D digital images compared to the Euclidean distance are obtained with a 5x5 neighbourhood mask denoted Chamfer 5-7-11. We used the 3D extension with respect to the background. The following layers (toward the inside of the stomach) are located at a distance that is an optimal approximation to the Euclidean distance according to [START_REF] Borgefors | Distance transformations in digital images[END_REF]. Indeed, after normalisation, the 2 nd layer is located at the distance 7/5 = 1.4 from the stomach border which is a good estimate for √ 2 1.41. The 3 rd layer is located at the distance 9/5 = 1.8 from the stomach border which is a good estimate for √ 3 1.73. Layers are associated with respect to their distance to the background. An example is given for a layer that corresponds to the inner layer of the stomach wall at distance 2.4 in Figure 5.7(b). The assumption that are made concerning the layer matching will also be discussed in the in Section 5.5.

Then, it is possible to compare all time-activity curves of all layers in order to observe some trends that could suggest some specic kinetic according to the distance. Figure 5.7(c) shows an example for the same subject. The analysis is dicult since there are too many layers. However, it seems that the activity corresponding to stomach wall increases faster at the beginning but then slower, while the activity in the stomach cavity layers is almost nonexistent but then increases faster up to an equivalent or greater average activity than in the wall.

Results

Stomach segmentation

Here, we present the results for the segmentation process for some interesting cases. For each case, we present one transverse slice centered in the stomach cavity where the original image in gray level is superimposed by the resulting segmentation (green line). Generally, the segmentation performed well for all subjects, segmenting what was expected. However, we notice some interesting results.

The segmentation for subject #0366 performed well as presented in Figure 5.8. This acquisition present clearly both challenges: a large stomach deformation and a progressive radioactivity uptake along the time. This large deformation results in a large variation in stomach size. Figure 5.9 shows the resulting segmentation for the subject #0370 performed well. The process went well although at T22, 111.5 min after injection the stomach shape abruptly changed. This was due by the fact that the mouse woke up during the acquisition sequence and moved. This induced a spatial shift but the stomach was still segmented and the stomach volume remained coherent. In the end, the method is eective in the slices where the stomach appear like a balloon. This shows a progressive activity uptake in the cavity that lled up to the margo plicatus (opening on the upper right corner). In the end, the method is less eective in the slices where the stomach appear like a bucket.

Despite the numerous segmentations that performed well, the process failed for some images acquired at early time after injection (less than 10 min). The activity uptake in the stomach was not high enough to compete with the low signal to noise ratio, or with the presence of the site of subcuta- neous injection in the stomach eld of view. However, after several minutes, the activity uptake in the stomach wall increased while it decreased in the site of subcutaneous injection, and the stomach was segmented again. This diculty occurs with all subjects #0396, #0397, #0398, #0402, #0420 that were injected subcutaneously and imaged at short after the injection. As an example, in Figure 5.11, the stomach wall can be distinguished at any on the selected slice by adjusting the intensity window but for the rst instant, it was not segmented. In Figure 5.12, it is clear that the presence of remaining activity in the site of injection aects the segmentation as it dened.

Finally, Figure 5.13 shows the temporal variation of the stomach volume for some cases that have been described. As expected, for the subject #0366 (see Figure 5.13(a)), the large variation in the stomach segmentation yields to a large variation in the volume as well. Figure 5.13(b) shows that the stomach volume slowly but regularly decreases with respect to the time. Only 99m Tc-pertechnetate uptake in stomach: biodistribution study the third to last volume measures presents a small discrepancy comparing to the other, but it corresponds to the time when the mouse woke up and moved.

Figure 5.13(c) is similar to the previous one, except for a small variation in the rst instants that occurs very early. However, the rest of the variation shows a slow decrease of the volume. Figure 5.13(d) shows that the stomach volume is much bigger for the rst acquisition than after. Actually, this is due to the fact that the volume corresponds here to the volume of the segmented activity.

And the segmented activity for the rst time point was the subcutaneous injection site instead of the stomach itself. The same goes for the two rst instants for the subject #0402 Figure 5.13(e), although this time the wrong segmented region was smaller. This rst volume (resp. two rst volumes) is (are) then not relevant and should not be considered for the following steps.

For the next volume measurements, volume slowly but regularly decreases with respect to the time, as it is shown in Figure 5.13(b). with respect to the stomach original image and its segmentation for the subject #0370. On the original image slices presented in Figure 5.14(a), we can see that the intensity window enlarge with the time, suggesting that the stomach activity uptake increases proportionally. At the same time, it seems that the stomach is slightly shrinking as suggested by the Chamfer distance 3D map presented in Figure 5.14(c). To be mentioned, the 5 th image corresponding to an acquisition 111.5 min after injection of 99m Tc-pertechnetate reveals that the animal motion during the sequence of acquisition. Still, despite the animal motion, the layer-based decomposition allows the tracking of the progressive 99m Tc-pertechnetate uptake in the stomach cavity.

Temporal and spatial activity curve

Spatial activity curves Figure 5.15 shows activity curves at dierent times. Each subgure shows the activity depending on the distance from the border of the stomach. More precisely, each point of the curve corresponds to the average activity observed on the layer dened by its distance to the border. Boxplots represent the range of activity in the layer, they show maximum and minimum, and rst and third quartile. We rst observe that the activity is more and more important along the acquisition but also that the activity does not evolve in the same way in all the layers. In Figure 5.15(a), at 8.73 min, the activity is quite low in all the layers. However, the external layers (approximately between 0 and 4.5 voxels) expose a higher activity, we assume that these layers correspond to the stomach wall.

The internal layers that have a lower activity correspond to the stomach cavity.

The Figure 5.15(b) and Figure 5.15(c) show the activity after 45.73 min and 82.75 min. At these instants, the activity is constantly growing, and the stomach wall still has more activity than the stomach cavity. Note that the point where the activity decreases is always around 4.5 voxels, this corroborates our hypothesis that the frontier between stomach wall and cavity is around this position.

In Figure 5.15(d) and Figure 5.15(e), at 112.62 and 134.92 min, the growth of the activity slows down, and the activity of the stomach wall stabilizes.

However, the activity of the stomach cavity still slowly increases and nally reaches approximately the same activity as the inside of the stomach wall.

To summarize, along time the activity shown by the stomach increases: at the beginning this activity is accumulated mostly in the stomach wall, and then, after a delay, the activity diuses towards the stomach cavity. global activity of the whole stomach, that has the same global trend. This corroborates the precise analysis that we made in the previous paragraph. always have a small activity may be lost in the segmentation after some time.

Temporal activity curves In

However, we showed in Figure 5.13(c), Section 5.4.4.1, that for this same mouse, the variation of the stomach volume was very low, but decreasing.

The conclusion is that the automatic threshold detection is not perfect and most probably looses some of the external layers of the stomach along time.

However, we consider in the following that, as the volume variation is quite low, the quality of the automatic segmentation is sucient enough for our study. In addition to time activity curves for stomach layers, we extracted the blood time activity curves for all acquisitions that have the heart in the eld of view. This time activity curve will be used to estimate the input function in the compartmental analysis that we will present in the next Chapter. They were obtained by measuring the average activity in a 3x3x3 mm 3 sphere lo- cated in the left ventricule. Figure 5.18 shows the time-activity curve (TAC) for all subjects, for each type of injection: subcutaneous (subgure (a)) and intraperitoneal (subgure (b)). As the injection can be considered as an impulse, the blood TAC is expected to be a bolus. As we can see for the subject #0530, in the subcutaneous case, the activity increases quickly, reaches its maximum around 20 min after injection and then decreases slowly. This shape can be considered as closed to a bolus. In the intraperitoneal case, the curve is only decreasing and if there was a peak it was in the rst instants, before the rst measure. We would probably need more data right after the injection in order to observe a peak and use this kind of injection method.

We choose to perform our compartmental analysis on the subject #0530 and with the subcutaneous injection because we have a very precise sampling in this case, and the blood temporal activity curve has the expected shape, i.e., it is closed to a bolus.

Discussion

The main contribution in this chapter is a region-based analysis for the study of iodide biodistribution in SPECT images of mice stomach. This simple approach allows the observation and the description of the progressive 99m Tcpertechnetate uptake in the stomach. We also described the challenges that had to be addressed in our analysis: progressive radioactive uptake and stomach deformation. The region-based approach has been chosen instead of a voxel-based. Voxel-based analysis would have allowed a much more precise study if the stomach was not deforming, but the movement of the stomach and the evolution of the activity makes this approach impossible. Our regionbased approach is based on layers, it relies on the following observations:

The stomach has a bean shape. It can be decomposed into layers.

In a sequence of images, layers can be associated with respect to the distance to the border.

The strengths of this layer-based approach are the following:

Our approach is adapted to stomach deformation.

The Chamfer distance transform is simple, ecient, and provides the best approximate for Euclidean distance.

The layer-based approach makes up for the irregularities of the stomach shape.

However, the region-based approach has some weaknesses induced by some assumptions and choices that were done. In the following of this section, we review these weaknesses, explain them, and describe how we could improve the method in the future.

Stomach segmentation short after the injection Intensity thresholding with Otsu method nds the optimal threshold that divides the data into two groups of voxel by minimising the standard deviation in both groups. As it is only based on the intensity, the process failed when the intensity in the object (e.g. the stomach wall) was too low compared to the noise or the site of injection. In practice, in these rst instants, it is often the injection site that is segmented. The consequence of this wrong segmentation is that all the activity measured in the layer will lead to wrong interpretation and aect further work such as the extraction of the time activity curves for the compartmental analysis that will be presented in Chapter 6. A solution for this would be to guide and constrain the segmentation in early acquisition using the 99m Tc-pertechnetate uptake in stomach: biodistribution study result of the segmentation at later instants. This way, the region segmented in the rst instants would be more precise and lead to a better analysis.

Approximation due to the Gaussian ltering The fact that we apply a Gaussian lter smooths the irregularities of the stomach but it also tends to enlarge the segmented area for the stomach, especially when the activity is low (i.e., short after the injection). However, we consider that the benets of the ltering are more important than this slight imprecision.

Stomach volume decrease As mentioned in Section 5.4.4.3, when the activity increases, the size of the stomach decreases because the external layers of the stomach may be lost in the automatic segmentation process. The main consequence is a slight layer mismatch in our analysis:

the rst layer at a given instant might disappear in the next instant.

However, we have shown that the decrease of the stomach volume is quite low, and we consider that this layer mismatch is negligible for our study. In the future, more complex segmentation techniques might be considered in order to limit the layer mismatch.

Morphological closing is not adapted to the shape of the stomach.

The stomach of the mouse is in fact made of two parts: the glandular wall and the cavity. The margo plicatus forms the limit between the two areas; it behaves a bit like a valve which prevent the food to go back to the forestomach. In the forestomach and in the margo plicatus, there is few or no NIS expressing cells. Consequently, the shape of the area of the stomach is not fully closed by the stomach wall. Indeed the region corresponding to the margo plicatus is very thin, and the layer approach will consider the margo plicatus and the part of the stomach cavity nearby as part of the stomach wall. However, this region of the stomach is quite small and we consider that the imprecision induced by this wrong classication as negligible. Detecting more precisely this area would be very complex and would probably not improve the quality of our study.

Chapter 6. Dynamic SPECT analysis of 99m Tc-pertechnetate uptake in stomach: biodistribution modelling

Introduction

In the previous chapter, we dened the biological context of this study and we showed that it is possible to quantify the 99m Tc-pertechnetate uptake in the murine stomach with dynamic single photon emission computed tomography (SPECT). The working hypotheses in terms of biology are: the 99m Tc-pertechnetate uptake in stomach is cyclic with a reabsoption phenomenon in the pars pylorica region of the stomach (see Figure 5.3)

there is a 99m Tc-pertechnetate retention in the glandular cells of the gastric mucosa (see Figure 5.2)

These hypotheses are done according to observations done on stomach immunohistochemistry (IHC) images. Knowing all of this, we do want to investigate further the mechanisms of this uptake. We would like to set up a model that explains these observed 99m Tc-pertechnetate dynamic uptake in the stomach. Then, we would like to estimate the model parameters that best t the observations. Among others, with this descriptive model, we would like to conrm that the sodium/iodide symporter (NIS) expression in the apical membrane of pars piloric cells contribute to the cyclic 99m Tc-pertechnetate uptake mechanism. Moreover, we would like to answer the biologists main question: Is the 99m Tc-pertechnetate accumulation in the gastric glandular cells superior to the secretion into the stomach cavity? and if yes, what is the limiting step? We remind that the role of gastric cells expressing NIS on their basolateral membrane is the accumulate and the secretion of 99m Tcpertechnetate. However, the functional role of gastric cells expressing NIS on their apical membrane has to be conrmed: Does it maintain or control the 99m Tc-pertechnetate concentration in the cavity?

Observations Indeed, the 99m Tc-pertechnetate uptake from the blood to the stomach wall is clearly associated to the expression of NIS proteins in the glandular cells. Then, it is also known that the 99m Tc-pertechnetate is also secreted by these glandular cells in the stomach cavity. Both of these observations have been assessed by SPECT imaging. However, since there is no NIS protein expression in the apical regions of some stomach wall cells, there is no proof of reabsoption of 99m Tc-pertechnetate by the stomach. Still, dynamic SPECT analysis showed that, after some time, the activity in both stomach wall and cavity tends to stabilise (i.e. the speed of activity uptake is slightly reduced). This observation suggested that there might be a 99m Tcpertechnetate reabsoption phenomenon by the stomach wall.

Objectives Consequently, the objective of this chapter is rst to dene a descriptive model of the 99m Tc-pertechnetate uptake in the stomach, then to identify the parameters of this model using the observations and nally to conrm the existence of this reabsoption according to the parameters that have been observed.

Dierent types of modelling The general objective of modelling is to nd some general laws that will describe a phenomenon. The model is good if it the real behaviour exactly follows the laws that have been discovered.

However, reality is often too complex, then simplied models provide a precise enough approximation of the phenomenon under certain assumptions. One could classify modelling in two main categories according to their objectives: the descriptive modelling that aims at describing a phenomenon according to some observations. It is used to explain how and why the studied phenomenon occurs. The eld of application of such models is wide, they are frequently used for understanding consequences or behaviour of some new products (pharmaceuticals).

the predictive modelling goes a step further. It aims at deducing laws that not only describe the present but also predict what will happen in new situations according to a wide set of data.

In this study, our primary objective is to conrm the existence of a phenomenon. This should be addressed by descriptive modelling, based on the data we presented in the previous chapter. The data we currently have and the current status of our analysis prevents us from doing any predictive analysis at this point. However, a possible follow-up of our work could be to use the model we design in order to predict the value of the parameters of the phenomenon for a new subject, or in a dierent context. This inverse problem can be solved in two dierent ways: a direct approach that uses a semi-analytic resolution, or an indirect approach that estimate the set of transfer parameters with optimisation approaches.

Compartmental analysis

The purposes and applications for compartmental analysis are multiple.

In the literature, a few biodistribution works have been conducted with compartmental analysis approaches on other organs than stomach using ET imag- Our approach In this work, we will focus on the study of the kinetics of one substance, the iodide analog 99m Tc-pertechnetate, and in one single organ -the stomach -since the general idea is to design a model for iodide uptake kinetic in the stomach.

This chapter presents the following contributions:

A compartmental model of the stomach (Section 6.2) that distinguishes the stomach wall from the stomach cavity;

An analytic and numerical resolution methodology for determining the time-activity curves (TACs) for each compartment of the stomach (Section 6.3) based on the activity of each layer obtained as described in the previous chapter;

The results: the parameters that characterise 99m Tc-pertechnetate uptake in the stomach such that those parameters t with the observations made in the previous chapter (Section 6.4). 99m Tc-pertechnetate biodistribution in mice, proposed in [START_REF] Lathrop | Biologic behavior of 99m Tc from 99m Tc-pertechnetate ion[END_REF] biodistribution in mice. For example, Figure 6.1 shows a schematic representation of the 99m Tc-pertechnetate biodistribution proposed in [START_REF] Lathrop | Biologic behavior of 99m Tc from 99m Tc-pertechnetate ion[END_REF]. This representation also gives temporal estimation of uptake or secretion of 99m Tcpertechnetate after an intravenous injection. However 99m Tc-pertechnetate uptake and kinetic in stomach has been mentioned in a few papers [START_REF] Ng Hofmeyr | Stomach scanning after intravenous 99mTc administration. A preliminary report[END_REF].

Compartmental model of the stomach

The diculty with the stomach is that it is a soft organ, with non predictable deformation due to alimentary bolus combined with respiratory motion and subject movement.

In this section, we will a schematic representation of the 99m Tcpertechnetate biodistribution in the stomach that corresponds to the block gastric mucosa + juice in Figure 6.1. The gastric juice is a digestive uid that is secreted by the gastric mucosa, i.e. various glands in the mucous membrane of the stomach (i.e. the gastric mucosa). This representation will be used for the kinetic study with compartmental analysis. kinetic of iodide uptake in the murine stomach. It is focused on the organ of interest: the stomach. Each box refers to a compartment, i.e. an organ or a set of organs, which is a macroscopic subunit containing a concentration x i of 99m Tc-pertechnetate. The observations (or measures) of this activity are the system outputs of the model, and are denoted y i . Arrows between boxes represent the transport of 99m Tc-pertechnetate between compartments, they are labeled with transfer parameters k ij ; some arrows have no label because they will not be used in the simplied model that we will describe in the next section.

Q Q Q Q Q Q Q Q Q Q s Q Q Q Q Q Q Q Q Q k k Dig2 k 2Dig k 0Dig k Dig0

Complete model

The input function u 0 (t) refers to the 99m Tc-pertechnetate injection in the mouse. Two types of injection can be performed: an impulse (bolus) that is a one-time injection of all the activity in the site of injection before the imaging sequence, or a continuous injection that is a permanent injection during the imaging sequence. In this study, the 99m Tc-pertechnetate injection was always an impulse. It is formally dened by u 0 (t) = D 0 δ(t).

Then, according to the site of injection, the blood activity function x 0 (t) is dierent. This blood activity function is experimentally measured in the left ventricle of the heart. In the case of an intravenous injection, all the activity is directly injected into the blood system, which generally induces the blood activity function to be of the form x 0 (t) = D 0 e -kt . However, intravenous injection are dicult to perform. There are two other possible sites for injection:

in the peritoneum (intraperitoneal injection), or under the skin (subcutaneous injection). For these two cases, the blood activity function can be modelled as of a sum of two exponential functions such as: x 0 (t) = φ 1 e -θ 1 t + φ 2 e -θ 2 t . The rst term represents the slow activity uptake in the blood system and the second term represents an even slower distribution of the activity through all the body. However, the kinetic is slightly dierent, since in the case of intraperitoneal injection, there is an activity diusion phenomenon in the peritoneum simultaneously to the progressive activity uptake by the blood system.

Once in the blood system, 99m Tc-pertechnetate is carried through all the body and is going to accumulate in several glandular organs such as the thyroid, salivary glands, stomach, or to be eliminated by kidneys.

As we focus on the interaction at the stomach level, Figure 6.2 details the exchange between the blood, the stomach wall, and the stomach cavity.

From the stomach wall point of view, 99m Tc-pertechnetate can be exchanged only either with the blood or with the stomach cavity. From the stomach cavity point of view, 99m Tc-pertechnetate is mainly secreted by the stomach wall, but exchanges might exist with the rest of the digestive track. Indeed, there is a progressive diusion of 99m Tc-pertechnetate in the cavity up to the forestomach cavity. Some 99m Tc-pertechnetate might be transported from the salivary glands to the stomach cavity in the saliva of from the stomach cavity to the duodenum at the same time as the alimentary bolus. 99m Tc-pertechnetate uptake in stomach: biodistribution modelling since the space that is occupied by blood is not a closed space and is hard to dene. As a consequence, the blood activity function cannot be deduced from the knowledge of the type of 99m Tc-pertechnetate injection in the body since the bioavailability fraction of activity in the blood is not known beforehand. Fortunately, the blood activity function can be easily measured in the left ventricle in SPECT images. Thus, we simplify our model by removing the two rst blocks: input function u 0 and the injection site. We now consider the input function of our compartmental model as equal to the blood activity x 0 as shown on the left of Figure 6.3. This is not an impulse function anymore. It can be approximated by sum of two exponential functions, as we will show in Section 6.3.

Simplied compartmental model

Diusion in the stomach cavity In the previous chapter, it has been observed that 99m Tc-pertechnetate is diusing in the stomach cavity from the stomach wall toward the middle of the cavity. Here, we rst put aside the diusion phenomenon inside the cavity, it might be studied in later improvements of the model. Consequently, in this rst study, we won't consider the diusion phenomenon in the cavity, but we will consider the stomach cavity as a whole.

99m Tc-pertechnetate in other organs and in digestive system We assume that the blood is carrying 99m Tc-pertechnetate through the whole body and that 99m Tc-pertechnetate can accumulate in organs that are deserved by blood. If it is accumulated in those organs, we assume that it cannot leave them, it is trapped like the model of the thyroid trap( [START_REF] Grant T Gullberg | Dynamic single photon emission computed tomography-basic principles and cardiac applications[END_REF]).

Then, the only glandular organ of the digestive system located before the stomach is the salivary glands, however according to Figure 6.1

(and [START_REF] Lathrop | Biologic behavior of 99m Tc from 99m Tc-pertechnetate ion[END_REF]), these glands only contribute to 5% of the activity after one hour. Additionally, in our case the mouse is sleeping, thus is assumed to produce less saliva that would go to the stomach. We thus consider that the 99m Tc-pertechnetate transport from the digestive system directly into the stomach is negligible.

In other words, from the stomach point of view, the 99m Tc-pertechnetate uptake can be considered as an autonomous cyclic phenomenon. As a consequence, we assume that there is no 99m Tc-pertechnetate transport parameters from the other organs into the blood. Additionally, we assume that there is no transport from the digestive system to the stomach cavity, and also that there might be leaks from stomach wall and cavity to the outside of the system, with transport parameters denoted k 01 and k 02 .

Transport between blood and stomach cavity As there is no direct blood injection into the stomach cavity, we remove from the model the transfer parameter k 20 of 99m Tc-pertechnetate from blood to the stomach cavity. However, we still consider there is a 99m Tc-pertechnetate transport k ext2 from the cavity to the outside, until it is re-transfered into blood. It is a non direct transport to the blood. In the next, by abuse of notation, we will denote it k 02 .

Mathematical representation of the model

The nal simplied model presented in Figure 6.3 can be studied using compartmental analysis similarly to the two-compartment model of Section B.

           ẋ1 (t) = -(k 21 + k 01 )x 1 (t) + k 12 x 2 (t) + k 10 x 0 (t) ẋ2 (t) = k 21 x 1 (t) -(k 12 + k 02 )x 2 (t) y 1 (t) = x 1 (t) y 2 (t) = x 2 (t) x 0 (t) = φ 1 e -θ 1 t
+ φ 2 e -θ 2 t (6.1)

In the Laplace equations, the input matrix U(s) is now expressed as follows:

U(s) = X 0 (s) 0 = φ 1 s+θ 1 + φ 2 s+θ 2 0 (6.2)
The transfer matrix G(s) remains unchanged (see Section B.3) but not the observation matrix Y(s) which depends on U(s).

Y(s) = G(s)U(s) = (s+(k 12 +k 02 ))k 10 (s-λ 1 )(s-λ 2 ) X 0 (s) k 21 k 10 (s-λ 1 )(s-λ 2 ) X 0 (s) = (s+(k 12 +k 02 ))k 10 (s-λ 1 )(s-λ 2 ) ( φ 1 s+θ 1 + φ 2 s+θ 2 ) k 21 k 10 (s-λ 1 )(s-λ 2 ) ( φ 1 s+θ 1 + φ 2 s+θ 2 ) (6.
3)

The denominator of the two terms of the observation matrix Y(s) are polynomials of degree 4 in s. Thus, each term can be decomposed in a sum of four rational fraction with denominator of degree 1 in s. Then, the inverse 99m Tc-pertechnetate uptake in stomach: biodistribution modelling Laplace transform ( L -1 ) of each term results in a sum of four exponential functions such as:

y 1 (t) = α 1,1 e -β 1,1 t + α 1,2 e -β 1,2 t + α 1,3 e -β 1,3 t + α 1,4 e -β 1,4 t y 2 (t) = α 2,1 e -β 2,1 t + α 2,2 e -β
α i , β i , φ 1 , φ 2 , θ 1 , θ 2 .
However, the analytic equations are too complex to be shown here. In practice, the system is solved both analytically and numerically with a mathematical software such as Maple.

Hypotheses according to the available data

The simplied model is a two compartmental model, with one compartment for the stomach wall and one for the stomach cavity. Our compartmental analysis relies on the following hypotheses:

Each SPECT image voxel contributes to each compartment with a percentage that can be calculated The rst challenge is to extract TACs for each compartment: the two compartments are physically close, the stomach wall is surrounding the stomach cavity. Thus the extraction of TAC for each compartment is hard since each voxel is a mix of each compartments and blood. Indeed, a SPECT image voxel measures 0.33×0.33×0.33 mm 3 which is more or less the same scale as the thickness of the stomach wall, thus most voxels contribute partly to the stomach wall and partly to the stomach cavity. The challenge here is to discover what percentage of each voxel corresponds to the wall and to the cavity, we call this percentage the mixing coecient. These coecients allows for transforming the layer-based activity measure into TACs for each compartment.

We can consider that the activity inside a compartment is homogeneous One of the hypotheses of compartmental analysis is that each compartment is homogeneous. Unfortunately, neither the stomach cavity nor the stomach wall is homogeneous. There are a progressive activity uptake in stomach wall coupled with secretion in the mucus and then a diusion phenomenon and all of these happen along the dierent layers. Consequently, the simplication hypothesis we do here is the homogeneity of each compartment; we consider the activity of the wall as the sum of the activities of the layers forming the wall.

The input is approximated by the activity of the left ventricle, expressed as a sum of two exponential functions We approximate the input function by a sum of two exponential functions; the TAC for this input function is estimated from the measure of the average activity in the left ventricle (see Section 5.4.4.4).

Each TAC can be approximated by a multiexponential function We suppose also that each TAC can be approximated by a sum of exponential functions.

Problem statement

We described in Section 6.2.2 the compartmental model that describes the 99m Tc-pertechnetate uptake kinetics in the murine stomach. The objective now is to determine the transfer parameters k ij of the model. These transfer parameters k ij represent the speed at which the substance ows from compartment j to compartment i, or from the outside (j = 0) to compartment i. In our case, these parameters are the two transfer parameters k 21 and k 12 between the two compartments, stomach wall and stomach cavity, the two transfer parameters k 01 and k 02 from each compartment to the outside, and k 10 the transfer parameter from the input to the rst compartment.

Knowing the TAC in the left ventricle, we will see in the next section that it is possible to t the input function with a sum of two exponential functions. In that case, we have shown in the previous section that the observations of each compartments should be of the form of a sum of four exponential functions.

There are two dierent ways to determine these parameters. Both ways consists in determining the kinetic model (i.e. the k ij constants) of a compartmental system from its observations (y).

A rst way consists in solving the inverse problem analytically [RW71, GRS + 10]. As we already have the input function expressed as a sum of two exponential functions, if we can approximate the output functions as a sum of four exponential functions, then the transfer parameters will be given by the analytic solution to the inverse problem.

A second way consists in solving the inverse problem numerically. Knowing the input functions, we test dierent sets of values for the transfer parameters, then we obtain the output functions of the system, i.e. compartment functions, as sum of (four) exponential functions. Finally, we compare these computed functions to the compartment observation and we optimise iteratively the transfer parameter in order to minimise the dierences between the computed and the observed compartment activities.

In the next sections, we will build the TACs for the dierent observations and try both the analytic and the numerical resolution methods to characterise the parameters of the system.

99m Tc-pertechnetate uptake in stomach: biodistribution modelling

Time-activity curve (TAC) acquisitions

First of all, before any other work, we now need to obtain the time-activity curves (TACs) for each of the involved compartments, i.e., the blood, the stomach wall, and the stomach cavity. These TACs are obtained by measuring the activity in volume of interests (VOIs) of dynamic SPECT images. In our case, we work with dynamic images that have been corrected for the radioactive decay and the cylinder factor, normalised with respect to both the animal weight and the injected radioactivity. Therefore the values of the image after conversion have no unit. By abuse of notation, we express these activity value in standardized uptake value (SUV) as it is usually done for positron emission tomography (PET) imaging.

Blood: model input function

First the blood time-activity curve is the measure of the activity variation in the left ventricle of the dynamic SPECT images. This curve is easily obtained.

Then, for the analytic resolution of the compartmental system, it is better to provide not only a time-activity curve but its approximation by functions such as impulse or exponential functions. However, the measured time-activity curve is neither a simple impulse nor a single exponential function. Though, it can be approximated easily by a sum of two exponential functions.

Approximation with sum of exponential functions

One assumes that the observations y of compartmental systems can be approached by continuous functions, that are sum of exponential functions [START_REF] Foss | A method of exponential curve tting by numerical integration[END_REF][START_REF] Jl Fresen | A note on Foss's method of obtaining initial estimates for exponential curve tting by numerical integration[END_REF]. Here we shortly describe how to approximate the set of observations, e.g. the blood time-activity curve u 0 , with a sum of two exponential functions y 0 (t) such as:

y 0 (t) = φ 1 e -θ 1 t + φ 2 e -θ 2 t
The idea is rst to dierentiate the theoretical expression of y 0 (t) twice and then re-integrate twice in order to obtain a formal representation such as:

y 0 (t n ) = a + bt n + cF (t n ) + dG(t n ) where                a = φ 1 + φ 2 b = φ 2 θ 1 + φ 1 θ 2 c = -(θ 1 + θ 2 ) d = -(θ 1 * θ 2 ) F (t n ) = tn 0 y 0 (λ)dλ G(t n ) = tn 0 λ 0 y 0 (ξ)dξdλ = tn 0 F (λ)dλ
Then, it is possible to estimate F (t n ), G(t n ), and y 0 (t n ) by integrating by parts on each segment delimited by the observations. The integration on each segment can be done according to the trapezoidal rule. This results in the following expressions. b, c, andd with least square minimisation algorithm between the observation y 0 and its estimate ŷ0 . Finally, it is possible to retrospectively determine the unknown constants φ 1 , φ 2 , θ 1 , and θ 2 that are combinations of the previously determined constants a, b, c and d. 99m Tc-pertechnetate uptake in stomach: biodistribution modelling

     F (t n ) = 1 2 n-1 k=0 (t k+1 -t k )(y 0 (t k+1 ) + y 0 (t k )) Ĝ(t n ) = 1 2 n-1 k=0 (t k+1 -t k )(F (t k+1 ) + F (t k )) ŷ0 (t n ) = a + bt n + c F (t n ) + d Ĝ(t n )

Stomach wall and cavity: model output

As mentioned in Section 6.2, if our TAC input model can be approximated by a sum of two exponential functions, we expect the observations in the two compartments of our model to be the sum of four exponential functions. In this section, we explain how we obtain the TAC for each compartment. These

TACs are expected to be approximated by a sum of four exponential functions using a similar approach to the one used for the blood function, but adapted to four-exponential functions. The rst goal is to obtain the TAC of the two compartments of the stomach from the activity in stomach layers. We have at our disposal the TAC for each layer of the stomach as presented in Figure 6.5. Finding TACs for the two compartments from these observations is challenging. On the one hand, compartmental analysis requires one TAC per studied compartment. On the other hand, we have a set of TACs, one per layer of the segmented stomach.

Principles

However, each layer is not homogeneous but a mix of compartment. It is assumed that each layer composition is time invariant. The question is then:

how much each compartment contribute to each layer?

Then, we assume that the activity concentration in each layer can be expressed as a linear combination of the activity concentration in the two compartments. Consequently, in order to build the TAC of the two compartments of the stomach, we need to determine which percentage of each layer belongs to the stomach wall, or the stomach cavity. We call these percentages mixing coecients. Those mixing coecient are dened as follows. Let m be the number of layers of the stomach model considered, and n the number of time frames. We suppose that the activity follows the rule:

L = K • C
where L denotes the activity concentration in layers, expressed as a m × n matrix;

K denotes the time-invariant mixing coecients, is it as m × 2 matrix and C denotes the activity concentration, expressed as a 2 × n matrix in compartments.

In other terms, we have:

     a l,1 (t) a l,1 (t 2 ) • • • a l,1 (t n ) a l,2 (t) a l,2 (t 2 ) • • • a l,2 (t n ) . . . a l,m (t) a l,m (t 2 ) • • • a l,m (t n )      =      κ 1,1 κ 1,2 κ 2,1 κ 2,2 . . . . . . κ m,1 κ m,2      • a c,1 (t) a c,1 (t 2 ) • • • a c,1 (t n ) a c,2 (t) a c,2 (t 2 ) • • • a c,2 (t n )
Assuming K is known, we can deduce the TAC for the two compartments by computing the pseudo-inverse of K [START_REF] Penrose | A generalized inverse for matrices[END_REF][START_REF] Sondermann | Stalistische Hefte Statistical Papers Best approximate solutions to matrix equations under rank restrictions[END_REF]:

C = ( t K • K) -1 • t K • L
The pseudo-inverse computes the best solution to a system of linear equations.

Here, we compute the concentration in each compartment that produces activity curves in each layer as close as possible to the observed ones 1 . The best t produced by the pseudo-inverse is in terms of least-squares minimisation. This solution will provide us the TACs for the two compartments of the stomach.

Identication of mixing coecients

The idea described above is pretty simple, however in practice, the diculty lies in the determination of the mixing coecient κ i,j .

First, we assume that the whole Stomach (S) is a combination of four dierent compartments identied by their 99m Tc-pertechnetate uptake properties: 1 The solution cannot be exact as we have m equations for only two parameters.

99m Tc-pertechnetate uptake in stomach: biodistribution modelling the gastric mucosa or NIS expressing cells, the gastric juice, the extra-cellular matrix, and the non diusible. The gastric mucosa is mainly composed by NIS expressing cells. By abuse of language, this compartment will be called stomach Wall (W). The gastric juice result from the secretion of the glandular cells. By abuse of language, it will be called stomach Cavity (C). The extra-cellular matrix is assumed to be of the same scale of magnitude in terms of activity than the blood. By abuse of language, this compartment will be called Blood (B). The Non Diusible is composed by everything that does not uptake 99m Tc-pertechnetate such as non glandular gastric wall and alimentary bolus. This compartment occupies most of the space, but as it does not uptake 99m Tc-pertechnetate, the contribution in terms of activity is null.

This can be represented analytically with:

S = W ∪ C ∪ B ∪ ND
In a rst time, we assume that the volume ratio of each compartment composing the layer can be used as the activity ratio for each compartment that contributes. In other words, we suppose that the activity in the layer is the combination of the activity of each compartment multiplied by a ratio that is the proportion of the volume of the layer that belongs to the compartment. More formally, denoting V W , V C , V B , and V N D the volume of the compartments and V l the volume of the layer, the activity of the layer l is:

T AC l = T AC W × V W + T AC C × V C + T AC B × V B + T AC N D × V N D V l = T AC W × V W V l + T AC C × V C V l + T AC B × V B V l + T AC N D × V N D V l
Note that we only need the fraction of the volume layer that is occupied by each compartment, not the volume itself.

For example, supposing that

V W = 10 100 V l ; V C = 5 100 V l ; V B = 1 100 V l ; V N D = 84 100 V l
We have

T AC l = 10 100 T AC W + 5 100 T AC C + 1 100 T AC B + 84 100 T AC N D
Consequently, in order to estimate the mixing coecients, we need to estimate the percentage of volume of each layer that is a part of each compartment (wall, cavity, blood, or non diusible). In order to obtain the approximation of the stomach wall and cavity TACs ( TAC W and TAC C ), the mixing coecient matrix K has been set up. As previously explained, one assume that the mixing coecients are equal to the volume ratio of each compartment in each layer. A rst rough estimation of these ratio (see Table 6.1), expressed in %, has been done according to the 99m Tc-pertechnetate uptake in stomach: biodistribution modelling anatomical spatial distribution of each compartment in a chosen stomach IHC image (see Figure 5.2(a)). In practice, this rough estimation of the ratio was done by measuring the space occupied by each compartment in a 1250x1250 pixel 2 neighbourhood at dierent position in the stomach wall and cavity, in the 2D IHC slice. For example, with this manual method we estimated that the six rst layers are composed by around 10% of NIS expressing cells and no cavity, the rest being the non diusible space and the extra-cellular matrix. Then, from the seventh to the tenth layers, around 25% of the space is occupied by NIS expressing cells with around 20% of cavity. Next, from the eleventh to the fourteenth layers, we estimated that around 5% of the space was occupied by NIS expressing cells and 80% of the space by the cavity.

Then, in the inner layers, only cavity and non diusible material remain. The proportion is more variable since in the stomach cavity, the non diusible material corresponds to the alimentary bolus and it is hard to quantify its occupation on IHC slices.

Then the mixing coecient matrix consists of only the ratio of the three compartments (Wall, Cavity and Blood) that contribute in term of activity to each of the rst 14 layers. K is of size 14 × 3.

Then the approximated compartment TACs, Ĉ, is computed according to the relation:

Ĉ = ( t K • K) -1 • t K • L
The accuracy of the mixing coecient matrix K is checked by recomputing L = K • Ĉ. The mixing coecient of the rst 14 layers were adjusted in order to minimise the sum of squared dierences (SSD ) between L and L. Table 6.1 presents the nal volume estimated contribution, in %, from which K was extracted. K was then used for all the following steps.

In practice, we chose to only consider the rst 14 layers of the stomach.

The 10 last layers are ignored on all the following of this chapter. The reason for this simplication is that the inner part of the stomach cavity is subject to the diusion phenomenon, that we do not study here. As we decided in Section 6.2.2 to neglect this phenomenon, we decided to ignore the last layers of the stomach cavity. The mixing coecient matrix that has been used is composed by the rst 14 rows (from L1 to L14) and the rst 3 columns (W, C and B).

The mixing coecient matrix has been used to obtain the TACs of the stomach compartments. for most time points (after time point 2) and for all layers. Only the rst 99m Tc-pertechnetate uptake in stomach: biodistribution modelling instant present a relative error above 25% but the activity measured in the layers is really low (smaller than 5 SUV) with an absolute error between 0.5 and 1 SUV.

Overall, the nal SSD summed over all the layer and all the time frames is 633.09, and the root mean square (RMS ) per layer and per time frame is 1.54, which represents around 2% of error for the latest time frames. We can thus consider that the approximated activity found by this method ts well with the observations, and as a result, the mixing coecient matrix is satisfying for this rst approach.

Results

In this section, we describe the dierent steps that conducted to a the rst estimation of the ve transfer parameters (k 01 , k 02 , k 12 , k 21 , and k 10 ) that best describe the 99m Tc-pertechnetate uptake phenomenon.

First, we present the analytic resolution of the inverse problem. The inverse problem approach is a straightforward as soon as we are able to approximate all observations with simple functions such as sum of negative exponential functions. If the input functions and compartment observations are known, we t these with multiexponential functions and we nd analytically the k ij . However, the estimation of TAC W and TAC C with multiexponential functions failed. This is partially due to the fact that there are not enough observations along the time. The acquisition sequence stopped 140 min after activity injection, before any activity equilibrium has been reach. Besides, we have no conrmation of the existence of any equilibrium. Thus, we can not make any assumption about the evolution of the uptake after 140 min, such

an assumption could lead to the convergence of a solution. As a consequence, it was not possible to estimate the transfer parameters of our compartmental model by this analytic way.

That is why, in a second time, we decided to solve numerically the inverse problem. This approach consists in determining numerically the transfer pa- In order to solve the inverse problem, we rst need to approximate the TACs by sums of four exponential functions. The analytic solution we presented in Section 6.2.3 will then allow us to compute the transfer coecients.

To approximate the TACs, we rst tried the method proposed by Fresen and Juritz [START_REF] Jl Fresen | A note on Foss's method of obtaining initial estimates for exponential curve tting by numerical integration[END_REF]. The result is shown on the two sub-gures at the top of 99m Tc-pertechnetate uptake in stomach: biodistribution modelling This is due to the fact that the computed TACs do not contain enough data;

in particular it is impossible to know when the activity will stop growing, or if it decreases, only from the TACs data. If we try to add a point further, e.g. an activity curve closed to zero after a long period of time, the result is slightly better but still unsatisfactory.

To verify that this failure was not only due to the method used, we briey experimented with Maple Fit function in order to try to t the TACs with multiexponential functions. Here again, the software has diculties to nd an approximation, and the approximation found is still unsatisfactory as shown by the two approximation on the bottom of Figure 6.7(c) and (d).

In both cases, not only the approximated curves are far from the original ones, but also the method is not robust and the approximated function varies a lot for a small change in the TACs, e.g. a small change in the mixing coecients.

The problem of tting a function with a multiexponential approximation is known to be a dicult one (see for example [START_REF] Forman S Acton | Numerical methods that usually work[END_REF] p. 253). Here, additionally to the basic complexity of the problem, we also have to face another diculty: there is not enough data to compute a good approximated multiexponential function. Indeed, from the TACs computed by the previous step, it is impossible to predict the evolution of the activity (i.e., whether and when it decreases). Consequently, we did not investigate further the resolution of the problem with analytic resolution.

Numerical resolution 6.4.2.1 Exploration of transfer parameters space

As a rst step, we tested a large number of combination for the set of transfer parameters. Instead of a regular sampling, we try a sampling at a logarithmic scale in order to determine the rough size of the parameters. We choose the following sampling for k ij : {5 × 10 -3 , 10 -2 , 5 × 10 -1 , 1, 5}. The number of possible values in this set are n = 5 and the number of parameters to estimate is z = 5. The number of combination of set of transfer parameter values that were tested is thus n z = 3125.

These sets of parameters were evaluated with respect to the sum of squared dierences (SSD ) in the stomach wall TAC (SSD W ), in the stomach cavity (SSD C ), and in both compartments (sum of the two). In order to nd the best candidate for each transfer parameter, the evaluation function that is minimised is the sum: SSD W +SSD C . Thus, the sets of parameters were sorted with respect to this sum. We additionally computed the RMS error in SUV.

Then, to further evaluate the precision of our approximation for transfer parameters, we consider the matrix C of TACs computed from the approxi-mated transfer parameters for both stomach wall and cavity. This matrix is multiplied by the mixing coecient matrix K in order to obtain the estimated TACs for all layers L. We compare L and L, the matrix for observed TAC in each layer. To make this comparison, we compute SSD L , the sum of squared dierences on the rst 14 layers (between L and L) and its RMS error.

The rst lines of Table 6.2 presents the six combinations of parameters that performed the best with respect to SSD W + SSD C (case (A)), as well as the sets of parameters that performed the best with respect to SSD W (case(B)) and with respect to SSD C (case (C)). This parameter space exploration step has been performed in order to determine initial set of values for the transfer parameters that would permit to t the compartment TACs with sums of four exponential functions.

Case (A): The six sets of parameters are slightly dierent however the RMS error on the TAC, around 60.00, is almost the same for all six sets. This suggests that these sets of parameters could be good candidates as initialisation of an optimisation process.

For the rst set, the RMS error on the compartment TACs is below 60 that represents around 30% of error for the latest time frame for TAC W and Case (B): The set of parameters gives a better estimation for TAC W than the sets of cases (A). However, the RMS error on the wall and cavity TACs is around 145 that represents around 73% of error for the TAC W and around 207% of error for the TAC C . This is not satisfying.

Case (C): The set of parameters gives a better estimation for TAC C than any other set of cases (A) and case (B). The RMS error on the TAC W and TAC C is around 85 that represents around 43% of error for the TAC W and around 120% of error for the TAC C . This is better than case (B) but still not satisfying. Surprisingly, the RMS error on the layer TACs is smaller than for the best case (A).

Indeed, the relatively good results for case (B) and for case (C) with respect to SSD W or SSD C reveal that when the priority is given to the one or the other compartment, the general results on the sum over the two compartment SSD are worse. However, the two compartments are linked together with the parameters k 12 and k 21 . Thus it is not relevant to just compare with respect to one or the other.

One desires to nd a set of parameters that results in a best tting of both wall and cavity TAC. Therefore, one chooses to select the best parameters with respect to the sum SSD W + SSD C . This sum will be used as the 99m Tc-pertechnetate uptake in stomach: biodistribution modelling evaluation function for the tness of the approximated solution.

6.4.2.2 Optimisation of transfer parameters with Nelder-Mead simplex algorithm

In order to optimise the set of transfer parameters, we use the Nelder-Mead (or downhill) simplex algorithm [START_REF] Nelder | A Simplex Method for Function Minimization[END_REF][START_REF] Wh Press | Downhill Simplex Method in Multidimensions, chapter Numerical recipes in C: the art of scientic computing[END_REF]. This downhill simplex algorithm is an iterative process that allows the determination of a local optimum for a multidimensional point which minimises a specied evaluation function.

The space dimension is equal to the number of parameter to estimate z = 5.

A simplex is composed by z + 1 = 6 points, that is one point more than the number of parameters to estimate.

This algorithm is executed twice. A rst execution is done with a simplex composed by six points determined by another mean (here, the parameter space exploration). Then a second execution is performed using as initialisation a simplex composed by six random perturbations of the nal point obtained at the end of the rst execution of the algorithm. As suggested in the previous paragraph, the chosen evaluation function is the sum of SSD of the two compartments TACs.

First execution of simplex algorithm The simplex was initialised according to the result of the parameter space exploration. The six sets of parameters that yield to the smallest SSD W + SSD C have been chosen to initiate the rst execution of the downhill simplex algorithm. As some of the parameters (k 10 and k 01 ) were of the same scale for all sets, we perturbed randomly all parameter values between +/ -50% of their initial value. The main reason for this perturbation is that it allows the algorithm to also investigate the possible variations on k 10 , this would not be the case if the value of k 10 were identical for all the initial points. As this is an iterative process, the stopping condition was the dierence between the result of the evaluation function at the centroid of the simplex and each point of the simplex. When the dierence was lower than ε, the process stopped. We choose ε = 1000 for this rst execution since it is a hundredth of the scale of the result for SSD W + SSD C in case (A). The set of parameters after stopping condition that led to the smallest result of the evaluation function is presented in Table 6.2 (case (D)). The resulting estimated TACs are also compared visually to the approximated ones in Figure 6.8(b). Qualitatively, the estimated TAC W (blue dashed line) is somehow closer to the approximated TAC W than it is the case for the TAC W (green continuous line). This is probably due to the fact that TAC W is around 2.5 greater than TAC W and the resulting error has more weight in the evaluation function.

Quantitatively, the process stopped after 32 iterations and the set of parameters results in an RMS error on the compartment TACs around 32, that Quantitative results of the parameter space exploration (see Section 6.4.2.1) and after downhill optimisation with the simplex algorithm (see Section 6.4.2.2) for subject #0530. This table presents several sets of transfer parameters and the errors measured between the observations ( Ĉ and L), and the estimations ( C and L).

The errors are presented in columns (1) to (6) The results of the parameter space exploration are presented in cases (A), (B) and (C). The results of the optimisation with downhill simplex algorithm are presented in cases (D) and (E). (A1-A6): The results for the best six sets of transfer parameters with respect to column (3) are presented since these six sets will be used as initialisation of the downhill simplex algorithm that requires an initialisation with a number of sets equal to the number of parameters + 1. (a) Best result among the 3125 combinations tested for the parameter space exploration, case (A1) in Table 6.2 (b) Best result after the rst execution of the downhill simplex algorithm, ε=1000, case (D) in Table 6.2 (c) Best result after the second execution of the downhill simplex algorithm, ε=10, case (E) in Table 6.2 is equivalent to the half of RMS for the best case of (A). This is equivalent to around 16% error on the TAC W for the latest time and around 45% error on the TAC C .

Second execution of simplex algorithm The rst execution of the simplex algorithm did not result in a satisfying estimation of the transfer parameters according to the resulting estimated TAC. Thus, a second execution of the simplex algorithm has been performed. This time, the stopping condition has been divided by 100 to reach ε = 10, and the simplex has been reinitialised with a set of six perturbations of the resulting set of parameters of case (D). The perturbation was randomly chosen between -50/ + 200% of the parameter values.

The set of parameters after stopping condition that led to the smallest result of the evaluation function is given in Table 6.2 (case (E)). The resulting estimated TACs are also visualised and compared with the approximated ones in Figure 6.8(c).

Qualitatively, the estimated TAC W (blue dashed line) ts well the approximated TAC W and so does the estimated TAC C (green dashed line) with TAC C . Thus, from a qualitatively point of view, the estimation of the transfer parameters is satisfying.

Quantitatively, the process stopped after 83 iterations and the set of parameters results in a RMS error on the compartment TACs around 3.3, that is equivalent to the tenth of the RMS error in case (D). This is also equivalent to around 1.6% error on the TAC W for the latest time and around 4.7% error on the TAC C .

Finally, the RMS error on the layer TACs, 12.32, is equivalent to around 30% less than the RMS error in case (A). In order to visualise the distribution of the error, we present respectively on Figure 6.9(a) and Figure 6.9(b) the signed error in SUV between the observed layer TACs and the estimated layer TACs L -L, and the relative error |L-L|

L

, at each time frame and for the rst 14 layers. One can see that for most of the time frames and layers, the error does not exceed 4 SUV in absolute. The error is larger for the three last time frames (17 to 19) and for the seven inner layers (layer 8 to 14). Indeed, this is conrmed by the relative error that is greater to 10% for these instants and layers. On the other hand, expected for the rst time point, the relative error does not exceed 10%-15% for most of layers and time points. At this step, one can conclude that the proposed approach is excellent. biology. As k 01 < 0, one can assume that k 01 = 0, this will mean that no 99m Tc-pertechnetate is secreted by the NIS expressing cells outside the con-sidered system. This is a reasonable hypothesis according to the biologists.

Biological interpretation of the transfer parameter values

However, the three other parameter values pose a problem in terms of biology. According to the biologists knowledge, k 12 should be small comparing to the others and it is not the case. The overestimation of k 12 might be a consequence of the short circuit of both k 02 and k 10 . This has to be conrm by additional experiments and a rened model of the system. Consequently, we cannot claim yet that the results are biologically relevant.

Discussion and Future works

In this chapter, we established a methodological approach based on compartment analysis for the modelling of the 99m Tc-pertechnetate uptake in the murine stomach. Our objective was also to study the feasibility of this approach. In order to simplify the model and reduce the number of parameters of the model, several assumptions have been made. The results (transfer parameters) that have been obtained at the end are satisfying knowing all the assumptions.

From this feasibility study, one can conclude that our modelling, together with the methodology proposed in this chapter, allows the estimation of parameters, here the transfer parameters between compartments.

The originality of this feasibility study is to take into account information that was available on some functional units such as the presence of NIS expressing cells in the glandular part of the stomach wall.

In this section we review each choice and assumption done in this work and discuss them. We conclude this section by proposing future works; these research directions present the limitations of our approach and how we could improve it, the next steps that should be undertaken to exploit the methodology we developed, and the alternative directions and choices that should be investigated.

General approach

Our objective in this work is to study the iodide uptake in the murine. For this study, we rst provided, in Chapter 5 qualitative results on this evolution.

The data available, and the experimental conditions conducted us to a rst choice concerning our approach.

Region-based approach (Choice 1) Iodide is mainly accumulated in two regions of the stomach, the wall and the cavity. Thus, one would like to consider only two regions for the uptake study. However, these two spatial regions are not homogeneous. Then, the regions should be spatialised in order to take into account the diusion phenomenon in the cavity. A voxel-based approach is not possible to study iodide uptake 99m Tc-pertechnetate uptake in stomach: biodistribution modelling in the murine stomach because of the movement and deformation of the stomach (see Chapter 5). Thus we considered an approach that permits us to gather the voxels in homogeneous regions. Once the regions dened for each image, we associate each layer along the time. We made the assumption that the layers also refer to the same anatomical region. This kind of approach was already advocated by Zang et al. for functional MRI [ZJL + 04], where voxels are gathered into homogeneous regions.

Regional Homogeneity (ReHo) is a voxel-based measure of brain activity which evaluates the similarity or synchronization between the time series of a given voxel and its nearest neighbors [START_REF] Zang | Regional homogeneity approach to fMRI data analysis[END_REF].

The diculty here, that diers from parcel-based approaches with parceling/regional homogeneity in functional MRI [ZJL + 04], is that our regions are interlocked layers: we base our analysis on a measure of the activity inside each layer of the stomach, each layer being composed of a mix of the dierent stomach regions.

Model proposal: compartmental modelling

To provide a model for the phenomenon observed and presented in Chapter 5,

we proposed a rst model of the stomach based on compartmental analysis, and then two steps of simplication. Those simplications were needed in order to deal with the equations of the model both analytically and numerically.

The general assumptions we made are related to compartmental analysis.

The system follows compartment modelling prerequisites (Assumption 1)

This assumption can be split into two hypotheses.

Temporal invariance of transfer functions The transfer functions between any two compartments are constants.

Nature of transfer equations We supposed that the transfer equations between compartments are ordinary dierential equations (ODEs), which are linear and of rst order.

The fact that our analysis concludes with relatively stable results shows that these two hypotheses were reasonable, at least in a rst time.

Observability (Assumption 2) We can observe the concentration of 99m Tc-pertechnetate in each compartment, and use these observation as system outputs. This hypothesis is somehow justied by the observations of Chapter 5 and Choice 1, but building the activity in each compartment from the layered observation is not trivial, as explained below.

In this work, we decided to study our subject with compartmental analysis but other approaches could have been considered. For example, another kind of approach used in medical imaging is the factor analysis of dynamic structures [START_REF] Frouin | Image Sequence Processing Using Factor Analysis And Compartmental Modelling[END_REF]. Its aim is to estimate the kinetic of dierent physiological mechanisms.

Computing TACs

To build the system input and outputs, we need to construct time-activity curves (TACs) for the input function, i.e. the activity in the blood, and each of the observed output, i.e., the activity in the stomach wall and cavity. The activity in the blood could be observed from the activity variation in the left ventricle of the dynamic SPECT images. For obtaining the TAC in the two other compartments, we make the following assumptions.

Time dependency of composition of each layer (Assumption 3) We made the assumption that the composition of each layer, in terms of the volume occupied by each compartment, does not change with respect to the time. This assumption does not take into account that, although the animal is anaesthetized, its stomach might be subject to local deformation, such as contraction or relaxation. More precisely, we made the assumption that the n th layer had the same composition over the time, but if the stomach contracts the layer number n will get closer to the stomach cavity and have a dierent composition. To summarise, the composition of each layer might slightly vary over the time but we neglected this variation in a rst time. In the future more complex segmentation techniques might be considered in order to limit the layer mismatch.

From these assumption and from an estimation of the mixing coecients that express the activity in each layer depending on the activity in each compartment, we approximated the TAC of each compartment. Then, we veried that the TACs of all the layers reconstructed from the TAC in the compartment was precise enough.

At the moment, the mixing coecients have been estimated from a simple observation of the IHC images, but a more precise estimation of these coecient could be performed (see Section 6.5.6 below).

Computation of the transfer parameters

Once the TACs and the compartmental model are known, the objective is to nd the transfer parameters in order to fully characterise the system. This characterisation can be done by solving the inverse problem. The inverse 99m Tc-pertechnetate uptake in stomach: biodistribution modelling problem consists in estimating the transfer functions from longitudinal observations. This inverse problem can be solved in two dierent ways: a direct approach that uses a semi-analytic resolution, or an indirect approach that estimate the set of transfer parameters with optimisation approaches. Here, we tested both ways but the chosen resolution is the numerical one with optimisation steps.

Numerical resolution (Choice 2) The attempt of determination of the compartment TACs with a sum of four exponential functions failed with a mixed analytic/numerical approach. The reasons for this failure are the diculty to t a function with a sum of four exponential functions, and the lack of data. Although, this kind of multiexponential function tting consideration can be improved with additional optimisation method. However we did not investigated further this direction, this might be a perspective of our work. Consequently, instead of trying to t sum of exponential functions, we build TACs from dierent transfer coecients, and t them with the TACs obtained by the preceding step.

This last step allowed us to nd transfer parameters by combining exploration and optimisation steps. The optimisation was performed with a simplex algorithm. We nally obtain a set of transfer parameters that t relatively well with the observed phenomenon.

The numerical resolution has been possible thanks to several simplication of the model toward a classical model (two-compartment system). Although, the solution that has been found is relevant for the system. In order to better match with the biological functional units, one should go back to a more complex model. Then, another kind of approach might be considered for the numerical resolution since there will be more transfer parameters to determine.

Indeed, the dimension of space of values to explore will increase by 1 for each new transfer parameter to nd. Consequently, the exploration of the space of values will be more complex and the time of computation will increase.

Evaluation of the method

We did not conduct complete evaluation of the robustness and adequacy of our method, however, our rst results are satisfying because the method allows us to compute transfer parameter, and the computation seems relatively stable.

Indeed, rst experiments where we tried to perturbate the solution or to nd other optimums for the transfer coecients led us to similar coecients (same order of magnitude). We can thus conclude for the moment that this analysis proved the feasibility of the study, and that our methodology is appropriate.

However, further developments and probably renement of the methodology, will be necessary to improve the robustness of the results and the precision of the stomach model we propose. It is also necessary to validate and exploit our results in terms of biological interpretation.

Future works

The rest of this section will detail the further developments that are, to my mind, necessary to overcome the main limitations of our approach and to improve our method and results.

Mixing coecient describing the composition of each layer

In a rst step, we estimated from observation of the IHC images the mixing coecients. We consider this rst estimation as satisfactory as it leads to reasonable results and shows the feasibility of our study. However, we wonder if this estimation from a single observation of a single IHC image is relevant. Investigations are conducted on the reconstruction of a three (spatial) dimensions (3D) image of the murine stomach from a set IHC images 2 . This 3D image will be considered as an atlas for the stomach and will be registered to the SPECT images. It will provide information on the spatial distribution of the dierent functional units that have been considered. Considering the number of approximations that will be performed such as the resampling required for the IHC, the nal atlas should not be used as a ground truth. The resulting 3D stomach atlas could be used to assess the mixing coecients that have been determined in Section 6.3.2.2.

Robustness of the estimation of transfer parameters k ij

While we checked that small variations of the starting values for the transfer parameters led to similar results, we did not perform any systematic robustness test on our obtained results. Here are some steps that we plan to do in the future in order to demonstrate the robustness and reproducibility of the estimation of the transfer parameters.

Exploration of the parameter space. We started from 5 values for each of the 5 parameters. While extending much this step costs a lot of computation time due to the combinatorial explosion of the explored space, as the current exploration is performed in a few minutes, we could slightly extend the range of parameters investigated.

Direct comparison with L (observed layer TACs) instead of Ĉ (approximated compartment TACs). At the moment, we approximate the TACs of the compartments from the ones of the layers, 2 Thanks a lot to Audrey Lamit, Julien Guglielmi and Philippe Pognonec for their time, help, support and involvement in the planning, experimental and acquisition process. Chapter 6. Dynamic SPECT analysis of 99m Tc-pertechnetate uptake in stomach: biodistribution modelling and optimise the dierence between the TACs obtained from the layers and the ones obtained from the compartmental analysis (i.e. from the transfer parameters). Instead, we could apply the mixing coecients on the TACs obtained from the compartmental analysis and obtain approximated TACs for all the layers; then we would directly optimise the dierence between the observed and the estimated TACs for the layers.

In other words, this approach consists in optimising directly the last two columns of Table 6.2. A rst experiment seems to provide similar results with the two methods, but this should be further investigated.

Modication of the initialisation of the two execution of downhill simplex algorithm. Several alternative ways of choosing the initial points of the downhill simplex and of perturbing this initial points should also be investigated.

Among these investigations, one promising direction is to run the second execution of the downhill simplex on a large number of perturbated initial set of parameters in order to have statistics (average/standard deviation...) and then check the reproducibility of our results.

Increase the precision of the downhill simplex algorithm. This consists in reducing the stopping condition ε, and increasing the maximal number of iterations. This requires more computation time and will not change much the results, but it should give us results with an improved precision.

Improve the speed the convergence To improve the speed of convergence of the simplex algorithm, we could try to minimise the dierence of normalised compartment or layer TACs, instead of the simple dierence.

Constrain transfer parameter search to positive values For the moment, we also explore the negative values for the transfer parameters while the biological signicance of such a negative value is not clear. We could reduce the exploration space to positive values for transfer parameters, but this could perturbate the eciency of the simplex algorithm.

Though it could prevent us from detecting unreasonable results, especially if the negative value is not close to zero. In our experiments, we have reached a negative value for one transfer parameter, in the last iteration, but this value is quite small. Consequently, constraining transfer parameter search to positive values might not provide a signicant gain.

Extending the dataset

At the moment, this work has been conducted on only one animal and one sequence of SPECT images. In order to fully validate the approach, it would be interesting to perform a group study. The rst outcome of such additional tests would be to know if the transfer parameters are really dierent, or if they are similar. Having parameters relatively similar for dierent animals would also convince us of the robustness of the results. However, to use another set of data we need it to satisfy some constraints.

These constraints are linked to the conditions of the imaging procedure. First, the image eld of view should include the blood in order to be able to estimate the blood input function. Then the type of injection should be the same.

Finally, the number of frames and the time lapse between each frame should be sucient. In other words, the total sequence of acquisition should contain not less than 19 time frame acquired over at least 2 hours. Among the set of data available (see Table 5.1), only four sets match this prerequisites: the subject #0530 studied in this chapter and the subjects #0366, #0370, and #0371. However, on the remaining 4 sets, only the subject #0530 has the heart visible in the eld of view.

Biological interpretation and validity

Last but not least, a rst interpretation of the rst result has been provided by the biologists. However, one result is not sucient to conclude on the relevance of the approach and on the biological interpretation. Indeed, some transfer parameter values that have been determined do not reect the biological reality. Consequently, it will be more relevant if we check rst the robustness and reproducibility of the method, and if we increase the complexity of the compartmental model. In this thesis we address the study of longitudinal phenomena with application to small animal SPECT imaging. Among the wide range of processes that can be studied, we focus on two particular subjects with biological challenges that led to the development of dedicated methodologies. Overall, several challenges related to the work with animal imaging were addressed.

First, the elaboration of a relevant biological acquisition protocol is a demanding task which requires a lot of interaction with all the contributors (biologists, computer scientists), and a good understanding of all the variables that inuence the biological and imaging outcome. The protocols had to deal with biological variability since the studied biological phenomena were not well known. For example, a new in vivo tumour model was considered for which the in vivo injection protocol was still on trial (cell lineage, animal lineage, sensitivity to radioactivity dose or Xray dose). The site of injection and the site of tumour xation were also very sensitive to change, hard to predict and yield to a dataset with a lot of variability. The setting up of longitudinal acquisition also required organisation in order to estimate the speed and duration of evolution of the observed phenomenon. Indeed, the study of longitudinal phenomena suggests quantication and/or modelling of the phenomenon over time, and the time scale of a new phenomenon might be dicult to predict.

Second, both phenomena that are studied here are observed in the abdominal region, which lead us to two additional sub-challenges. The rst one is the lack of imaging that would serve as detailed anatomical reference, indeed a simple computed tomography (CT) scan does not oer a good contrast in the soft tissues. And the use of additional contrast agent without interference with the protocol is not guaranteed. One could consider approaches such as the ones presented in [START_REF] Maroy | Quantitative organ time activity curve extraction from rodent PET images without anatomical prior[END_REF]. This method succeeds in extracting TAC from rodent PET images without the anatomical information.

The second one but not the least is the intrinsic motion of the considered area. Indeed, the subjects were anaesthetized mice, thus there were still breathing and their stomach and surrounding organs were moving. In both studies we thus had to take into account the movement of the organs.

Contributions

The major contributions were in terms of methodology. The rst method addresses the problem of reconstruction of SPECT images without respiratory motion artefacts. The second method proposes a model for the kinetic of 99m Tc-pertechnetate uptake in the murine stomach.

Simulated breath-hold reconstruction in micro-SPECT In order to avoid respiratory motion artefacts in SPECT images, we proposed a method that selects retrospectively and maximises the data used for a motionless reconstruction. This could be seen as a breath-hold like acquisition or a smart respiratory gated acquisition. This method guarantees a high signal-to-noise ratio (SNR) while eliminating respiratory motion artifacts. Moreover this method is adapted to breathing rhythms particular to anaesthetized animal (gasps) with high variability in terms of amplitude and length. This method is based on the extraction of intrinsic information; it is robust to the apparent asynchronism between the true motion in the images and the respiratory signal externally acquired. The quality and accuracy of reconstructed images are comparable to a non-gated reconstruction in terms of noise measure and SNR, and intermediary between gated and non-gated reconstruction for lesion measurements (SUV peak , SUV mean and lesion volume). Thanks to this method, new and more challenging studies concerning organs or tumour modelling in abdominal region can be addressed with SPECT imaging.

This contribution led to the development of a rst version of a software implementing this method. This software with a graphical interface allows biologists to obtain a simulated breath-hold reconstruction from a SPECT acquisition.

Study of 99m Tc-pertechnetate uptake in the stomach In order to help biologists to understand the process of iodide uptake and secretion in the stomach, we propose a method that permits us to build a rst compartmental model describing the iodide analog 99m Tc-pertechnetate uptake in the murine stomach. This conducted to the design of a compartmental model together with a rst estimation of hidden parameters, i.e., transfer parameters between compartments. The method relies on an original layer-based description of the stomach that permits us to obtain a relevant anatomical information of the stomach, and that makes up for the progressive 99m Tc-pertechnetate uptake and irregularities in the stomach shape. Additionally, this layer description is robust to global stomach deformation such as contraction or relaxation. Still, several challenges had to be taken into account. A rst challenge concerns the intensity-based stomach segmentation process without anatomical prior information. This segmentation step failed in the early time frame (just after activity injection) since the site of injection that was in the eld of view, was segmented instead of the stomach. A second challenge was that each layer does not correspond to only one compartment. It is a mix of several compartment and in a certain sense, even if each compartment is supposed to be homogeneous, a layer is not. We thus considered the set of layers as a mix of compartment, each mixing coecient being representative of the composition of each layer. A rst determination of these mixing coecients was proposed and yields really good results relatively to the method that was used to determine them. A few other challenges were only partially addressed, like the robustness of the method, or the non-uniform composition of each layer. Also, the method makes several simplifying hypotheses. For example, we do not take into account the local temporal stomach deformation that induces a slight layer mismatch between the outer layers. The rst layer at a given instant might disappear in the next instant. Despite these points that can be improved, the method yields relevant results in terms of estimation of time-activity curves.

Future works and perspectives

The development of these two methods open a wide range of possibilities.

First we present several direct improvements of our methods and results that could be envisioned. Then we propose a few longer term research directions. Validation of the breath-hold like reconstruction method on simulated data This validation step consists in simulating the SPECT imaging of breathing mice in order to obtain a large set of data. These numerical mouse phantoms should include small lesions randomly located in the abdomen and a variable breathing rhythm. More precisely, one would have to animate a numerical mouse phantom with the same particular respiratory pattern that has been observed (gasps followed by a rest period). However, such a model has not been developed yet. For example, the MOBY phantom [STF + 04] implemented a respiratory motion similar to humans, and would not be adequate.

Additionally, such a work requires the elaboration of a numerical model for our camera in order to simulate the acquisition and reconstruction processes.

The disadvantage of this validation process is that it involves a set of data that will be perfectly produced by the simulation software without reecting the reality of the manipulation and acquisition of living animals. However, such a simulation will provide a large set of data that will be valuable for the validation of our method.

Iodide uptake study in stomach Most of the future works related to this subject have been detailed in Section 6.5.6. Among them, the main points are the study of the robustness of the approach; the assessment of mixing coecients from IHC slices; and the biological interpretation of the phenomenon.

Long term research directions

I present below a couple of longer term research objectives, they consist in improving the methods in the same direction of each approach I presented, but also in combining both approaches in order to study new research directions.

Abdominal tumour growth quantication and modelling Studying the abdominal tumours, and quantifying and modelling their growth was the initial goal of this thesis. However, this objective required preliminary works on the respiratory motion correction. This preliminary work turned to be more challenging than expected and has been addressed in Chapters 3,4 and in Appendix A of this thesis. From this, the tracking and accurate quantication of tumours in the abdominal region will rely on strong basis. A tumour matching algorithm [WGR + 01] should be developed in order to track the tumours along the time by taking into account their growth, the possible fusion of spatially close tumours, or the apparition of a new lesions. This kind of approach should include registration methods such as block matching [START_REF] Ourselin | Block Matching: A General Framework to Improve Robustness of Rigid Registration of Medical Images[END_REF]. Then, the tumour growth in the peritoneal cavity can be studied and will set up a computational model for tumour growth [START_REF] Konukoglu | Tumor growth modeling in oncological image analysis[END_REF] in a deformable environment. Such a model should greatly help the biologists in the understanding of the growth of this kind of tumours.

Layer-based approach applied to the tumour heterogeneousness

The original layer-based description of the stomach presented in Chapters 5 and 6 could be extended to other organs but also to tumours. Indeed, some tumours can be observed until a late stage. At late stage, the segmentation of tumours based on SPECT images reveals only the active part of the tumour, that is constituted by vascular cells. However, it does not reveal the necrosis part of the tumour. In the case of subcutaneous tumours, this can also be observed on CT imaging since the necrosis part sinks but not for other tumours.

First, in the case of NIS expressing tumorous cells, the necrosis of the tumour can be predicted by a decrease of the uptake in the internal layers of the tumour in a preceding imaging.

More generally, it has been also mentioned in [UTT + 12] that the tumours interiors is never homogeneous, for the reasons mentioned above. In this article, authors presents the results of the rst investigations done with high resolution SPECT imaging of small animal: they assess the heterogeneousness of tumours.

More interestingly, the idea of applying a layered-based approach to tumours is to better understand the structure and evolution of the tumour.

Indeed, we could study the tumour not as an homogeneous system, where we calculate the average activity on the whole tumour. Instead, tumours could be studied as a set of subregions such as layers. We would calculate the average activity in each subregions. This would allow the identication of dierent parts of the tumour: necrosis, proliferating cells, but also some other intermediate parts. We hope that such an approach, with a region-based analysis of the tumour, would provide a much better understanding and modelling of the tumour growth.

A.2 Preliminary investigations on the respiratory signals We conducted extensive experiments but we could not nd a law expressing what the length or amplitude of breathing depends on. This could be expected from the initially observed high variability of the respiratory cycle, but it was crucial to try to nd such a law because its existence would have helped the setting of the threshold during the acquisition, and could also have been exploited in our motion correction method. We tested:

the inter individual variability on dierent days and on the same day of acquisition;

the intra individual variability on follow-up acquisitions.

All these observations were done on a set of data acquired with the same respiratory monitoring constraints and recorded with the respiratory monitoring system. For all acquisitions, the operator adjusted the anaesthetic gas rate to keep the animal breathing at the same rhythm. An identical pressure threshold on rising edge was chosen for all acquisitions, this threshold was used to delimit each cycle. An average respiratory cycle was computed from all cycles for each acquisition. At this preliminary stage, the average respiratory cycle is computed by summing up all the extracted cycles as they are (without resampling) and by dividing it by the number of cycles. As all the cycles were not of the same length we added to each cycle parts of the previous and next cycles in order to have all cycles length at 1400 ms: 400 ms before the beginning of each cycles and up to 400 ms after each cycles, so that all the intervals were of the same 1400 ms-length and could be summed up. Later, i.e. in all experiments shown in the following sections, each cycle will be resampled to the average cycle length before the summing since the variability of cycle length is high. To conclude, the hypothesis of a rule dependent to the subject does not hold either.

Third, for the same two rst lineages (C57Bl/6 and NOD-SCID), the same animals were not only imaged several times but the acquisitions were done over two or three weeks. We compared then the average cycles with respect to the In this section, we investigate the similarity between the respiratory signal and the motion information intrinsic to the image. On the one hand, the average respiratory cycle is computed from the respiratory signal, as explained First, the graphics suggest that the lesions are moving periodically, with the same frequency as the respiratory motion.

Second, the graphics show that the lesion motions are temporally shifted with respect to the moving phase of the average respiratory cycle. This suggests that there is a time shift between the lesion motion recorded in the listmode acquisition and the respiratory signal recorded with the pressure sensor. However, in the introduction chapters, we made the assumption that the respiratory motion could be represented by the pressure signal that is recorded with the respiratory monitoring system. With this observation, it seems that this assumption should be considered with caution. The respiratory signal record, combined with a rising edge threshold detection, allows for a good cycle splitting since we observed that the movement of the lesion was following the same law as the average respiratory cycle concerning the length of the cycles. However, the respiratory signal can not be used to directly determine with accuracy the phase of the cycle in the listmode since we observed a slight time shift between the lesion motion and the signal variation that is smaller than the time shift observed in Figure A.8 between the respiratory signal and the intrinsic measure from 4D CT. We investigated the possibility that this time shift could be dependent on the imaging modality used; this is the reason why we also studied 4D CT images.

CT intrinsic signal In 4D CT images, the motion signal can be characterised either by the volume variation of the lungs or by the motion of some particular points of the lungs. The part of the lungs that is the most aected by the respiration in mice are the bases, just above the diaphragmatic dome.

Here, both lung volume and diaphragmatic dome motion have been calculated. More precisely, two distinctive points of the lungs were localised, at the top and at the bottom of the diaphragmatic dome.

A.3.1 Respiratory cycles duration variations

First, we studied the variability of the length of the respiratory cycles within each of the eight acquisitions considered. As the breathing rhythm was definitely not stable during the whole acquisition, the rst question concerned the distribution of the cycle length, and if this distribution was Gaussian or not.

We would like to suppress the data recorded during the outlier cycles (cycles that are much longer or shorter than the average) from the image reconstruction process. The reason for discarding these outlier cycles is that these much longer or shorter cycles suggest a default in the acquisition (missing data, bad detection of rising edge, ...). In addition, we would like to be able to estimate how much data will be discarded at the cycle selection step. If the distribution is Gaussian, we could x the threshold at the same value (typically a multiple of the standard deviation) for all acquisitions. This threshold would have to be chosen in order to select only the cycles within a reasonable range and remove only the outlier cycles. Moreover, the number of discarded cycles could be estimated depending on the threshold that has been used. reconstruction in micro-SPECT: Supplementary data the skewness: that is the symmetry of the distribution with respect to the average length;

the kurtosis: that is the attening of the distribution.

The visual examination of the gures already suggested that the cycle length distributions are not Gaussian although some of them look relatively close. The additional test described above conrmed that none of the acquisition could be safely approximated by a Gaussian distribution.

As our objective here is to discard the outlier cycles, i.e. the cycles that are longer or shorter than a given length range, we decided to perform a selection cycle step on the respiratory signal according to a cycle length. One needed to dene the cycle length window that preserve the majority of the cycles.

As we observed, according to the width of the length window, the ratio of remaining cycles after the selection outliers cannot be described by a general law. However, for all acquisitions, the length window [µ -3 * σ; µ + 3 * σ]

includes most relevant cycles. Thus, this window has been used to initialise the cycle selection step of the reconstruction process.

Another observation can be done from A second improvement that is suggested by this variability is that all cycles should be resampled to the average cycle length. This cycle resampling step will increase the chance for a temporal window to correspond to the same phase in each cycle. We assume that the length of each respiratory cycle phase is proportional to the cycle length. Resampling will then be an indirect manner to work with relative temporal window. However, the reconstruction system does not allow for the denition of temporal window proportional to the cycle length. The system only accepts absolute temporal window. As a consequence, a certain number of events are misused (see Table 4.2). These misused events justify a resampling step. In other words, cycle length resampling should be applied in order to assign each event to the right phase. This will guarantee that an equivalent number of event is used for the reconstruction of each frame of the 4D image. Indeed, if resampling is applied, then the instant corresponding to the beginning of the cycle plus N milliseconds will belong to the same phase in every cycle. On the contrary, in a non-resampled signal, the beginning of the cycle plus N milliseconds could be for example in the rst half of a longer cycle and in the second half of a shorter cycle.

After this resampling step, we will notice a slight variability in the amplitude of the cycles. Indeed the variability in amplitude shown in As we saw in the previous section, the cycles were not all of the same length before resampling. The sequence of cycle lengths in one acquisition can be seen as a pattern. This pattern is assumed to be identical in the two triggering signals. This pattern is used here to nd the shift in cycle number (cycle shift) that registers the two triggering signals.

A rst approach is based on the assumption that each cycle has the same length in both triggering signals. We naively expected that if we try all possible cycle shifts, there will be one for which each cycle ts the length of a cycle in the other. However none of the shifts succeeded in nding a null dier- A second approach is based on the assumption that there is a cycle shift We thus assumed that for some cycles the dierence of length is so large that the least square minimisation failed. A new assumption is done which considered that most cycles have a 1 ms-dierence in length but it might happen that few outlier cycles have a larger dierence in length for unknown reasons. A third approach was tested based on the minimisation of the number of cycles that have a length dierence greater than 1 ms in absolute value (see the third column of Table A.1). Cycle shifts identical to the second approach have been found for all acquisitions except for the acquisition 1839.

A new shift was determined. This dierence is decreasing over the cycles for both acquisitions. It suggests that a lot of cycles from the respiratory signal are missing some milliseconds.

In total, this represents around 600 missing milliseconds for the acquisition 0182 and around 3500 missing milliseconds for the acquisition 1839. The two graphics show a cumulative dierence between the cycle length (respiratory signal minus listmode) that reaches up to around two cycle lengths over the whole acquisition. The slope of the left hand gure suggests that the cumulative dierence is more or less linear. However, the graphic for acquisition 1839

shows a high variation of this cumulative dierence at the end of acquisition, after around 1300 cycles. Something strange happened during the record of either the respiratory motion or the listmode. This is why we also applied the method on a truncated signal of both acquisition 1839 and 1841, which leaded us to better results. As the cumulative dierence seemed to be linear, we expected that the loss of milliseconds was regular. Instead of remaining identical, the number of missing milliseconds randomly varies around 3 or 4 per 10 cycles, with some exceptions. For the acquisition 1839, we also observed the number of missing milliseconds but only on the rst 1200 cycles (before the abnormal variations). This number also randomly varies around 3 or 4 ms per 10 cycles and up to 31 ms of dierence for one sample of 10 cycles. These observations can be completed by additional measures presented in the forth column of Table A.1. On the right column, the numbers of cycles that have a dierence of -1 ms, 0 ms, 1 ms, and an absolute dierence greater than 1 ms are given. For all the acquisitions, around 4 ms are lost every 10 cycles; this ratio seems to be relatively uniform. However reconstruction in micro-SPECT: Supplementary data between them is that the gating signal has been triggered either on the falling edge or on the rising edge, with or without an additional temporal delay.

This additional delay is shifting the gating index in order to synchronise data in the middle of the motionless phase. A motion detection histogram was superimposed to the average respiratory cycle. Figures show that in any case, the phases that received the more votes were the same, with respect to the peak of respiratory motion. We should mentioned the fact that for the case of rising edge detection with a delay (see Figure A.17(d)), the average respiratory cycle is slightly more noisy than for the other methods. This is due a higher variability of the breathing rhythm during this acquisition featured a higher variability. This result strengthen the wide application to this method. In practice, the rising edge detection of the motion, without delay, was chosen by default for the acquisition protocol. However, this choice does not aect the result of the motion detection method. Any other parameters could have been chosen.

A.5 Conclusion

In this appendix, we presented additional results, which both explain the reasons that leaded us to the method presented in Chapter 4, and show that this method is eective.

The rst analyses we conducted leaded us to the conclusion that it is very dicult to obtain acquisitions with stable breathing either in frequency or in amplitude. Consequently our method must be able to deal with such highly variable breathing rhythms and amplitude.

Second, we observed that the movement recorded by the pressure balloon was indeed representing the respiratory motion of the mouse. However, we also observed a temporal shift between the balloon observations and the lesion movement in the images; this temporal shift is dependent on the position of the lesion and we were unable to nd a systematic law for predicting the temporal shift. Overall, these observations make the detection of motionless phases more dicult because they imply that it is impossible to predict the characteristics of the breathing movement, neither for the movement observed by the pressure sensor, nor for the movement of the observed lesions when the pressure sensor movement is known. Consequently, it is necessary to design a new method for detecting motionless phases only from the information contained in the listmode, and tolerating a high variability in the respiratory rhythm. We also realised that there was no perfect match between the respiratory signal and the gating information stored in the listmode, which leaded us to rely only on the information intrinsically available in the listmode.

Finally, this appendix also provides additional arguments to show the efciency of our method, and to justify the dierent steps we designed for the reconstruction of a motionless 3D SPECT image. We showed that the lost events for the G4D reconstruction were coming from outlier cycles that were impairing the reconstruction because they were introducing a wrong count of events; some events being counted several times, and some other omitted. We also have shown in more details that our reconstruction methodology provides a much more reliable compromise between SNR and accuracy of the results that the existing classical methods (G4D and NG3D). We concluded this appendix by a short study on the sensitivity to noise; this study has shown that our algorithm is relatively reliable and that it is appropriate to use it for the kind of images provided by our camera.

ies only one substance, each compartment refers to a dierent space; and the interactions between compartments are ow of substance. When two substances or more are studied, two compartments can overlap and occupy the same physical space. In that case, the interactions between compartments are transformations from a substance to another.

The interactions between compartments are described by ordinary dierential equations (ODEs). The resolution of the ODEs will give the functions that describe the quantity or concentration of each substance in each compartment.

To summarise, the compartmental modelling implies the denition of the set of compartments based on physiological and anatomical information, and the determination and the resolution of ODEs.

B.2 Principles

The denition of a model requires a certain number of steps. First, we need to dene the compartmental system by its elements and under some assumptions.

B.2.1 System denition

State variable, x i , with i ∈ 1..n, and n the number of compartments, System inputs, u i = external input on compartment i, System outputs, y = observation of compartment(s), Unknown parameters, parameters characterising the transfer functions between the dierent compartments, but also with the system inputs and outputs, they are represented by three matrices A, B, C.

More precisely, each element k ij of the matrix A represents the speed at which the substance transits from compartment j to compartment i. Similarly, each element b ij of the matrix B represents the speed at which the substance transits from input j to compartment i. Overall, the derivative 1 of the concentration of the substance x in compartment i depends on the substance x i in compartment i that ows to other compartments j with a factor k ji , on the substance x j in each compartment j that ows to compartment i with a factor k ij , and of the input u q that injects substance into the compartment i with a factor b iq , more formally: ẋi = k ji x i + k ij x j + b iq u q 1 We denote by ḟ the derivative of f with respect to time.

The matrix C represents the relation between the system output functions y (e.g. observations) of each compartment i of the system and the substance concentration x in the dierent compartment. This represents the fact that the activity measures that are done on each compartment represent a fraction of the activity or a mix of the several compartment activity that is eectively there. Here, we consider that C is the identity matrix, more formally:

y i = c ij x j
All these elements can be summarised and represented analytically.

x =     x 1 x 2 . . . x n     , u =     u 1 u 2 . . . u n     , y =     y 1 y 2 . . . y n     , ẋ = Ax + Bu, y = Cx B.2.

Assumptions

The assumptions are essential in compartmental analysis since they will yield to simplication of both the system and the resolution of the ODE. Those assumptions are: the transfer functions between any two compartments are constants, denoted = k ji (or k j←i ) from compartment i to compartment j. This assumption is very frequent in a rst time modelling; it is called the temporal invariance.

the transfer equations between compartments are ODEs, which are linear and of rst order.

The system outputs are the state variables, in other words, we can observe the concentration of 99m Tc-pertechnetate in each compartment, and use these observation as system outputs.

B.2.3 Resolution of the inverse problem

The inverse problem [START_REF] Bellman | On structural identiability[END_REF][START_REF] Rubinow | Compartment analysis: an inverse problem[END_REF] consists in determining the transfer constants, k ji , from the observations (u and y). The transfer constants are the speed of transfer, transport or reaction between two compartments.

The resolution of ODEs usually yields to sum of two exponential functions if the input function of the system is an impulse. If the system is solved analytically, then we obtain a sum of exponential functions parameterised by the constants of the system. Thus, if we approximate the observations by multiexponential functions, this will allow us to deduce the value of the where the dierent matrices refer to: X(s) = X 1 (s)

k 21 D 1 λ 1 -λ 2 ( 1 s-λ 1 -1 (s-λ 2 ) )    (B.11)
Finally, it is possible to obtain back the output temporal functions y 1 (t) and y 2 (t) with the inverse Laplace ( L -1 ). Nous élaborons un modèle compartimental du transport de l'iodure dans l'estomac.

Nous dénissons les coecients de mélange qui nous permettent de calculer l'activité dans chaque compartiment à partir de l'activité mesurée dans chaque couche.

Nous mettons en place une méthodologie combinant des étapes de résolution analytiques et numériques permettant la caractérisation complète du cycle de prise de iodure.

Enn, nous calculons les coecients de transfert pour notre modèle compartimental. Ces coecients sont les paramètres du modèle compartimental. Les résultats sont très satisfaisants pour une première étude.

E.5 Organisation du manuscrit

Ce manuscrit s'articule autour de deux parties montrant les deux contributions principales:

1. La première partie traite du mouvement respiratoire et de la reconstruction d'images TEMP du petit animal sans artefact de mouvement. Plus intéressant encore, l'idée d'appliquer une d¢omposition en couches

1 6 1. 1

 161 Imagerie TEMP 4D du petit animal -Estimation du Mouvement Respiratoire et de la Biodistribution de l'Iode Résumé: L'objectif de cette thèse est d'étudier temporellement des phénomènes évolutifs à l'aide de la tomographie d'émission monophotonique (TEMP). La première partie de cette thèse traite le problème du mouvement respiratoire dans les images TEMP de souris. Nous présentons ici une méthode permettant de détecter ce mouvement respiratoire dans les images TEMP 4D, d'extraire un signal respiratoire intrinsèque, et de déterminer la phase du cycle respiratoire sans mouvement la plus large possible. Les données enregistrées durant ces phases sans mouvement sont alors utilisées pour reconstruire une seule image TEMP 3D sans artefact de mouvement par acquisition. Les images ainsi reconstruites présentent un bon compromis en terme de statistiques et de précision des mesures par rapport aux images TEMP 3D de base et TEMP 4D. Dans la deuxième partie, nous étudions la cinétique d'accumulation de l'iode dans l'estomac de souris par le biais de l'étude de la biodistribution du pertechnetate de technetium ( 99m Tc-pertechnetate) avec des images TEMP 4D. An de comprendre le rôle biologique de cette accumulation dans l'estomac, nous avons modélisé le phénomène par une approche d'analyse compartimentale avec un modèle simpliée à deux compartiments (paroi et cavité stomacale) et une entrée (sang). Les courbes temps -activité (TAC) de chaque compartiment sont déduites des observations et une première estimation des paramètres a été obtenue. Mots-clés: TEMP du petit animal; Images dynamiques; Mouvement respiratoire; Synchronisation respiratoire; Context and Motivations: Longitudinal studies and applications to small animal imaging . . . . . . . . . . . . . . . . . . . 3 1.2 Material and experiments . . . . . . . . . . . . . . . . . . . . . 4 1.3 Respiratory motion . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.4 Dynamic distribution of iodide in stomach . . . . . . . . . . . 5 1.5 Manuscript organisation . . . . . . . . . . . . . . . . . . . . . . Context and Motivations: Longitudinal studies and applications to small animal imaging

Figure 2 . 1 :

 21 Figure 2.1: Presentation of pre-clinical imaging system eXplore speCZT CT 120: (a) Description of the assembly of animal table SPECT and CT imaging systems, (b) Description of the pinhole collimator operating, (c) Schema describing CT acquisition process, (d) Description of SPECT acquisition process (pinhole collimation).

  Contrary to adult humans, animals cannot stand still in a specic position for a certain period of time (acquisition), except died or asleep. When imaging alive subjects, anaesthesia is used for articially send animals to sleep. It avoids animal motion and allows to control the breathing rhythm. Awake healthy wild type mice breathing rate ranges from 84 to 230 breath per min (equivalent to 0.26 to 0.71 seconds per cycle) while for immuno-depressed strains, respiratory rate ranges from 177 to 420 breath per min (equivalent to respiratory cycle length ranges from 0.14 to 0.34 sec) [SZI00]. Martiniova et al. [MSL + 10] distinguished dierent anaesthetic alternatives: intra-peritoneal or intra-muscular injection of mixture of ketamine and other drug, or isourane inhalation. When animal are anaesthetised, their respiratory rhythm generally slows down.

(

  at the time of the acquisition) [NER + 03, VLB + 07]; motion compensation, which consists in incorporating a motion model in the image reconstruction algorithm [RMK + 07]. A variant of the motion gating approaches, called motion correction, has been characterised as a third technique. Motion correction works on dynamic or three spatial + 1 temporal dimensionss (4Ds) image. It consists in regis-

4

 4 Image statistics . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 4.3.5 Qualitative comparison . . . . . . . . . . . . . . . . . . . . . . . 54 4.3.6 Quantitative comparison . . . . . . . . . . . . . . . . . . . . . . . 55 4.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 4.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 This chapter presents in its entirety the following journal paper: M. Breuilly, G. Malandain, J. Guglielmi, R. Marsault, T. Pourcher, P.R. Franken and J. Darcourt. Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT. Physics in Medicine and Biology, 2013. Chapter 4. Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT

  instance, Kuntz et al. used phase-based gating [KDZ + 10]. However the lit- erature is less abundant than in the clinical eld and mostly concerns µ-

4. 2 . 4 . 2

 242 Gated 4D reconstructionGated 4D reconstruction (G4D) consists in reconstructing a temporal sequence of N 3D images. The number of phases N is set to 15 in the experiments.Using the trigger set in pressure signal, we compute the respiratory cycle average duration, denoted by C, over the whole acquisition. The duration of each reconstruction window is then D = C/N . Therefore the n th gating window of the i th cycle is dened by the temporal interval [t i + (n -1)D, t i + nD[ where t i is the time of the i th trigger and n ∈ [1, N ]. A 3D image for each gating reconstruction in micro-SPECT window is then reconstructed from the data of the corresponding window.

Figure 4 . 1 :

 41 Figure 4.1: Histogram of cycle duration This gure presents the distribution of cycle duration for a pressure signal record, with its Gaussian curve t. Vertical red dotted lines indicate the average cycle duration C (middle) and the two limits C -3σ and C + 3σ that have been dened for the cycle selection. Cycle with a duration outside the limits are discarded for following reconstruction.

Figure 4 .

 4 Figure 4.1 presents a typical example of the distribution of cycle durations for an acquisition. The acceptance set to [C -3σ;C +3σ] will guarantee to keep around 99.7% of the cycle if the distribution is Gaussian, while it will discard the outlier cycles.

gure 4 .

 4 2(b).

Figure 4 . 2 :

 42 Figure 4.2: Reconstruction windows for G4D reconstruction The two graphs represent two samples of the same pressure signal. Vertical lines delimit the reconstruction windows for each phase of the G4D image. Sample of pressure signal with cycles shorter than the average duration C f highlights the overlapping of phases from consecutive cycles (a); Sample of pressure signal with cycles longer than C f highlights the fact that data at the end of some cycle are not used for reconstruction (b). Moreover, to reconstruct a phase of the G4D image (for example phase 15), G4D reconstruction use data that does not correspond to the same time in the respiratory cycle. Sometimes it corresponds to the end of exhalation (a), sometimes it corresponds to the end of inhalation (b).

Figure 4 . 3 :

 43 Figure 4.3: Singularity of anaesthetized mice pressure signal with motion phase detection Pressure signal recorded with the monitoring system (a); and normalized average pressure signal (continuous curve) and normalized histogram resulting from the image-based motion phase detection (dashed curve), the horizontal line depicts the threshold determined by the Otsu method (b). The phases getting more votes (above the threshold) are associated with the motion phase and the ones with less votes (below the threshold) are associated with the motionless phase.

First, all the

  reconstructed images (3D images, as BH3D, or 3D phase images from 4D images, as G4D, G4DSR) were normalised to NG3D. The values of the reconstructed voxels depend directly on the number of counts used in the reconstruction process. NG3D is the 3D image reconstructed with all the available counts, while all other images only uses a fraction of them. Let us denote by N c(I) the number of counts used for the reconstruction of image I. By mutiplying the values of image I after reconstruction by the ratio N c(N G3D)/N c(I), the resulting values are comparable with the ones of NG3D. Second, the values are converted into Standardized Uptake Values (SUV). Thanks to the previous normalisation, values of all reconstructed images are comparable and the conversion into SUV is the same for all images. Image values are rst converted into uptake values by calibration using an external known 99mTc source acquired simultaneously with the image. SUV values are obtained by normalising the uptake values by the total injected activity divided by the animal weight.4.2.7 Image measurementsA number of measures, either image-based or lesion-based, were performed for each reconstructed image.A volume of interest was manually delineated in the liver with the help of both anatomical (CT) and SPECT (NG3D) images avoiding any abnormal lesion, yielding a volume of homogeneous SUV. Calculating the standard deviation (SD) over this region allow the noise of the reconstructed image to be estimated.

Figure 4 .

 4 Figure 4.4 presents the volume rendering of fused CT and SPECT (NG3D) images of peritoneal carcinoma lesions after 3 weeks of growth.Figure 4.4

Figure 4 . 4

 44

Figure 4 . 4 :

 44 Figure 4.4: CT/SPECT volume rendering highlighting peritoneal carcinoma lesions. Volume rendering of fused CT and SPECT (NG3D) highlighting peritoneal carcinoma lesions after three weeks of growth, salivary glands and stomach (S).

Figure 4 .

 4 Figure 4.3 (a) presents a sample of pressure signal of free-breathing mice. Table 4.1 presents the results of the pressure signal analysis on the 8 acquisitions. The average duration of the respiratory cycle varies from 658 ms to 766 ms while its standard deviation ranges from 44 to 153 ms. This signal is characterised by a still phase at the end of exhalation. The image-based motionless phase detection demonstrated that the shorter phase associated with motion state lasts 5 or 6 phases out of 15 (see table 4.1) which corroborates our initial observation of a duration of around 1/3 of the respiratory cycle. Moreover, it has been detected that the motion state is time-shifted away from the pressure signal. According to the image-based motionless phase detection, the time shift varies from 3 to 5 phases.

Figure 4 .

 4 Figure 4.3 (b) presents the normalised average respiratory cycle (continuous curve) for one acquisition, after the cycle selection and the cycle resam-

Figure 4 .

 4 6 illustrates these reconstructions with sagittal views showing 2 peritoneal lesions. The peritoneal location can be visualised with the CT image (a). NG3D image (b) shows the underestimation of the lesion uptake and the overestimation of size compared to G4D images, with a signicantly improved signal-to-noise ratio. The zoomed-in images of the G4D end-of-exhalation (c) and end-of-inhalation (d) images show the high lesion uptake and the respiratory induced lesion displacement, but are corrupted by a high noise level. On the BH3D reconstructed image (e), the obtained result is a trade-o between the previous methods:

Figure 4 . 6 :

 46 Figure 4.6: Sagittal views from sub-diaphragmatic peritoneal lesions comparing NG3D, G4D (15 phases), BH3D SPECT reconstruction methods. Top row: anatomical reference CT image (a), NG3D reconstruction (b), end-of-exhalation phase from G4D reconstruction (c), endof-inhalation phase from G4D reconstruction (d), BH3D reconstruction (e); bottom row: zoomed-in on lesions of interest of the SPECT images. The SUV colour map has been set up on SUV max values in G4D images for the selected lesion.

Figure 4 . 7 :

 47 Figure 4.7: Quantitative results for lesions and image characteristics Quantitative results for lesions SUV mean (a), lesions volumes (b), SNR (c), and noise (d). Box plots show the median, upper and lower quartiles and range; + indicate outliers. Statistical signicance: ** for p < 0.01, *** for p < 0.001, p-value is given when non signicant (p > 0.05).

Contents 5 . 1

 51 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.1.1 Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.1.2 Dynamic SPECT acquisition and analysis . . . . . . . . . . . . . 5.1.3 Challenges . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.1.4 Region-based analysis of the 99m Tc-pertechnetate biodistribution 5.2 Biological context . . . . . . . . . . . . . . . . . . . . . . . . . . 5.2.1 Biological objectives . . . . . . . . . . . . . . . . . . . . . . . . . 5.2.2 Biological knowledge . . . . . . . . . . . . . . . . . . . . . . . . . 5.3 Material . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.3.1 Animal model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.3.2 Dynamic SPECT imaging protocol . . . . . . . . . . . . . . . . . 5.4 Method and qualitative results . . . . . . . . . . . . . . . . . . 5.4.1 Segmentation of dynamic SPECT images of stomach . . . . . . . 5.4.2 Layer-based decomposition of the stomach . . . . . . . . . . . . . 5.4.3 Kinetic evolution of activity in stomach layers . . . . . . . . . . . 5.4.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .5.1 IntroductionIn this second part of the thesis, we study the biodistribution of 99m Tcpertechnetate in murine stomach. The idea is to automatically segment an organ in a temporal sequence of three (spatial) dimensions (3D) images. Then, a layer-based decomposition of the segmented organ will get rid of the animal motion. In the information provided by the layer decomposition, we try to integrate notions of functional compartment. Chapter 5. Dynamic SPECT analysis of 99m Tc-pertechnetate uptake in stomach: biodistribution study 5.1.1 Objectives The main objective of this work is to demonstrate the feasibility of the modelling of the 99m Tc-pertechnetate biodistribution with a very simple model. For this purpose, we rst study the 99m Tc-pertechnetate uptake kinetic with three spatial + 1 temporal dimensions (4D) single photon emission computed tomography (SPECT), and then we propose a simplied compartmental model that describes this progressive uptake. The compartments of the model corresponds to real compartment that are biological functional units: the sodium/iodide symporter (NIS) expressing cell in the glandular wall and the mucus secreted in the stomach cavity. The model will help biologists to complement their knowledge about the role of 99m Tc-pertechnetate in extrathyroidal tissues.

Figure 5 .Figure 5 . 1 :

 551 Figure 5.1: Anatomy of a murine stomach [KRK11]. 1 = distal esophagus, 2 = forestomach, 3 = margo plicatus, 4 = pars fundica, 5 = pars pylorica, 6 = proximal duodenum. The dotted line highlight the position of the margo plicatus.

99mFigure 5 . 2 :Figure 5 . 3 :

 5253 Figure 5.2: IHC image of murine stomach.Figure (a) shows a section of the stomach with nucleus (blue) and NIS protein (brown) marking. Arrows indicate the forestomach (1), the esophagus (2), the duodenum (3), the glandular stomach (4) and the margo plicatus (5,6). Figure (b) shows a zoomed-in image on the mucosa. Figure (c) shows a zoomed-in on a transverse section of a follicle of the glandular stomach mucosa. Figure (d) describes the follicle made of 13 epithelial cells: extracellular space (1), apical region (2), cellular space (3), stomach cavity (4), NIS protein marking on basolateral region(5)

Figure 5 .

 5 Figure 5.1 show an ex vivo stomach where we can identify each parts. The nonglandular forestomach has a thin wall and is transparent. The glandular stomach has instead a a thicker wall as shown on the IHC images Figure 5.2(a).These two regions are separated by a limiting ridge or margo plicatus. The glandular stomach does present a thick epithelium. It is divided into three parts: the small cardia adjacent to the limiting ridge (margo plicatus); the fundus, where the mucosa is folded with rugae; and the antrum or pars pylorica, which has relatively smooth mucosa. On Figure5.1, the limiting ridge is emphasized by the white dotted line but it is clearly visible by the dierence of tissue color. On Figure5.2(a), the limiting ridge is the visible, slightly raised division at the left hand side junction with esophagus.

99m

  Tc-pertechnetate uptake in stomach: biodistribution study biodistribution analysis. Although biodistribution are done on whole-body images, here we will only consider the cardiac and abdominal regions.

5. 4 . 1 Figure 5 . 4 :

 4154 Figure 5.4: Structuring element for the 2D morphological closing in order to segment both stomach wall and cavity, the 3D structuring element can be easily deduced: each voxel at an Euclidean distance less or equal to 2 voxels is associated with the value 1.

Figure 5 . 5 :

 55 Figure 5.5: Distance transformation in 3D in a 5x5x5 neighbourhood. These gures illustrate the distance computation in a 5x5x5 neighbourhood. The distance is given with respect to the inner voxel at level 0. Values are given for the planes at level 0, +1 and +2. Distance in planes -1 and -2 are obtained by symmetry. (a) presents the Euclidean distance transformation from the inner voxel. (b) presents the 5-7-11 Chamfer distance transformation which is the a 3D integer approximation for the Euclidean distance before division by 5.

for the Chamfer 5 - 7 -Figure 5 . 6 :

 5756 Figure 5.6: Stomach for the subject 0371 at early time after 99m Tcpertechnetate injection (a), and the resulting distance transform (b). Background (outside of the stomach) has a zero distance value. Foreground (stomach) distance increase with the distance to the background.

The 4 th

 4 layer is located at the distance 11/5 = 2.2 from the stomach border which is a good estimate for √ 5 2.24. The distances can thus be expressed in voxel unit, that represent the distance between voxel centers.An example for the Chamfer distance transformation is given in Figure5.6.

Figure 5 .

 5 Figure 5.6(a) shows the image at the rst time point, 11.5 min after 99m Tcpertechnetate injection, of the 4D image for subject #0371 centered in the cavity of the stomach. The stomach can be distinguished by the layer of activity that surrounds the cavity without acivity at that time. After segmentation of the wall and a morphological closing step, the Chamfer distance transformation is computed.Figure 5.6(b) shows the resulting distance 3D

Figure 5 . 4 ?

 54 6(b) shows the resulting distance 3D map. The distance zero is attributed to all voxel from the outside of the stomach (background) and starting from the outer layer of the segmented stomach, the distance increases according to the Chamfer distance transformation.99m Tc-pertechnetate uptake in stomach: biodistribution study5.4.3 Kinetic evolution of activity in stomach layersThe spatial and temporal evolution of the stomach activity uptake are simply observed with respect to the layers dened by the Chamfer distance transformation previously described. right 2.(a) Spatial activity boxplot at rst time acquisition: T0 + 12.5 min. Axes: Distance from stomach border (left) to inner stomach (right) in abscissa, image intensity in (%Id/cc) in ordinate (b) Time activity boxplot for layer at distance 2.4 from stomach border. Axes: Time after injection (min) in abscissa, image intensity in (%Id/cc) in ordinate (c) Time activity curves for each layer. Axes: Time after injection (min) in abscissa, average SUV in ordinate

Figure 5 . 7 :

 57 Figure 5.7: Example of spatial-activity and time-activity observations for subject #0371. Distance-activity curves The spatial activity curves are observed for each time point of the dynamic images. Each layer that has been segmented, either in the stomach wall or cavity, is considered. An example is given for the subject 0371 in Figure 5.7(a). The spatial variation conrms the visual assessment: the activity increases from the border up to small distances (around 2.4 voxels) and then decrease up to a nearly zero activity in the cavity. (See Section 5.4.2 for details on the distance computation.) According to the Chamfer distance transformation Figure 5.5, distances larger than 4.5 voxels should all be gathered in a single thick layer associated to stomach cavity. This question will be discussed in Section 5.5.

Figure 5 .

 5 Figure 5.10 exhibit the interesting segmentation for the subject #0371.

Figure 5 . 8 :

 58 Figure 5.8: Stomach segmentation for subject #0366 at several time point, same transverse slice

Figure 5 . 9 :Figure 5 . 10 :

 59510 Figure 5.9: Stomach segmentation for subject #0370 at several time point, same transverse slice

Figure 5 .Figure 5 . 12 :Figure 5 .

 55125 Figure 5.11: Stomach segmentation for subject #0396 at several time point, same transverse slice

5. 4 . 4 . 2

 442 Stomach segmentation and Chamfer distance transformation (a) Original transverse slice (b) Result of the segmentation superimposed over the normalised image after Gaussian ltering (σ = 1.00) (c) Result of the Chamfer distance transformation

Figure 5 .

 5 Figure 5.14: Temporal evolution of iodide uptake after intraperitoneal injection of 99m Tc of subject 0370. From left to right: T01: T0 + 16.5 min, T07: T0 + 43.5 min, T14: T0 + 74.5 min, T21: T0 + 106.5 min, T22: T0 + 111.5 min.

Figure 5 .

 5 Figure 5.14 presents the result of the Chamfer distance transformation

Figure 5 .

 5 Figure 5.15: Activity boxplots with respect to the distance from stomach wall to stomach cavity for mouse #0530. Axes: Distance from stomach border (left) to inner stomach (right) in abscissa, image intensity in (%Id/cc) in ordinate.

Figure 5 .

 5 Figure 5.17 shows the variation of the threshold that was automatically computed for segmenting the whole stomach. The sub-gure on the left shows that the threshold grows along time, similarly to the global activity. The subgure on the right conrms the trend that appears on the left hand side: the variation of the automatic threshold grows more or less at the same rate as the maximum intensity. The conclusion is that some of the external layers may be lost along time, because the global activity grows. The external layers that

Figure 5 .

 5 Figure 5.16: Mouse #0530: time activity curves per layer and in the whole stomach

Figure 5 .

 5 Figure 5.17: Mouse #0530: segmentation threshold variation through the time.

99mFigure 5 .

 5 Figure 5.18: Plasma temporal activity curves. Axes: Time after activity injection (min) in abscissa, average SUV in ordinate.

  In this chapter, a mathematical model based on compartmental analysis approach is proposed. Compartmental analysis is frequently used in pharmacology to estimate drug exchanges in organs, or in biomedical emission tomography (ET) to study dynamically the biodistribution of markers. So far, most biodistribution studies considered organs as a whole but only a few of them considered organs as a set of dierent functional subunits (subregions or tissues). Here both the size of the stomach and the camera resolution were suitable for the distinction between of the stomach wall and the stomach cavity.More generally, compartmental analysis consists in describing phenomena as compound exchanges between compartments. Each compartment is supposed to have a homogeneous concentration of the studied compound but they distinguishes from each other from the uptake kinetic. Exchanges of compound between compartments are described by transfer functions. The ordinary dierential equations that describe the transfer functions are linear 99m Tc-pertechnetate uptake in stomach: biodistribution modelling functions with time-invariant parameters. The bases of compartmental analysis are summarised in Appendix B. We also remind the major denitions and principles that can be found in the literature (ref. to[START_REF] Godfrey | Compartmental models and their application[END_REF][START_REF] Jacquez | Compartmental analysis in biology and medicine[END_REF] for further information). The reading of this appendix might be useful for the understanding of the mathematical description of the compartmental analysis.The compartment analysis used for the characterisation of a biodistribution model can be considered in two ways. Either by solving a direct problem that consists in simulating the kinetic evolution from the knowledge of the transfer functions of a theoretical model. Or, as it is the case here, by solving an inverse problem that consists in estimating the transfer functions from longitudinal observations and characterise the model parameters[START_REF] Bellman | On structural identiability[END_REF][START_REF] Rubinow | Compartment analysis: an inverse problem[END_REF].

  ing and dierent radiopharmaceuticals: on the thyroid [Hay78, DCH + 08, FGV + 10], on the heart [RGH98, HRZG98, RGH00, KG01, RGH02, GRS + 10, ZGK10, NYK12] and on the brain [JST + 12, NSG + 13].

A

  Figure 6.1:

Figure 6 . 2 :

 62 Figure 6.2: Complete compartmental model describing the studied problem:

Figure 6 .

 6 Figure 6.2 presents a rst proposition for the compartmental system that describes the dynamic distribution of 99m Tc-pertechnetate in the murine body.

12 Figure 6 Figure 6 . 3

 12663 Figure 6.3: Compartmental system after second simplications: twocompartment catenary system

Figure 6 . 4 :

 64 Figure 6.4: Blood time-activity curve (red continuous line) and its approximate function (blue dotted line) for subject #0530. The next step is the determination of the unknown constants a,b, c, and dwith least square minimisation algorithm between the observation y 0 and its estimate ŷ0 . Finally, it is possible to retrospectively determine the unknown constants φ 1 , φ 2 , θ 1 , and θ 2 that are combinations of the previously determined constants a, b, c and d.

Figure 6 .

 6 Figure 6.4 presents the observed blood time-activity curve for the subject #0530 (cf Table 5.1), and its multiexponential approximate function. The resulting coecients are φ 1 ≈ -1.93, φ 2 ≈ 2.24, θ 1 ≈ 0.123, and θ 2 ≈ 0.0103. The squared residuals between the observations and the estimate function at the observed time is approximately equal to 0.0826. The root mean square (RMS ) is equal to RMS = SSD 19 = 0.0659.

Figure 6 . 5 :

 65 Figure 6.5: Stomach layer TAC for all the distances from the outer layer of the segmented stomach (d = 1.00) to the inner layer for subject #0530.

3

 3 TACs for the stomach wall and cavityNow, we propose an estimation of the mixing coecient matrix K, where each coecient has been chosen according to the estimated volume ratio of each compartment in each layer.

Figure 6 .

 6 6(a) presents the resulting TACs for the two compartments: stomach wall TAC W and stomach cavity TAC C compared to the blood TAC B .

Figure 6 .

 6 6(b) presents the observed TACs for each layer in order to give a reference for the amplitude of activity variation.

Figure 6 .Figure 6 . 6 :

 666 Figure 6.6: Stomach wall and cavity TACs for subject #0530.

Figure 6 . 7 :

 67 Figure 6.7: Stomach wall and cavity TACs multiexponential approximation.

Figure 6 .

 6 Figure 6.7(a) and (b). The results is clearly very far from the original curves.

  85% of error for the latest time frame for TAC C . Moreover, the RMS error on the layer TACs is around 40 that represents around 57% of error for the latest time frame. The estimated TAC W and TAC C are also qualitatively compared to the approximated ones. As shown on Figure 6.8(a), the estimated TAC W and TAC C obtained with the rst set of parameters of case (A) does not t very well the approximated TAC W and TAC C .

  nb of layers×nb of frames = equivalent SUV error/layer/time frame.

Figure 6 . 8 (= 10 )(see Figure 6

 68106 a) presents the graphic for case (A1); (B): Best transfer parameters with respect to column (1) and results; (C): Best transfer parameters with respect to column (2) and results. (D): Best parameters and results with respect to column (3) after the rst execution of the downhill simplex algorithm (ε = 1000) (see Figure 6.8(b)); (E): Best parameters and results with respect to column (3) after the second execution of the downhill simplex algorithm (ε

Figure 6 . 8 :

 68 Figure 6.8: Compartments TACs and their approximates functions. These gures present the compartment TACs:TAC W for the wall (blue con-

Figure 6 . 9 :

 69 Figure 6.9: Dierence (signed error) and relative error between observed and estimated layer TAC for the rst 14 layers and for all time frame for subject #0530.

  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 7.2 Future works and perspectives . . . . . . . . . . . . . . . . . . 129 7.2.1 Direct improvements of the contributions . . . . . . . . . . . . . 129 7.2.2 Long term research directions . . . . . . . . . . . . . . . . . . . . 130

7. 2 . 1

 21 Direct improvements of the contributions Software for breath-hold like image reconstruction The software dedicated to the reconstruction of SPECT images without respiratory motion artefacts should be further developed. Although it has been developed for biologists, its current version is just a beta version and is not yet adapted to a routine use.

A. 2 . 1

 21 Study of the reproducibility of the respiratory signal A preliminary simple test was the observation of the dierent average respiratory cycles according to dierent criteria: animal lineage, disease, subject, in order to to see if there is a rule that describes the variability from one acquisition to another (Figure A.1, Figure A.2, Figure A.3, Figure A.4, and Figure A.5).

First, we compared

  Figure A.1: Comparison of average respiratory cycle according to the mouse lineage, for several animal and for several days of acquisition.

Figure A. 2 :

 2 Figure A.2: Comparison of average pressure cycle obtained on repeated acquisitions of the same animal (C57Bl/6 lineage).

Figure A. 4 :

 4 Figure A.4: Comparison of average pressure cycle according to the day of acquisition (NOD-SCID lineage).

Figure A. 6 :

 6 Figure A.6: Comparison between center of mass motion (along each axes: leftright, dorsoventral and craniocaudal, and the Euclidean (or 2-norm) distance) and average respiratory cycle for one lesion in mouse #022, with respect to the rst 50 ms phase.

(

  see Figure A.6 and Figure A.7). The graphics show that the lesion is not moving along the left-right axis (Figure A.6(a)), while it is slightly moving along the dorsoventral axis (Figure A.6(b)) with a shift above 1 voxel, and signicantly moving along the craniocaudal axis with a cumulative shift above 4 voxels. The Euclidean distance of the lesion motion shows a shift of 4 voxels from the initial position.

Figure A. 7

 7 Figure A.7 presents the Euclidean distance illustrating the motion of several lesions. It suggests that all the lesions are animated by a similar motion, with a range of motion up to 4 voxels from the initial position. Two observations should be noted from these graphics.

Figure A. 8

 8 Figure A.8 presents the average respiratory cycle, the lung volume variation, and the movement of the two points identied above for three distinctive acquisitions. First, the graphics show that the three intrinsic signals have a similar variation. This suggests that the lung volume variation is related to the diaphragm motion, and that the diaphragm motion only could be used as the

Figure A. 9 Figure A. 9 :

 99 Figure A.9 shows the variability in the distribution of cycle lengths for each acquisition together with the estimated Gaussian distribution law (estimated with three dierent methods). Some of the acquisitions present two modes in the distribution, meaning that the cycle lengths did vary during the acquisition around two dierent average lengths as we can see in Figure A.9(a), Figure A.9(b), and Figure A.9(c). For one acquisition in particular (Figure A.9(f )), the distribution of cycle length was particularly narrowed around the average. Concerning the other acquisitions, the distribution follows more or less a Gaussian distribution. Though the rst estimate (red curves in gures) include the distributions more or less nicely. Additionally to this visual observation, the cycle length distributions were also numerically

Figure A. 9 .

 9 We wanted to nd a law that describes nicely the distribution of cycle lengths, but actually, if the mice are breathing at the same rhythm at all times, the distribution should be an impulse function at a specic cycle length. On the contrary, the cycle lengths vary. In the previous paragraph, we suggested to perform a cycle selection step. However, the selected range (±3 * σ) still leads to a relatively high variability. This can be explained visually with the observation of the respiratory signal. Each cycle is distinguished from the previous and next one by trigger indices. This trigger is recorded simultaneously in the respiratory signal and in the listmode. All the respiratory cycles can be superimposed on one gure per acquisition by aligning them according to the trigger indices at the beginning of cycles. Figure A.10 presents the superposition of all cycles before and after the cycle selection for one acquisition of the mouse #033. These gures show the variability of both cycle duration and amplitude. When all cycles (around 1400) are added, the resulting signal is very noisy due to cycles heterogeneity (Figure A.10(a)). After outlier cycles rejection, the resulting signal becomes much (Figure A.10(b)) To summarise, the benet of the cycle selection step is illustrated by the fact that respiratory cycles vary similarly in the same window of amplitude. All breathing phases are positioned around the same instant at the beginning of the cycle with respect to the trigger.

A. 3 . 2

 32 Figure A.10(b)is unchanged when resampling relatively to cycle length. The triggers have been determined from an absolute threshold that has been xed for the whole acquisition. It is not adapted to the amplitude variation. Consequently, each trigger does not correspond roughly to the same phase of the cycle; xing the trigger at 20% of the amplitude variation between top and bottom amplitude of the cycle, if it was possible, would have given a more accurate result. In other words, we could have worked with a relative amplitude threshold. Unfortunately, setting a relative trigger is not possible with our system. Since both record (respiratory signal and listmode) are synchronised thanks to the common triggers, an additional improvement could be considered: the a posteriori redenition of the rising edge detection threshold on the respiratory signal. Instead of a xed value that is identical for all cycles, the threshold would be adapted to each cycle. This threshold would be a relative position on the excursion of pressure between inhalation and exhalation. Once redened, a new triggering signal could be regenerated and modied in both respiratory signal record and listmode. However,the inuence of the adjustment of the trigger according to the amplitude variability is quite low since the trigger will be shifted from only a few ms only. Considering the changes required to modify the triggering method and the low impact this aspect should have on the quality of the results, we kept the xed trigger method for delimiting cycles. Discrepancy between triggering signals In this section, we investigate the discrepancy between the two triggering signals that are recorded during each acquisition. Triggering signals are simultaneously recorded in the respiratory signal and in the listmode. Even though, the acquisition of the two signals is done independently, we assume reconstruction in micro-SPECT: Supplementary data

  (a) Comparison of cycle length of matching cycles. Cycle shift has been determined with LS minimisation (b) Comparison of cycle length of matching cycles. Cycle shift has been determined with outlier minimisation

Figure A. 11 :

 11 Figure A.11: Comparison between cycle lengths of triggering signals embedded in the respiratory signal and the listmode after cycle shifting, for acquisition #0182.

  (a) Cumulative dierence on cycle length between the two signals after cycle shifting (b) Number of millisecond of dierence between respiratory signal and listmode per 10 cycles after cycle shifting.

Figure A. 13 :

 13 Figure A.13: Dierence between cycle lengths after registration with outlier minimisation method of the triggering signals, for two acquisitions 0182

  Figure A.11(b) and Figure A.12(b) show the pairwise comparison of cycle length for the same two acquisitions, with the new shifts. For acquisition 1839, most cycle lengths now match. The cumulative pairwise dierence (respiratory signal -listmode) for both acquisitions 0182 and 1839 is additionally presented in Figure A.13(a) and Figure A.14(a).

  Figure A.14(b) show the sum of pairwise dierence per sample of 10 cycles.

  Figure A.17: Comparison between the average pressure cycle and motion detection (detailed method in Chapter 4 Section 4.2). The comparison was conducted on four consecutive acquisitions with dierent gating parameters: gating on falling edge (a), gating on falling edge and delay (b), gating on rising edge (c), gating on rising edge and delay (d)

X 2 5 ) 1 = 1 ∆ 12 (B. 8 ) 2 ) where λ 1 λ 2 are dened by: λ 1 , λ 2 = 1 2 -(B. 9 )D 1 ( 1 s-λ 1 + 1 s-λ 2 )

 25111282212911112 (s) , Y(s) = Y 1 (s) Y 2 (s) , U(s) = U 1 (s)The observation matrix Y(s) can be also expressed as:Y(s) = G(s)U(s) (B.6)Since both B is assumed to be identity, we have:G(s) = (sI -A) -1 B = (sI -A) -1 = s + (k 21 + k 01 ) -k 12 -k 21 s + (k 12 + k 02 ) -(s) s + (k 12 + k 02 ) k 12 k 21 s + (k 21 + k 01 ) (B.7) ∆(s) = |sI -A| = s + (k 21 + k 01 ) -k 12 -k 21 s + (k 12 + k 02 ) = s 2 + (k 21 + k 01 + k 12 + k 02 )s + (k 21 + k 01 )(k 12 + k 02 ) -k 21 k It is possible to factorise the determinant ∆(s) as (s -λ 1 )(s -λ (k 21 + k 01 + k 12 + k 02 ) ± (k 21 + k 01 + k 12 + k 02 ) 2 -4(k 21 k 02 + k 01 k 12 + k 01 k 02 )At last, general expression for G(s) and Y(s) are given:G(s) = s+(k 12 +k 02 ) (s-λ 1 )(s-λ 2 ) k 12 (s-λ 1 )(s-λ 2 ) k 21 (s-λ 1 )(s-λ 2 ) s+(k 21 +k 01 ) (s-λ 1 )(s-λ 2 ) (B.10) Y(s) = G(s)U(s) = (s+(k 12 +k 02 )) (s-λ 1 )(s-λ 2 ) D 1 k 21(s-λ 1 )(s-λ 2 ) D 1 λ 1 +(k 12 +k 02 ) λ 1 -λ 2 λ 2 +(k 12 +k 02 ) λ 2 -λ 1

1 1 . 1 B. 4

 1114 Figure B.2 presents the dierent functions (inputs and outputs) for this example.

  Figure C.2: Sagittal (a) and coronal (b) planes

2 .

 2 La deuxième partie décrit la modélisation de la biodistribution de l'iodure dans l'estomac murin. Le chapitre 2 décrit les contextes biologique et technique de cette thèse. Il décrit aussi les contraintes inhérentes au matériel et aux modèles animaux. La partie I se concentre sur le mouvement respiratoire et se compose de deux chapitres. Le chapitre 3 présente un état de l'art sur les méthodes de correction du mouvement respiratoire. Le chapitre 4 reprend in extenso l'article [BMG + 13] dans lequel nous proposons un procédé dédié à la recon- struction optimale sans artefact de mouvement d'images TEMP 3D. La partie II porte sur la biodistribution d'iodure dans l'estomac murin et se compose de deux chapitres. Dans le chapitre 5, nous décrivons de manière qualitative la biodistribution du pertechnétate de technétium 99m Tc, analogue de l'iodure, dans l'estomac, à partir de l'observation d'images TEMP 4D. Dans le chapitre 6, nous caractérisons le phénomène d'absorption de l'iodure à l'aide d'un modèle compartimental simplié. Enn, le chapitre 7 conclut cette thèse avec la liste des contributions et présente diérentes perspectives possibles pour ce travail. Une liste des publications est fournie page 133. Les annexes les plus utiles de cette thèse sont les suivantes. L'annexe A donne des résultats supplémentaires à propos du procédé de reconstruction d'images TEMP 3D sans mouvement. L'annexe B décrit les principes de l'analyse compartimentale. région abdominale. Cette contribution a conduit à l'élaboration d'une première version d'un logiciel mettant en oeuvre cette méthode. Ce logiciel avec interface graphique permet aux biologistes d'obtenir une reconstruction d'images TEMP en quasi apnée. Etude de l'absorption du pertechnétate de technétium 99m Tc dans l'estomac An d'aider les biologistes à mieux comprendre le processus d'absorption et la sécrétion de l'iodure dans l'estomac, nous proposons une méthode qui permet de construire un premier modèle compartimental décrivant le cycle de la prise de pertechnétate de technétium 99m Tc, analogue du iodure, dans l'estomac murin. Ceci a conduit à la conception d'un modèle compartimental, avec une première estimation des coecients de transfert entre les compartiments. La méthode repose sur une décomposition originale de l'estomac en couches qui nous permet d'obtenir une information anatomique pertinente de l'estomac, et qui compense l'absorption progressive du pertechnétate de technétium 99m Tc et les irrégularités de forme de l'estomac. En outre, cette description en couches est robuste à la déformation générale de l'estomac. Néanmoins, plusieurs dés ont dû être résolus. Un premier dé concerne le processus de segmentation de l'estomac basée sur l'intensité a priori sans information anatomique. Cette étape de segmentation a échoué pour les images aux premiers instants (juste après injection d'activité) puisque le site d'injection, qui était dans le champ de vue de l'acquisition a été segmenté à la place de l'estomac. Un deuxième dé est que chaque couche ne correspond pas seulement à un compartiment. Chaque couche est un mélange de plusieurs compartiments, et même si chaque compartiment est supposé être homogène, une couche ne l'est pas. Nous avons donc considéré l'ensemble des couches comme un mélange de compartiments, chaque coecient de mélange étant représentatif de la composition de chaque couche. Une première détermination de ces coecients de mélange a été proposée. Cela a conduit à de très bons résultats au vu de la méthode utilisée pour les déterminer. Quelques autres dés n'ont été que partiellement résolus, comme la robustesse de la méthode, ou la composition non-uniforme de chaque couche. En outre, la méthode est basée sur plusieurs hypothèses simplicatrices. Par exemple, nous ne prenons pas en compte les déformations locales de l'estomac au cours du temps, ce qui engendre un léger décalage entre les couche externes. La première couche dans une image d'une séquence donnée peut disparaître dans l'image suivante de la même séquence. Malgré ces points qui peuvent être améliorés, la méthode donne des résultats pertinents en termes d'estimation des courbes temps-activité.

F. 2 F. 2 . 1

 221 Travaux futurs et perspectives Le développement de ces deux méthodes ouvre un large éventail de possibilités. Nous présentons d'abord plusieurs améliorations directes de nos méthodes et les résultats qui pourraient être envisagés. Ensuite, nous proposons quelques axes de recherche à long terme. Améliorations directes des contributions Logiciel pour la reconstruction d'image en quasi-apnée Le logiciel dédié à la reconstruction d'images TEMP sans artefact de mouvement respiratoire doit être davantage développé. Bien qu'il ait été développé pour les biologistes, sa version actuelle est juste une version bêta et n'est pas encore adaptée à une utilisation en routine. Validation de la méthode de reconstruction d'images TEMP en quasi-apnée sur des données simulées Cette étape de validation consiste à simuler des images TEMP à partir d'un fantôme de souris incluant un modèle de la respiration, an d'obtenir un plus grand nombre de données. Ces fantômes numériques de souris devraient présenter de petites lésions situées aléatoirement dans l'abdomen et un rythme de respiration variable. Plus précisément, il faudrait animer un fantôme numérique de la souris avec le même motif de la respiration que celui qui a été observée (halètements suivis d'une période d'apnée). Cependant, un tel modèle n'a pas encore été développé. Par exemple, le fantôme MOBY [STF + 04] a été implémenté avec un mouvement respiratoire semblable à celui de l'homme, et ne serait pas susant. En outre, un tel travail exige l'élaboration d'un modèle numérique de notre caméra an de simuler les processus d'acquisition et de reconstruction. L'inconvénient de cette étape de validation est qu'elle implique un ensemble de données qui seront parfaitement produites par le logiciel de simulation sans reéter la réalité de la manipulation et de l'acquisition à partir d'animaux vivants. Cependant, une telle simulation orira un vaste ensemble de données qui pourra être utile pour la validation de notre méthode.Etude de l'absorption d'iodure dans l'estomac La plupart des travaux futurs liés à ce sujet ont été détaillées dans Section 6.5.6. Parmi eux, les points principaux sont l'étude de la robustesse de l'approche; l'évaluation des coecients de mélange de coupes d'immunohistomchimie (IHC), et l'interprétation biologique du phénomène.

F. 2 . 2

 22 Directions de recherche à long termeJe vous présente ci-dessous quelques objectifs de recherche à long terme. Ceuxci consistent à améliorer les méthodes en suivant la même direction que chacune des contributions présentées, mais également en combinant les deux approches an d'étudier de nouvelles directions de recherche. Quantication et modélisation de la croissance de tumeurs abdominales L'étude des tumeurs abdominales, notamment la quantication et la modélisation de leur croissance a été le premier objectif de cette thèse. Cependant, cet objectif a nécessité des travaux préliminaires sur la correction du mouvement respiratoire. Ce travail préliminaire s'est avéré être plus dicile que prévu et a été abordé dans les chapitres 3, 4 et dans l'annexe A de cette thèse. Grâce à cela, le suivi et la quantication précise des tumeurs dans la région abdominale reposera sur des bases solides. Un algorithme d'appariement de la tumeur [WGR + 01] pourrait être développé an de suivre les tumeurs au cours du temps en prenant en compte leur croissance, la possibilité que des tumeurs proches fusionnent, ou l'apparition de nouvelles lésions. Ce type d'approche devrait inclure des méthodes de recalage basé sur des algorithme d'appariement de blocs [ORPA00]. Ensuite, la croissance de tumeur dans la cavité péritonéale peut être étudié et permettre la mise en place d'un modèle informatique de la croissance de tumeur [KPCA08] dans un environnement déformable. Un tel modèle devrait grandement aider les biologistes dans la compréhension de la croissance de ce type de tumeur. Décomposition en couches appliquée à l'hétérogénéité des tumeurs La description originale de l'estomac présenté dans les chapitres 5 et 6 pourraient être étendue à d'autres organes, mais aussi à des tumeurs. En eet, certaines tumeurs ne peuvent être observées qu'à un stade avancé. A un stade tardif, la segmentation des tumeurs à partir d'images TEMP révèle seulement la partie active de la tumeur, constituée par des cellules vascularisées. Cependant, la segmentation ne permet pas de révéler la partie nécrosée de la tumeur. Dans le cas de tumeurs sous-cutanées, la partie nécrosée peut être observée avec l'imagerie TDM car celle-ci s'aaisse. Dans le cas de cellules tumorales exprimant la protéine NIS, la nécrose de la tumeur pourrait être prédite par une diminution de la prise de pertechnétate de technétium 99m Tc dans les couches internes de la tumeur au cours du temps. Plus généralement, il a été également mentionné dans [UTT + 12] que l'intérieur des tumeurs n'est jamais homogène, pour les raisons mentionnées ci-dessus. Dans cet article, les auteurs présentent les résultats de premières enquêtes eectuées avec une imagerie TEMP à grande résolution du petit animal: ils évaluent l'hétérogénéité des tumeurs.

Table 2 .

 2 1: Dierent mice lineage with its immune system characteristics compared to the ones of wild type mice (healthy lineage). = stands for a lym-

	phocyte activity equivalent to wild type,	(respectively	) for a reduced
	(respectively increased) activity of lymphocytes compared to wild type.

cells are trapped, and orthotopic tumour are going to develop. For biological interest and also because it is easier to proceed without surgeon intervention, an intraperitoneal injection has been chosen and peritoneal tumours are going to be studied.

Table 3

 3

	.1 summarises

  Reconstruction is then based on a selection of the acquired data. Such signals can be acquired by an

	external device and include electrocardiogram [LRS + 06], Real-time Position
	Management [GRS	+ 09] or Multidimensional Respiratory Gating [NHAQ + 11],
	pressure sensor (chest belt) [CCP + 10, vHJ + 11], and optical ber [BLW	+ 12].
	These signals can also be extracted directly from the data itself, by looking
	for time variation in raw projections in CT [KM11] or in ET [BMME	+ 08,
	SOBM09].

For periodic signals, phase-based gating consists in dividing the period into several phases of equal duration, with the assumption that the motion can be neglected within each phase. This method is used for both cardiac [BED + 10] and respiratory [BPD + 10, BMME + 08, GRS + 09, LRS + 06] gating for human imaging and allows the reconstruction of dynamic or 4D (3D+t) images. Using only a fraction of the acquired data deteriorates the signalto-noise ratio (SNR) of each 3D image. Therefore, Dawood et al. propose to co-register each 3D image and then to sum them up [DLJS06]. With this method, obtained 3D images have minimal motion artefacts and improved SNR. Furthermore, amplitude-based gating assumes that the signal amplitude is representative of the motion of interest. A range of amplitude values is selected, and the associated acquired data are used for reconstruction. For instance, images at end of exhalation or inhalation can be reconstructed by selecting either the lowest or the highest amplitude values. In human imaging, such a method has been used for cardiac motion in SPECT [KIK + 07] and for respiratory motion in PET [BPD + 10, CCP + 10, NHAQ + 11, SOBM09, vHJ + 11].

  After rejection of the outlier cycles, a new respiratory cycle average duration C s is calculated and leads to a new duration of gating window D s = C s /N . The n th gating window of the i th cycle is dened by [t i + (n -1)D s , t i + nD s [ where t i is the time of the i th trigger.

Table 4 .

 4 1: Respiratory signal analysis for G4D (15 phases) reconstructions.

	Acquisition	Cycle µ ± σ Motion duration Time shift
		(ms)	(ms/nb phases)	(ms/nb phases)
	0182	762.8 ± 112.9 255 (5)	153 (-3)
	0183	740.5 ± 152.7 294 (6)	147 (-3)
	0185	751.0 ± 48.9 250 (5)	200 (-4)
	0186	765.6 ± 43.6 255 (5)	204 (-4)
	0187	712.5 ± 44.6 240 (5)	192 (-4)
	0189	702.9 ± 100.4 282 (6)	188 (-4)
	0190	658.5 ± 52.1 220 (5)	220 (-5)
	0191	713.8 ± 59.1 240 (5)	192 (-4)
	4.2.6 Image value normalisation

Reconstructed images were converted into Standardized Uptake Values (SUV) for comparison purpose. The conversion is divided in two steps.

Table 4 .

 4 2: Percentages of counts used for the dierent SPECT images recon-

	struction schemes: NG3D, G4D (15 phases) and BH3D. Values are presented
	as mean ± SD.
		Number of counts used (%)
	NG3D	Total	100.00
	G4D	Used counts / phase	6.65 ± 0.19
		Total	99.76 ± 2.88
		Omitted counts	3.64 ± 3.14
		Counts used twice	3.51 ± 1.55
		Percentages of counts used twice per phase
		Phase 1	31.44 ± 7.11
		Phase 2	9.67 ± 7.94
		Phase 3	5.29 ± 7.02
		Phase 4	3.60 ± 4.87
		Phase 5	1.79 ± 2.48
		Phase 6	0.46 ± 0.68
		Phase 7	0.22 ± 0.47
		Phase 8	0.10 ± 0.26
		Phase 9	0.01 ± 0.04
		Phases 10 to 15	0
	G4DSR	Used counts / phase	6.45 ± 0.29
		Total	96.75 ± 4.36
	BH3D	Total	63.99 ± 4.58
	4.3.6 Quantitative comparison
	Here, G4D images are reconstructed with both 6 phases and 15 phases.
	Image noise, measured in a homogeneous liver area, for BH3D method (see
	gure 4.7 (d)) was signicantly reduced (-60.83% and -38.25%) in comparison
	to 15 phases-G4D and 6 phases-G4D (p < 0.001) and closer to NG3D mea-
	sured noise, although it was still signicantly higher (increase of +21.37%,
	p = 0.0012).
	The SUV peak and SUV mean for the BH3D method (see gure 4.7 (a)) show
	an increase of respectively +9.67% and +14.59% with respect to NG3D re-
	construction (p < 0.001), but remained signicantly lower (respectively -

NG3D, G4D and BH3D reconstruction methods were compared with respect to image-based and lesion-based quantitative measurements. The average values of the lesions' SUV peak , SUV mean and volumes, noise estimation in homogeneous liver area and signal-to-noise ratio (SNR) are presented in table 4.3.

5.29% and -11.99%) than those obtained by 15 phases-G4D reconstruction (p < 0.001).However, SUV peak and SUV mean were not signicantly dierent to the ones for 6 phases-G4D reconstruction (slight decrease of respectively reconstruction in micro-SPECT

Table 4 .

 4 3: Quantitative results for NG3D, G4D and BH3D reconstruction methods: SUV peak , SUV mean , lesion volume (threshold at 40% of the SUV max value), noise estimation in homogeneous liver area and SNR. Values are presented as mean ± SD. ± 40.02 19.81 ± 14.98 31.82 ± 22.50 45.61 ± 35.86

	Parameters	NG3D	G4D	a	BH3D
				(15 phases)	(6 phases)
	Noise	b	0.15 ± 0.04	0.47 ± 0.14	0.30 ± 0.08	0.18 ± 0.05
	SUV peak	7.26 ± 5.32	8.18 ± 5.56	7.87 ± 5.46	7.87 ± 5.47
	SUV mean	6.90 ± 5.26	8.63 ± 6.00	7.85 ± 5.67	7.80 ± 5.65
	Volume (mm 3 )	2.14 ± 1.29	1.60 ± 0.87	1.84 ± 1.00	1.79 ± 1.02
	SNR (SUV mean /Noise)	48.85

a G4D values are the average values of end-of-exhalation images correspond- ing to the motionless phases -i.e. nine or ten phases (respectively four)

Table 5 .

 5 1: Description of data used for the study of iodide uptake in the stomach.

	Tc-pertechnetate injection were
	99m

  3 in Appendix B. Let us adapt the equations to our current model. As mentioned in the previous section, the input blood activity function x 0 (t) is obtained from the measure of activity in the left ventricle. It is frequent to assume that the input function of the compartmental model are particular functions such as impulse or continuous function. In our study, this input function is neither an impulse nor a continuous function. However it can be approximated by a sum of two exponential functions, this will be shown in Section 6.3.

  2,2 t + α 2,3 e -β 2,3 t + α 2,4 e -β 2,4 t where each coecient α i and β i is a combination of k 01 , k 02 , k 21 , k 12 , k 10 , φ 1 , φ 2 , θ 1 , θ 2 . Conversely, the transfer parameters k 01 , k 10 , k 02 , k 21 , k 12 can be expressed and computed as combination of the

Table 6 .

 6

	1: Estimated compartment volume contribution (in % of layer vol-
	ume) of each compartment to each layer for subject #0530. These ratio were
	rst estimated from the observation of a chosen stomach IHC image. Then,
	the ratio for the rst 14 layers were manually optimised.
	Layer Nb W C B ND
	L1	15	0	1	84
		20	0	1	79
		22	0	1	77
		23	0	1	76
		24	0	1	75
	L6	25	0	1	74
		25	15	1	59
		20	25	1	54
		20	35	1	44
		15	45	1	39
	L11	12	55	1	32
		10	65	1	24
		10	75	1	14

Table 6 . 2 :

 62

Table A

 A

	.1: Discrepancy between respiratory and listmode triggering signals
	with two dierent minimisation criteria: Least square minimisation on the
	cycle length dierence, Minimisation of the number of outlier cycle (absolute
	dierence between cycle length has to be greater than 1 ms).
	Acquisition	Number of shifted cycles		Number of cycles
	number	at the beginning		with dierence
		of respiratory signal		between cycle length
		after signal registration		(respiratory signal
		after two minimisation		-listmode)
		methods:		equal to:
		LS	Outlier	-1	0	1	|di| > 1ms
		minimisation	minimisation
	182	370	370	569	870	0
	183	316	316	548	874	0
	186	6	6	593	936	0
	187	0	0	639	1138	0
	189	498	498	446	810	0
	190	57	57	485	970	0
	191	337	337	427	723	0
	1839	39	44	412	856	1
	1839 (truncated)	44	44 385 785 0	17
	1840	28	28	444	801	0
	1841	28	28	398	847	2
	1841 (truncated)	28	28 378 804 0	18
	1842	48	48	447	658	0
	that the two triggering signals are identical during the same temporal win-
	dow. Indeed, the respiratory signal record usually starts short before, and
	nishes short after the listmode, thus the beginning of the listmode triggering
	signal should be slightly shifted with respect to the respiratory signal one and
	entirely included into it.

Acknowledgements

We thank Audrey Lamit for the supply of biological material, the planning of experimentation, the manipulation of animals, and the support during acquisition.

Part II IODIDE DYNAMIC DISTRIBUTION IN STOMACH

A.1 Introduction

In this appendix, we present complementary observations, investigations, and results to the published paper entitled Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT presented in Chapter 4.

The overall aim of this appendix is twofold. reconstruction in micro-SPECT: Supplementary data First, we explain the reasons that led us to the approach presented in Chapter 4. We present preliminary observations that led us to deeply investigate the reasons of respiratory motion variations, and that initiated the proposed method for motionless three (spatial) dimensions (3D) single photon emission computed tomography (SPECT) reconstruction. Some variability of the respiratory rhythm (length and amplitude of pressure cycles) of animal are naturally observed during the acquisition procedures. The operator had to concentrate his/her attention in order to prevent high rhythmic variation by manually adjusting the rate of anaesthetic gas. However, this was not sucient to fully prevent variability of the respiratory rhythm. We have compared the average respiratory cycles extracted from the pressure sensor monitoring system for the set of acquisitions that we have and we tried to nd some common similarity and trend. These characteristics will be described in Section A.2.1. Then, we wanted to check that the abdominal lesions were animated with the same motion as the one describe by the average respiratory signal. In practice, from the preliminary comparison between the temporal motion of the lesions centre of mass (COM) and the respiratory cycle, we observed a temporal shift of the peak of motion. This comparison of observed motion and respiratory signal will be detailed in Section A.2.2.

Second, the message given by these two observations justify a deep study of two phenomena. In a rst time, we analysed the respiratory signal recorded by the monitoring system. In a second time, we studied the correlation between gating signals provided by the monitoring respiratory system and the one recorded by the camera (listmode). These investigations suggested to work exclusively with the gating signal (triggers) recorded in the listmode for postacquisition reconstructions. Section A.3 deals with all investigations done on the gating signals.

Finally, in Section A.4, we provide additional results showing the eciency of our approach. We show the performance of the successive improvement steps of our method for motionless 3D SPECT reconstruction. These improvements have been motivated by the investigations of the gating signal recorded in the listmode. In addition to the performance at each step, we show complementary results for the ultimate test of sensibility to the noise.

The results presented in the following were gathered in preliminary phases of our work and are intended to justify our approach. reconstruction in micro-SPECT: Supplementary data day of acquisition. Indeed, even if the same operator set up and controlled the animal for all the acquisitions, the placement of the animal in the bed and the positioning of the pressure sensor under the abdomen might be not exactly the same from one day to another. Moreover, the animal imaged on the same day were of the same lineage and age, and were injected with tumourous cells on the same day. We assume that the lesions were growing with the same speed and that the subject were aected by the lesions in a same manner.

Thus, we considered that all animal imaged on the same day were in a similar state. For all separate days of acquisition, there were no reproducibility with reconstruction in micro-SPECT: Supplementary data These observations suggest that the respiratory signal record might be failing. It cannot be used easily to give additional information about the cycles (amplitude variation associated with respiratory phase). As a consequence, it has been decided to give importance only to the information that is directly provided by the listmode in terms of cycle length. Though the listmode does not contain directly the respiratory signal, a respiratory motion information can be extracted intrinsically from the data, before or after image reconstruction. In our case, we chose to extract an intrinsic respiratory signal from the 4D image, after reconstruction.

A.4 Performance at each step of the method for motionless 3D SPECT reconstruction

In this section, we will comment in more details the benets of the proposed method with respect to the statistics at each step of its elaboration. This b Noise measure refers to the SD in a homogeneous liver area.

reconstruction in micro-SPECT: Supplementary data of lesion dened at 40% of the maximum SUV. As a consequence of the large increase in noise and smaller increase in SUV in G4D, the SNR also decreased by 60% with respect to NG3D image.

Second, G4D images after cycle selection (G4DS) have a slightly higher level of noise than in G4D images: noise increased by 2% in comparison to G4D and G4DS. This is comprehensible since around 2% of data were discarded.

The SUV measurements slightly increased comparing to G4D: SUV peak by 0.4% and SUV mean increased by 1%. The lesion volumes slightly decreased by 1.25%, with respect to G4D images.

Then, we perform a G4D reconstruction after cycle selection and resampling (G4DSR); we compare this G4DSR image to G4D and G4DS images.

The measurements of SUVs and lesion volume were improved: SUV peak and SUV mean were both increased by 2% with respect to G4DS image, and were increased by respectively 2.5% and 3% with respect to G4D image. The lesion volumes were decreased by around 3% with respect to G4DS image and by around 4% with respect to G4D image. These results are assessing the improvements of image accuracy as expected since the cycle resampling step of the listmode aimed at re-aecting the counts to their corresponding phase of the dynamic image. However, the noise increased again by around 2% comparing to G4DS image, that makes a 4% increase with respect to G4D and a 226% increase with respect to NG3D. The SNR also decreased by 0.5% comparing to G4DS, that makes a 1.5% decrease with respect to G4D and a 60% decrease with respect to NG3D.

As a result the quantitative accuracy of the images was improved by the successive modications of listmode (cycle selection and resampling) with respect to the lesions of interest. This suggests that each phase of the dynamic images were free of motion artifacts. However, the image quality was aected by the fact that less data were used for each phase image.

The modications on the cycles (selection resampling) slightly improved the accuracy. Additional work is necessary to improve the quality with respect to the dynamic images. The idea is then to distinguish the phases in motion from the motionless phases intrinsically to the dynamic images. According to the results with respect to the image accuracy,to extract the intrisic motion signal, we use the G4DSR images rather than the G4D or the G4DS images. Working on the G4DSR images increases the chance of determining the consecutive phases that best correspond to the motionless phase of the respiratory signal. Finally, the motionless images are simply a gated images reconstructed from data in the motionless phase according to the intrinsic motion signal. Motionless phase was determined during the study of intrinsic intensity variation along the time in region of interests (ROIs) of G4DSR images. It corresponds to two third of the respiratory cycle (around 10 phases over 15). These 10 phases were the one with the less variation in intensity with respect to time. Motionless images are also denoted as breath-hold like 3D

A.4. Performance at each step of the method for motionless 3D SPECT reconstruction 163 image (BH3D). In quantitative terms, the BH3D images have a better quality than G4DSR images with both a lower level of noise (63% decrease with respect to G4DSR) and a higher SNR (134% increase with respect to G4DSR).

BH3D images also have a better accuracy than NG3D image since BH3D images exhibit a 8% increase of the SUV peak , a 13% increase of the SUV mean , and a 16% decrease of the lesion volumes. Therefore, the improvement when compared to NG3D and G4D is mainly due to the motion correction with a minimal impact of cycle selection.

A.4.2 Sensibility of the image-based motion detection method with respect to noise

The test of sensibility to noise was described in Section 4.3.3 of Chapter 4. Visually the noisy lesion starts to be hard to distinguish with a noise equal to 422, that corresponds to a SNR around 2.4. with this method. It is not surprising since the lesion was barely visible on the images themselves and the method is based on the intensity variation of voxels along the time. A lower SNR means that the signal is drown in the noise and it makes it more dicult or impossible to distinguish the lesion from the background.

As a conclusion, the image-based motion detection method works on lesions with a minimal SNR of 3. Then, for images with SNR below, it will require additional constraints on the resulting histogram of votes to be able to automatically determine the motionless phases from the moving phases.

The method is then appropriate for images as the ones we used in Chapter 4

since the average SNR in the lesions is around 20, and greater than the limit of 3.

A.4.3 Sensibility of the detection method to the edge detection

In this section, we discuss the motion detection method performance with respect to the choice of detecting rising edge instead of falling edge on the pressure variation during the acquisition. Figure A.17 presents the results of the motion detection method of four dierent acquisitions. The dierence Appendix A. Amplitude-based data selection for optimal retrospective reconstruction in micro-SPECT: Supplementary data Appendix B

Compartmental modelling

The work presented in Chapter 6 relies on compartmental analysis. Compartmental analysis can be used for a wide range of applications. It is a frequent approach used in pharmacology, biology and medicine for studying dynamic biodistribution of compounds since it is quite simple to use and to understand.

In this section, we provide rst a description of the bases of compartment analysis (for additional information please refer to [START_REF] Godfrey | Compartmental models and their application[END_REF][START_REF] Jacquez | Compartmental analysis in biology and medicine[END_REF]. The idea is to dene all the elements required to build an appropriate and simple model for our problem.

B.1 Denitions

Some expression for compartmental analysis can be ambiguous, so we remind their denitions here. A compartmental system is a physiological or biochemical nite set of microscopic subunits located in an environment.

The substance is the studied material. Several substances can be studied simultaneously especially in pharmacology were substances react and yield to others.

The compartment is a macroscopic subunit of homogeneous and well-mixed substance. This means that in emission tomography (ET), the radioactivity concentration of the tracers is the same in the whole compartment (no gradient) or that the substance in this compartment is submitted to the same reaction or transformation function.

In ET, each subunit distinguishes from the others by dierent kinetic for the same substance. The interactions are characterised by the substance transfer or transport, or by substance transformation (between subunits or with the environment).

The space denotes the physical place where substances release or dilute.

A compartment does not always match with a specic closed space. Several compartments can be physically mixed. When a compartmental model stud-constants. In [START_REF] Foss | A method of exponential curve tting by numerical integration[END_REF], authors made these assumptions. They show that it is possible to determine the transfer coecients by approximating the discrete observations with sums of exponential functions.

The inverse problem can be solved analytically or numerically.

B.3 Example of a two-compartment system

For the sake of understanding the principles of compartmental analysis, we study below the case of a catenary two-compartment system. We assume that there is only one impulse input on compartment 1. The objective of this direct problem is to determine analytically x 1 (t) and x 2 (t), t ≥ 0 for the system represented in We detail below the dierent steps and equations of such a study in this simple case. Indeed, we will see in Section 6.2 that the compartment model of the stomach can be model by a two-compartment system, with To rephrase: similar transfer functions. However in the case of the stomach, the input of the system will be more complex.

Comp. 1 Here, we consider that (C) is the identity matrix. The system can be described by the following equations:

Laplace transformation gives: Appendix C

Pre-clinical image convention

Such as human radiological convention, there is a specic vertebrate radiological convention. This convention species anatomical directional terms and positioning convention. In general, the right and left sides are inverted that is the real left of the animal appears on the right of the image and the reverse.

The axis going from head to tale is named Craniocaudal axis, the axis going from left to right with respect to the spine is named Left-right axis, and the axis going from spinal column to abdomen is named Dorsoventral axis.