Julien Vibert

Véronique Thomas-Vaslin
email: veronique.thomas-vaslin@upmc.fr

Modelling T Cell Proliferation: Dynamics Heterogeneity Depending on Cell Differentiation, Age, and Genetic Background

Keywords: Bromo-deoxyuridine EdU: 5-ethynyl-2 cell proliferation, modelling immune system, aging, lymphocyte dynamics, cell differentiation, cell death, cell cycle, flow cytometry, mathematical model, state transition diagram

Cell proliferation is the common characteristic of all biological systems. The immune system insures the maintenance of body integrity on the basis of a continuous production of diversified T lymphocytes in the thymus. This involves processes of proliferation, differentiation, selection, death and migration of lymphocytes to peripheral tissues, where proliferation also occurs upon antigen recognition. Quantification of cell proliferation dynamics requires specific experimental methods and mathematical modelling. Here, we assess the impact of genetics and aging on the immune system by investigating the dynamics of proliferation of T lymphocytes across their differentiation through thymus and spleen in mice. Our investigation is based on single-cell multicolour flow cytometry analysis revealing the active incorporation of a thymidine analogue during S phase after pulse-chase-pulse experiments in vivo, versus cell DNA content. A generic mathematical model of state transition simulates through ODEs the evolution of single cell behaviour during various durations of labelling. It allows us to fit our data, to deduce proliferation rates and estimate cell cycle durations in sub-populations. Our model is simple and flexible and is validated with other durations of pulse/chase experiments.

Our results reveal that T cell proliferation is highly heterogeneous but with a specific "signature" that depends upon genetic origins, is specific to cell differentiation stages in thymus and spleen and is altered with age. In conclusion, our model allows us to infer proliferation rates and cell cycle phase durations from complex experimental 5-ethynyl-2'-deoxyuridine (EdU) data, revealing T cell proliferation heterogeneity and specific signatures.

Introduction

Cell division is a characteristic of biological systems, with variability of rates of proliferation and interphase duration according to the type of organism, organ, and period of life. In the immune system, cell proliferation is essential for the high turnover of diversified lymphocytes that insures dynamic cognitive function and maintenance of body integrity. The production of T cells occurs through processes of proliferation, differentiation, migration, selection and death in the thymus to finally export less than 6% of cells to peripheral tissues as naïve T cells [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF]. T cells circulate through the whole organism and interact with environmental and body antigens to proliferate and differentiate into effector/memory T cells upon antigenic recognition. Quantification and interpretation of lymphocyte dynamics remains a challenge in systems immunology [START_REF] Thomas-Vaslin | Lymphocyte Population Kinetics[END_REF][START_REF] Thomas-Vaslin | Quantifying Lymphocyte Division, Methods[END_REF]. Studies on T cell dynamics at various ages have questioned the mechanisms of thymic involution (4) beginning at puberty, followed by immuno-senescence and inflammation with aging [START_REF] Franceschi | Human immunosenescence: the prevailing of innate immunity, the failing of clonotypic immunity, and the filling of immunological space[END_REF], and physiological accumulation of effector/memory peripheral T cell clones, particularly of CD8 phenotype, in both aged humans [START_REF] Buchholz | CD8+ T cell differentiation in the aging immune system: until the last clone standing[END_REF] and mice [START_REF] Clambey | Non-malignant clonal expansions of CD8+ memory T cells in aged individuals[END_REF]. Although the number of T cells decreases with aging, promiscuous CD4 T cells accumulate [START_REF] Deshpande | Self-recognition drives the preferential accumulation of promiscuous CD4(+) T-cells in aged mice[END_REF].

Thymic hypoplasia or thymectomy in young individuals also affects the repertoire of peripheral T cells and immunocompetence [START_REF] Sauce | Altered thymic activity in early life: how does it affect the immune system in young adults? Current Opinion in Immunology[END_REF]. Thus, thymic production is primordial so as to maintain a homeostatic and dynamic equilibrium in terms of cell populations and repertoires, thereby allowing an effective and adaptive immune surveillance. Sequential T cell development and CD4/CD8 lineage decision models were previously described [START_REF] Penit | Sequential events in thymocyte differentiation and thymus regeneration revealed by a combination of bromodeoxyuridine DNA labeling and antimitotic drug treatment[END_REF][START_REF] Scollay | Thymic emigration: conveyor belts or lucky dips?[END_REF][START_REF] Germain | T-cell development and the CD4-CD8 lineage decision[END_REF][START_REF] Singer | Lineage fate and intense debate: myths, models and mechanisms of CD4-versus CD8-lineage choice[END_REF]. However, there remain open questions concerning the rates of thymocyte production that result from various processes, such as differentiation, proliferation, selection, and death of T cells, and that influence the homeostasis and life-span of naïve and effector/memory T cells in secondary lymphoid organs. Several theoretical approaches and experimental immunological protocols are available to investigate lymphocyte dynamics [START_REF] Thomas-Vaslin | Quantifying Lymphocyte Division, Methods[END_REF] and to model them [START_REF] Thomas-Vaslin | Lymphocytes dynamics and repertoires, modeling[END_REF]. A number of studies have tried [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF][START_REF] Von Boehmer | The life span of naive alpha/beta T cells in secondary lymphoid organs[END_REF][START_REF] Tough | Turnover of naive-and memory-phenotype T cells[END_REF], often with the help of modelling [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF][START_REF] Asquith | Lymphocyte kinetics: the interpretation of labelling data[END_REF][START_REF] Asquith | Lymphocyte kinetics in health and disease[END_REF][START_REF] Macallan | Measurement of proliferation and disappearance of rapid turnover cell populations in human studies using deuterium-labeled glucose[END_REF][START_REF] Boer | Quantifying T lymphocyte turnover[END_REF], to quantify thymocyte production, thymic output, and the processes involved in maintenance of a T lymphocyte dynamic equilibrium in the periphery in humans or mice [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF][START_REF] Callard | Modeling T-and B-cell growth and differentiation[END_REF][START_REF] Borghans | Quantification of T-cell dynamics: from telomeres to DNA labeling[END_REF][START_REF] Bains | Quantifying thymic export: combining models of naive T cell proliferation and TCR excision circle dynamics gives an explicit measure of thymic output[END_REF][START_REF] Bains | Quantifying the development of the peripheral naive CD4+ T-cell pool in humans[END_REF][START_REF] Bencheva | Computer Modeling of the Immune System Reconstruction after Peripheral Blood Stem Cell Transplantation[END_REF]. It was suggested that the naïve T cell pool is maintained by peripheral division in humans, at variance to thymic output in mice, and with significant differences in lifespan of CD8 peripheral T cells in young and old C57BL/6 mice [START_REF] Den Braber | Maintenance of Peripheral Naive T Cells Is Sustained by Thymus Output in Mice but Not Humans[END_REF]. Some findings also suggested different genetic regulations of thymic output according to strains of mice [START_REF] Dulude | The Magnitude of Thymic Output Is Genetically Determined through Controlled Intrathymic Precursor T Cell Proliferation[END_REF]. The flexibility of T cell behaviour has already been approached, in particular for peripheral cells [START_REF] Yamanaka | The dynamic lives of T cells: new approaches and themes[END_REF], but references concerning quantitative proliferation are still lacking.

Altogether, these studies were performed with various experimental methods, on various species, strains, ages, and cell populations. Consequently, no clear consensus has been reached concerning the influence of genetics and age on T cell dynamics. We have previously observed that both the genetic background and aging affect T cell dynamics and repertoires [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF]. Indeed, C57BL/6 and FVB mice differ in their T cell composition at steady state but also in their kinetics of T cell reconstitution after a transient immunosuppression. In FVB mice, there is an acceleration of thymic involution and of the immunological aging quantified by alterations in the homeostatic control of lymphocyte numbers and repertoire diversity in the periphery [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF].

To determine whether or not the differences in the dynamic behaviour of T cells previously observed between the FVB and B6 mice strains are related to the proliferation capacities of T cells during their differentiation, we assayed cell proliferation in vivo by active DNA labelling of cells in S phase. We designed an original in vivo cell labelling experiment with 5-ethynyl-2'-deoxyuridine (EdU) during a pulse-chase-pulse experiment, followed by multi-parameter flow cytometry analysis, on the basis of earlier methods that measured bivariate amounts of incorporated label versus DNA content in single cells [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF][START_REF] Dolbeare | Flow cytometric measurement of total DNA content and incorporated bromodeoxyuridine[END_REF]. To our knowledge, there is currently no model able to analyse the bi-dimensional dot plot flow cytometry data, in order to infer dynamic parameters of cell cycle phases after complex pulse-chase experiments. Thus, to interpret the complex cell dynamics during the experiment, we designed a generic mathematical model that describes the dynamics of cells that transit from G0/G1 to S and then to G2/M phase of the cell cycle according to variable periods of pulse-chase. Fitting of the model to the experimental bidimensional flow cytometry data allowed us to infer rates of proliferation, and to give an estimate of mean duration times of cell cycle phases. Our mathematical model is also able to fit data from other pulse/chase experiment protocols and by extension can potentially be used to investigate cell cycle kinetics in any cell type.

This study provides a detailed proliferation assessment and cell quantification in subpopulations of T cells, along their differentiation from thymus to spleen, and through aging. For the first time, we show that genetic origin and age drive a particular "signature" of proliferation and cycle phase durations, according to T cell differentiation stage and T cell lineage. FVB mice have significantly lower rates of T cell proliferation than C57BL/6 mice in both thymus and spleen. Proliferation decreases with age in both strains. These results are discussed in the perspective of T cell dynamics and proliferation of T cells, where aging and genetic peculiarities could be of importance.

Results

Assessment of in vivo cell proliferation by active DNA labelling

To determine proliferation properties of cell populations, active labelling of cells during S phase of the cell cycle was done with EdU, a thymidine analogue, during pulse/chase periods. In our specific case, two pulse periods were separated by a fourteenhour chase interval (Fig 1). Since the in vivo half-life of EdU is about one hour [START_REF] Cheraghali | Pharmacokinetics and bioavailability of 5-ethyl-2'-deoxyuridine and its novel (5R,6R)-5-bromo-6-ethoxy-5,6-dihydro prodrugs in mice[END_REF], this was equivalent to a pulse-chase-pulse experiment comprising an initial two-hour pulse followed by a fourteen-hour chase and then a second thirty-minute pulse. This protocol was chosen because it allows the labelling of cell populations with low proliferation rates, which would not be detected with only one EdU injection. Multicolour flow cytometry then allowed the identification of the cell cycle phases and the quantitative active EdU labelling status, using bi-dimensional EdU label/DNA content dot plots that represent single cell analysis. According to DNA content, cells can be classified in G0/G1 phase (DNA content equal to 1), G2/M phase (DNA content equal to 2), or S phase (DNA content between 1 and 2). According to the intensity of EdU label, one can distinguish between labelled cells that have progressed through S phase during the pulse, and nonlabelled cells. Thus, the dot plots represent the end snapshot of the cell population evolution through periods of pulse/chase. Three groups of cells can be delimited on our dot plots (Fig 1): G0/G1 unlabelled, G2/M unlabelled, and EdU labelled cells. Since the experiment ends with a pulse phase, there are no remaining unlabelled cells in S phase. Fig 1 . EdU pulse chase experiment for quantification of cell proliferation. Mice were injected intra-peritonealy with EdU at 0, 1, and 16 hours equivalent to a "pulse-chasepulse" experiment. Thirty minutes after the last pulse, the mice were sacrificed, lymphoid organs removed, and cell suspensions were analysed by single-cell multicolour flow cytometry. The bi-dimensional dot plot reveals the intensity of EdU label (log scale) and DNA content (linear scale arbitrary unit) in lymphocytes. Manual gating allows identification of unlabelled cells in either G0/G1 (G) or G2/M (M) that have not transited through S phase during pulse, and EdU labelled cells (G'+S'+M') that have transited through S phase.

Differential cell proliferation in T cells according to differentiation stage, age and strain

Our purpose was to assess the variation of T cell dynamics and heterogeneity in individuals, in time and location, by quantification of T cell proliferation in two different strains (FVB and B6), through aging (at 2 and 18 months) and during differentiation in lymphoid organs (from thymus to spleen). As exemplified in Fig 2A, in the thymus, the percentage of EdU labelled cells varies significantly according to strains of mice and the stage of thymocyte differentiation. Thus, comparing the FACS profiles between typical 2month-old FVB and B6 mice, one can observe that the proportion of EdU + cells is much lower in FVB mice, showing lower proliferation capacities in this strain. A given group of mice depicts a particular signature of proliferation according to strain and age (Fig 2B )    revealing the heterogeneity of proliferation across the stages of differentiation. The percentages of labelled cells after 16 hours are given in S2 Table . With aging, at 18 months, which corresponds to middle life in mice, the proliferation of thymocytes in B6 is reduced by about 1.5 fold in all populations. In 18-month-old FVB mice, the proliferative capacities are even more reduced in immature CD3 -CD4 -CD8 -TN cells and thereafter are quite comparable to the altered proliferation capacities already detected in young FVB mice. In conclusion, T-cell dynamics is very heterogeneous, depending on the differentiation stages as detected by the variation in the EdU uptake. Cell proliferation is altered with age, and genetic origin influences the proliferation capacities, which are depressed in FVB mice. Thus, a typical biological signature of thymocytes is revealed, allowing one to distinguish strains of mice and ages, although precise quantification of dynamic parameters is required.

Model of labelling experiment and fit to experimental data

Our aim was to quantitatively estimate proliferation parameters and mean life-span of cells (here restricted to inter-mitotic time), according to age and strain, in a standard way independent of pulse-chase experiment duration, and to express results as a proliferation rate/day. Thus we designed a mathematical model to interpret and fit the data obtained from the experimental cytometry bi-dimensional dot plots. To model system dynamics and the variation of "stocks" of cells as a function of time, equivalent representations can be proposed: the hydraulic metaphor, stock-flow/transition diagrams, integral equations, and differential equations [START_REF] Ossimitz | The basics of system dynamics: Discrete vs. continuous modelling of time[END_REF]. We chose to explain our model with a state-transition diagram (Fig 3), a visual language that is accessible to biologists and provides an equivalence to a system of ODEs [START_REF] Bersini | State-transition diagrams for biologists[END_REF][START_REF] Thomas-Vaslin | Dynamical and Mechanistic Reconstructive Approaches of T Lymphocyte Dynamics: Using Visual Modeling Languages to Bridge the Gap between Immunologists, Theoreticians, and Programmers[END_REF]. A state-transition diagram allows representation of the parallel processes occurring in cells and during the experiment.

Through time, the cells can progress across three stages, G0/G1, S, G2/M or exit the population (death, differentiation, migration). During the alternative phases of pulse and chase of EdU, the cells can exist in one of two stages, either EdU unlabelled (Edu -) or EdU labelled (EdU + ). Thus, by combination, the cells can progress in six different stages represented by six differential equations. Moreover, the cells progress through the cycle in the context of EdU bio-availability (pulse or chase). Altogether, the multiplicity of the states (3x2x2) results in a system of twelve equations (see Materials and Methods section and Tables 1 to 3 for the description of the mathematical model). Under our experimental conditions, the EdU labelling protocol is modelled by a first pulse phase (between t=0 and t=2 hours), a chase phase (between t=2 and t=16 hours) and a second pulse phase (between t=16 and t=16.5 hours). We defined an experimental result as a two-dimensional vector, representing the three percentages of G0/G1 unlabelled, G2/M unlabelled, and labelled cells which are equivalent to the G, M, and G'+S'+M' in the model. Fitting to experimental data was done by simulating the progression of cells through the cell cycle and labelling over a bi-dimensional range of parameters. To constrain parameters, we used the biological hypotheses detailed in Materials and Methods (Table 3).

Average proliferation rates and cell cycle durations in thymus and spleen from B6 and FVB mice

We first compared proliferation rates for whole thymus and whole spleen, abstracting the phenotypic heterogeneity of T cells, in young and old B6 and FVB mice (Fig 4 ,S1 Table). In the thymus, proliferation rates can vary up to three-fold, from 32% per day in young B6 mice to 9% per day in old FVB mice. Proliferation rates of thymocytes in young B6 mice are in general twice the value of those in young FVB mice (p<0.05). The proliferation rates of thymocytes from old B6 mice and young FVB mice are similar (around 17%/day). In both strains, there is a two-fold drop of proliferation rates in thymus with aging (p<0.05 in spleen and p<0.01 in thymus). In the spleen, in contrast to the thymus, the proliferation rates increase with age. From about 3% per day in both young strains, it doubles in old FVB (p<0.05) and quadruples in old B6 mice (p<0.01).

Fig. 4. Proliferation rates and cell cycle phase durations in thymus and spleen.

Proliferation rates (%/day), G0/G1 and G2/M phase durations (in days) in B6 and FVB mice (n=4 per group) according to age (2 and 18 months) in whole thymus (left) and whole spleen (right). Values are indicated in the box plot, the median is the 50th percentile, the lower and upper limits of the box are the 25th and 75th percentiles and the whiskers indicate 1.5x(interquartile range). *: p<0.05. **: p<0.01. (See Materials and Methods section for statistical tests used). Statistics are for proliferation rates only. Mean values and SD are given in S1 Table .     The duration of G0/G1 phase, S, and G2/M gives an indication of the mean time between two cell cycles (inter-mitotic time), which can be calculated as 1/proliferation rate. Again, these values are heterogeneous; the duration of G0/G1 phase increases with aging, G2/M duration varies from 2 hours in young B6 thymocytes to up to 63 hours in young FVB splenocytes.

Model validation: Simulation of another pulse-chase experiment from the literature

In order to consolidate our results and validate our mathematical model, we simulated a different pulse/chase experiment, as performed in 1990 by Baron and Penit. This in vivo experiment also assesses whole thymocyte kinetics from six-to eight-weekold B6 mice, using a thymidine analogue Bromo-deoxyuridine (BrdU) uptake, versus DNA content measurement [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF]. To simulate this experiment with our model, we set up a single injection of BrdU at the beginning of the experiment (1-hour pulse), followed by a nineteen-hour chase phase, during which within the BrdU-labelled cells the percentages of cells in G0/G1, S, and G2/M cells were quantified according to DNA content. We used cell cycle parameters retrieved from our own fitting procedure on the whole thymus of our 2-month old B6 mice (underlined in red in S1 Table ), and we simulated the experiment by setting the initial state to be 100% of cells in S phase. We also set α=0 and d G =a G . The comparison of the experimental results provided by Baron and Penit (15) and the simulation of their pulse/chase protocol with our parameters obtained from B6 data is shown in Fig 5 . The ability of our simulation to visibly approach Baron and Penit's experimental results confirms the plausibility of the biological parameters obtained in our young B6 mice and points to the fact that our model can be run with other experimental conditions of pulse/chase.

Fig. 5. Validation of the mathematical model and parameter values with pulse/chase experimental data from (15)

. Young B6 mice received a single BrdU pulse, and within the BrdU-labelled cells the percentages of cells in GO/G1, S or G2/M is recorded over a nineteen hour chase period. The points and dashed lines represent Baron experimental data. Simulation of the same pulse chase conditions using our mathematical model, with the mean parameter values (S1 Table) obtained from our 2-month-old B6 mice (n=4). The points and continuous lines represent the results of our simulation.

Signature of proliferation according to T cell differentiation stage, age, and genetics

The complexity and heterogeneity of cell populations that compose the thymus led us to analyse, more deeply, the proliferation rates of cells according to their differentiation stages, thus across time. Thymocytes were further decomposed in the earliest immature triple negative thymocyte CD4 -CD8 -CD3 -(TN), in DP (double positive thymocyte CD4+CD8+, sub-divided into DP CD3 lo and DP CD3 hi according to low or high expression of CD3), CD4, and CD8 mature CD3 + populations (S1 Figure upper panel). In the spleen, CD4 and CD8 T cells were separated as CD44 lo (naïve) and CD44 hi (antigen-experienced). CD4 + Foxp3 + (regulatory T cells, Treg) cells were also identified as a sub-population of CD4 T cells (S1 Fig lower panel). From all data collected in the thymus, we were able to establish a linear correlation between the observed percentages of EdU + cells and the fitted proliferation rates that is similar in the two strains of mice and is independent of the age (Fig 6). The equation y=-0.000884+0.0119x allows estimation of the proliferation rates, y, for further experiments using the same pulse/chase protocol duration, from x =%EdU + cells (r 2 =0.99, p<0.01). Estimated parameter values are given in S2 Table. [START_REF] Vasseur | Scheduled kinetics of cell proliferation and phenotypic changes during immature thymocyte generation[END_REF]. Hence in the thymus, while DN1 cells (the early thymic precursors) display minimal proliferation, in particular in old FVB mice (1.2%/day), CD8 mature thymocytes are the cells that display the highest proliferation rate (45%/day in young and a 2-fold decrease to 21%/day in old B6 mice) (S4 Table ). Thus, the intermitotic time for CD8 CD3 + mature thymocytes in young B6 are of about 2 days (S2 Table ), suggesting expansion of mature thymocytes before thymic export, while such proliferation is limited in FVB. The relative durations of G0/G1 and G2M/ phase and thus inter-mitotic time increase with aging in most thymic populations (S2 Figure). In the spleen aging perturbations are more irregular. Altogether these durations are again heterogeneous according to the stage of differentiation and T cell lineage.

Fig. 7. Quantification of dividing cells during thymocyte differentiation.

For each stage of thymocyte differentiation, the number of cells and the number of EdU+ cells (in millions) is quantified by FACS after pulse/chase. The correlation between these numbers is shown according to strains and ages of mice (upper panel). Box plots represent the quartiles and the lines represent the evolution of the median values through stage of differentiation thus across time, in four mice per group (lower panel). DN CD3 - cells are in yellow, DP cells are in red and decomposed into DPearly (CD4 hi CD8 hi ) and DP late (CD4 med CD8 med ). Then DP cells differentiate either into CD4 + CD3 + (blue) or into CD8 + CD3 + (green) mature T cells.

Numbers of cells per differentiation stage in the thymus are in accordance with the progression of thymocytes within 4 weeks as previously modelled in young FVB mice [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF]. These numbers are correlated to the numbers of EdU + cells (Fig 7 upper panel). Thus, an exponential growth of cells occurs during their residence in the thymus from DN2 where they actively divide up to late DP cells (Fig 7 lower panel). The minimal estimate of exponential cell divisions in B6 mice is 13 divisions in about 22 days, higher than the previous 9 divisions estimate [START_REF] Penit | Cell expansion and growth arrest phases during the transition from precursor (CD4-8-) to immature (CD4+8+) thymocytes in normal and genetically modified mice[END_REF]. This estimation does not take into account the cell death that occurs essentially at DN4 (1-3% of cells are dying) and at early DP stage for cells co-expressing the highest levels of CD4 and CD8 (death slightly increases from 13 to 15% in B6 mice and from 8 to 12% in FVB mice with ageing) during the thymic selection processes (Fig 7, S5 Table ). We calculated the ratio of cell number EdU + /death as a performance index of T cell expansion. In FVB, the EdU + /death ratio is lower than in B6 at all stages of differentiation (S5 Table, Fig 8). Thus, although the total number of cells reaches similar expansion, the balance between division and death appears quite different in both strains. The percentage of EdU + cells (red) and dead cells (blue) among each differentiation stage identified is given (lower panel). Quantification is done by FACS and given as the median (line) of 4 mice per group (points). In DN2 and DN3 cells in some B6 mice death is so rare that quantification is below 10 cells/thymus and sometimes up to 0 and cannot be represented on the log scale. The ratio of EdU + /dead cells gives a performance of cell expansion.

When cells migrate from the thymus to the spleen they reduce their proliferation. Naïve CD44 lo T cells reduce their proliferation rate to <1%/day in all groups (S1 Fig, S3 Table, S4 Table ). Thus a minimal inter-mitotic time of 226 days (a value similar to 218 days found in young B6 mice ( 27)), can be estimated for the more active CD8 naïve T cells in 2-month-old mice. In contrast, effector/memory T cells expressing CD44 hi after antigenic stimulation [START_REF] Budd | Distinction of virgin and memory T lymphocytes. Stable acquisition of the Pgp-1 glycoprotein concomitant with antigenic stimulation[END_REF], are actively labelled, suggesting inter-mitotic times of 12 to 16 days in CD4 and CD8 from 2-month-old B6 while it increases about two-fold (32 and 37 days) in FVB old mice. Regulatory T cells are the T cell population with maximal proliferation rates in young B6 (21%/d), dividing every 5 to 12 days, respectively, in 2and 18-month-old mice. In contrast, proliferation rates are as low as 1.4%/day suggesting that Treg cells divide only every 83 days in FVB mice, independently of the age.

Discussion

Investigation of cell proliferation and turnover in space and time currently remains a challenge in immunology, in order to quantify the parameters of cell population dynamics in vivo according to tissue/organ and cell population heterogeneity, as well as to identify peculiarities related to the state of individuals, such as age and genetic origin [START_REF] Thomas-Vaslin | Lymphocyte Population Kinetics[END_REF][START_REF] Thomas-Vaslin | Quantifying Lymphocyte Division, Methods[END_REF][START_REF] Thomas-Vaslin | Organism State, Lymphocyte[END_REF][START_REF] Thomas-Vaslin | Life Span, Turnover, Residence Time[END_REF][START_REF] Thomas-Vaslin | Modeling, Cell Division and Proliferation[END_REF].

Experimental approach and mathematical modelling

The model developed here allows estimation of proliferation rates and inter-mitotic time (with estimated parameters and confidence intervals for G0/G1 and G2/M cell cycle phase durations as described in Protocol S1 and Table S6) and allows quantitative comparison of the effect of age and genetics on lymphocyte dynamics from their early differentiation in thymus to their final maturation in spleen.

Heterogeneity of experimental approaches can be solved by modelling. Various

experimental approaches, such as active cell labelling during cell cycle with nucleotide analogues, deuterated water or glucose, passive DNA content measurement, depletion of dividing cells, assessment of cyclin activity, and measurement at the level of populations or single cells, were recently reviewed [START_REF] Thomas-Vaslin | Lymphocyte Population Kinetics[END_REF]. Mathematical or computer modelling is necessary to interpret complex data and obtain standardized values for parameters such as proliferation rates and cell cycle phase durations. This is the only way to compare results issuing from different experimental approaches (2) (P. Loap et al. in preparation). BrdU, or more recently EdU labelling pulse-chase experiments, with various durations and routes of administration were used to investigate lymphocyte turnover and life-span in mice. Only a few reports concerned bi-dimensional analysis of BrdU/DNA content [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF][START_REF] Rocha | Accumulation of bromodeoxyuridine-labeled cells in central and peripheral lymphoid organs: minimal estimates of production and turnover rates of mature lymphocytes[END_REF], and these investigations were limited to whole organ cell populations. Indeed, sensitivity of fluorochromes to reveal BrdU labelling, and the small number of parameter investigations limited studies using old cytometers. In the last few years, the introduction of Clickit EdU technology and multi-colour flow cytometry has allowed investigation of cell phenotype together with EdU and DNA content labelling, as set up in this paper.

Since complex experimental pulse chase periods change from one study to another, as in Baron's paper [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF], the only way to interpret the results and to infer quantitative parameter values of proliferation per day is to develop a mathematical model.

Our mathematical model, designed with visual language accessible to biologists, is easy to manipulate to quantify lymphocyte dynamics parameters and can be used to simulate all kinds of pulse/chase experimental designs for further investigations.

Moreover, it is flexible and can be easily modified to accommodate other phases of the cell cycle or different kinetics of label incorporation. This model was successfully applied to previous results from Penit group obtained with a different experimental pulse/chase protocol (cf Fig 5).

Cell heterogeneity and granularity of modelling. In the mathematical model and the process of fitting, it is assumed that the studied population has uniform cycling parameters. We make abstraction of certain important transition factors, such as differentiation, death, and migration between different compartments, as if the population was a closed system.

The estimated parameter values must be considered as averages over the whole population studied. Values for the total thymus are in fact an average over all subpopulations in the thymus, weighted by their relative importance in numbers of cells (DP cells are the most numerous intra-thymic population, and therefore results for the total thymus are similar to those for DP cells). A correct understanding of kinetic heterogeneity within the whole population therefore requires further identification into sub-populations so as to increase the resolution. However, precise rates of proliferation are difficult to calculate for heterogeneous populations during the differentiation process without integrating transition to cell death and transition between cell stages into a mathematic model. This is a general problem in the literature since to our knowledge there is currently no model that integrates such complex dynamical transitions. Indeed there is cell transition within populations (from rare DN, to DP to SP in thymus and from naïve to effector/memory in spleen), transition to cell death, and transition to another organ (from thymus to spleen) or recirculation. Here we can only quantify the numbers of dead cells at instant time (and not estimate rates/day) and death is not integrated in the current mathematical model. Inter-mitotic time and cell phase duration are only indicative but allow to compare the behaviour of various cell populations and conditions. While CD8 CD3 + mature thymocytes from young B6 mice have high proliferation rates (allowing for an inter-mitotic time of about 2.2 days) and possible expansion before thymic output, the CD8 CD3 + mature thymocytes from old FVB mice display a 100-fold lower proliferation rate and thus would be exported before having time to divide.

One should consider that between two pulses or within long pulse periods, cells may transit from differentiation stage or migrate to other tissues. However with a sixteenhour long pulse/chase/pulse experiment these transitions are rather limited. The observation of a linear relationship between the number of cells labelled by EdU and the proliferation rate (Fig 6) is probably explained by the fact that as the duration of the experiment is smaller on average than the duration of one full cell cycle, labelled cells have not had the time to transit to another compartment, so that the number of labelled cells is globally proportional to the rate of entry into S phase and the proliferation rate, as mice are in steady-state. Nonetheless if we were to model a longer experimental protocol, it would undoubtedly be necessary to take into account these processes. Moreover, the parameter values observed with our sixteen-hour pulse/chase/pulse are confirmed by simulating the context of a single pulse followed by a nineteen-hour chase, as in the Baron paper that assessed whole thymocyte proliferation.

Considering these limitations, we believe that our experimental design and modelling is able to give, if not exact absolute values, at least a correct understanding of the relative differences between the heterogeneous populations studied.

Overlap of phases with DNA content analysis. Since DNA content only allows one to distinguish three phases in the cell cycle, it is not possible to give separate estimates of durations for G2 and M phases or, more problematically, for G0 and G1 phases. Thus G0/G1 phase duration is an average between G0 and G1 durations, weighted by their relative importance in numbers of cells. New protocols and modelling introducing Ki-67 labelling to separate G0 from G1 cells are under investigation (P. Loap et al. in preparation).

Parameter constraint and fit. The relatively low information content of the experimental dot plots (three proportions, i.e., a two-dimensional vector) constrains our search for parameters. We notably restrain ourselves to the steady-state case (Hyp.4) since we believe that the in vivo conditions of the experiment in normal mice justify this hypothesis. We also fix loss of cells in S and G2/M to be null (Hyp.5), as confirmed by observations [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF]. Finally, we fix S phase duration to 6.5 hours (Hyp.7); this is a value that has been measured in the same type of cells, as simulated here [START_REF] Baron | Study of the thymocyte cell cycle by bivariate analysis of incorporated bromodeoxyuridine and DNA content[END_REF]. A stretched cell cycle duration model [START_REF] Witten | Modeling cellular aging: Variable cell cycle lengths[END_REF] was recently proposed for mature T and B cells [START_REF] Dowling | Stretched cell cycle model for proliferating lymphocytes[END_REF], but there is no indication this model could also occur in immature cells.

Heterogeneity of lymphocyte dynamics according to cell population, age and genetic origin

The seminal work of Penit has addressed the quantification of proliferation in the thymus of young B6 mice [START_REF] Vasseur | Scheduled kinetics of cell proliferation and phenotypic changes during immature thymocyte generation[END_REF]. However, to our knowledge we are the first to give proliferation rate/day at various granularities (from whole organ to rare populations) according to the time evolution of lymphocyte differentiation, from thymus to spleen, and in relation to the age and genetic origin of mice. FVB and B6 mice were chosen, because these strains belong to different phylogenetic clusters [START_REF] Witmer | The development of a highly informative mouse Simple Sequence Length Polymorphism (SSLP) marker set and construction of a mouse family tree using parsimony analysis[END_REF][START_REF] Goios | mtDNA phylogeny and evolution of laboratory mouse strains[END_REF] and display dissimilar T-Cell Receptor (TCR) repertoires as explained by chromosomal deletion of six Beta chain Variable region (BV) TCR families in FVB [START_REF] Osman | FVB/N (H2 q) mouse is resistant to arthritis induction and exhibits a genomic deletion of Tcell receptor V beta gene segments[END_REF]. B6 and FVB mice differ in their lymphocyte and T cell compositions [START_REF] Petkova | Genetic influence on immune phenotype revealed strain-specific variations in peripheral blood lineages[END_REF] and numbers in the various populations at steady-state [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF]. Moreover, FVB mice display accelerated thymic involution and aging of VB repertoire as compared to B6 mice [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF]. Thus, this suggests that variation in proliferation rates or cell death might explain these differences.

We have previously modelled thymocyte dynamics in young FVB mice and shown that 83% of the produced thymocytes die, with 35% of DP dying per day [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF]. Here, we quantify thymocyte cell death as apoptotic cells in sub G0/G1 and as cells with a DNA content >4, corresponding to macrophages engulfing dying cells. While in total thymus about 2-3% of cells are in death at an instant time, major cell death occurs between DN3 and DN4 stages if VDJ recombination fails (3-5% of instant death) and at the DP stage by negligence and negative selection processes, with up to 15% instant cell death at early DP stages before CD4 and CD8 modulation. During death by neglect [START_REF] Szondy | Thymocyte death by neglect: contribution of engulfing macrophages[END_REF] and negative selection macrophages remove dying auto-reactive thymocytes, before cell membrane permeabilisation and chromatin condensation [START_REF] Dzhagalov | Elimination of self-reactive T cells in the thymus: a timeline for negative selection[END_REF]. While the thymus of FVB young mice contains a higher number of cells than B6 mice (90 vs 81 millions of cells), we observe global defective proliferation rates as compared to B6 mice, in particular at early DN1 and DP stages but also for the proliferation of CD8 mature thymocytes (3%/day in FVB vs 40%/day in B6), and accumulation of cells at late DP stage. Conversely, in FVB mice, thymocyte death is generally higher. Thus, the death/EdU + ratio is higher at all stages except in DP. [START_REF] Giraud | Transient depletion of dividing T lymphocytes in mice induces the emergence of regulatory T cells and dominant tolerance to islet allografts[END_REF]. Our results show that genetics and age influence T cell dynamics, in particular increase in the time spent in G0/G1 in thymus with aging, related to the G0 elongation theory of aging. Aging thus increases the inter-mitotic time, and contributes to the decrease of the total number of thymocytes as the individual ages. Impaired thymocyte proliferation and higher cell death in FVB mice could be related to the immaturity of their thymic epithelial cells ( 52) that normally interact with thymocyte during differentiation, selection, and migration, as recently modelled in B6 mice [START_REF] Thomas-Vaslin | Dynamical and Mechanistic Reconstructive Approaches of T Lymphocyte Dynamics: Using Visual Modeling Languages to Bridge the Gap between Immunologists, Theoreticians, and Programmers[END_REF].

These defects in FVB mice may participate in the accelerated thymic involution and aging and the high CD4/CD8 ratio [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF], also observed here (CD4/CD8= 7 for FVB and 2 for B6). As modelled by Cohen's group in a reactive animation [START_REF] Efroni | Emergent dynamics of thymocyte development and lineage determination[END_REF], the cell competition for epithelial cells and decreased ratio of the dissociation rates of CD8 + and CD4 + thymocytes to 0.3, would induce an increase of the CD4/CD8 ratio to 7, thus possibly limiting the emergence of CD8 and their proliferation in FVB mice. As a consequence it could explain that while the recovery to steady state values following transient depletion of dividing cells is effective in young FVB mice [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF], this process is defective in old FVB mice, leading to CD8 clonal expansions in periphery in spleen [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF]. The defective Treg proliferation observed here in FVB mice with increased time passed in G0/G1 and G2/M even in young mice, may explain these TCR repertoire alterations in old mice. Indeed, in FVB mice, the transfer of Treg from young mice or the stimulation of their proliferation by low doses IL-2 treatment, initiated at 15 months of age, prevents repertoire alterations and CD8 clonal expansion at 2 years [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF].

The global proliferation rate of splenocytes increases with aging from 3 to 13% per day in B6 and from 3 to 7% per day in FVB. This 2-to 4-fold increase in proliferation with aging reflects, however, considerable heterogeneity according to naïve or effector/memory T cells. Our previous modelling using transient depletion of dividing cells (1) led to the conclusion that "recent thymic emigrants" in the spleen, which represent around 50% of all naïve T cells in young mice, make 1 to 2 cell divisions during the 9 days after thymic output, allowing to complete the thymic clonal expansion (of positively selected T cells) up to obtain at least 2 4 to 2 5 cells per naïve cell clone. The mean proliferation rate of effector/memory T cells is 14 times that of naïve T cells. This high level of proliferation in effector memory antigen-experimented T cells (inter-mitotic time varies from 12 to 19 days in both strains) is in accordance with our previous results

showing that immunological memory maintenance is dependent on active T cell division and that specific depletion of dividing T cells induces immunological amnesia while primary responses were intact [START_REF] Bellier | Turning immunological memory into amnesia by depletion of dividing T cells[END_REF]. However, with aging, the inter-mitotic time of CD44 hi cells increases. This suggests that the accumulation of CD44 hi cell number observed with aging is related to increased cell longevity rather than to active proliferation. This decrease in proliferation of effector memory T cells with aging could be related to hyper-glycosylation of T cell surface macromolecules, leading to altered T cell signalling [START_REF] Akha | Signal transduction in the aging immune system[END_REF]. Altered signalling pathways occur in T cells with aging, in particular as alteration of lipid raft polarisation, altering early steps of T cell activation, with increased activity of SHP-1 acting as a negative feedback on lymphocyte proliferation [START_REF] Goronzy | Goronzy signaling pathways aging T cells[END_REF][START_REF] Fulop | Cellular signaling in the aging immune system[END_REF]. This could explain the altered vaccine response in old individuals, and should be considered in systems biology approaches for vaccine design [START_REF] Six | Systems biology in vaccine design[END_REF].

The CD4 + Foxp3 + regulatory T cells, which represent less than 10% of CD4 splenic T cells, display the highest proliferation rate within splenic T cells. This is in accordance with their peculiar selection by medullary thymic epithelial cells [START_REF] Modigliani | Lymphocytes selected in allogeneic thymic epithelium mediate dominant tolerance toward tissue grafts of the thymic epithelium haplotype[END_REF] presenting promiscuous antigens because of Auto-Immune REgulator (AIRE) expression (60) of autoreactive Treg, then triggered in peripheral tissues to division by self-antigens cognition. We have previously shown that dominant tolerance to allogeneic grafts is related to antigen specific Treg cells selection in thymus [START_REF] Thomas-Vaslin | Thymic epithelium induces full tolerance to skin and heart but not to B lymphocyte grafts[END_REF]. In fact, a delicate balance between proliferation or suppressive function of Treg is of importance to insure the negative feedback loop, which controls the proliferation of other T cells. When Treg cells are under "suppressive mode", they are anergic, cannot proliferate [START_REF] Schwartz | T cell anergy[END_REF] and thus less sensitive to specific depletion of dividing T cells than other CD4 T cells. Thus, specifically targeting cells in proliferation can induce dominant tolerance [START_REF] Giraud | Transient depletion of dividing T lymphocytes in mice induces the emergence of regulatory T cells and dominant tolerance to islet allografts[END_REF].

Our general conclusion is that FVB mice have lower proliferation rates than B6 mice, and that with aging T lymphocytes have lower proliferation rates, although all mice are kept in the same environment. This certainly contributes to an accelerating aging. This reveals a specific "signature" of proliferation across differentiation stages in each of the two strains. This is in favour of the influence of genetic background and age on T cell dynamics and thus on the homeostatic equilibrium of the immune system. Therefore, the sensitivity of dividing cells -to depletion treatment by conditional immuno-pharmacogenetics that also reveal cell dynamics [START_REF] Thomas-Vaslin | Comprehensive assessment and mathematical modeling of T cell population dynamics and homeostasis[END_REF][START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF][START_REF] Giraud | Transient depletion of dividing T lymphocytes in mice induces the emergence of regulatory T cells and dominant tolerance to islet allografts[END_REF][START_REF] Bellier | Turning immunological memory into amnesia by depletion of dividing T cells[END_REF] or to gene therapy where retrovirus integration requires cell division (63) -could be very different with aging and genetic origin. Cell proliferation, age and genetics then influence the induction or control memory [START_REF] Bellier | Turning immunological memory into amnesia by depletion of dividing T cells[END_REF] and induction of tolerance [START_REF] Thomas-Vaslin | Immunodepression & Immunosuppression during aging[END_REF][START_REF] Giraud | Transient depletion of dividing T lymphocytes in mice induces the emergence of regulatory T cells and dominant tolerance to islet allografts[END_REF]. Thus, knowledge on lymphocyte dynamics and proliferation of sub-populations is of importance for clinical applications.

Materials and Methods

Mice

C57BL/6/N (B6) and FVB/N (FVB) mice were obtained from Charles River Laboratories maintained in SPF conditions and used at 2 and 18 months of age. Mice were manipulated according to European council directive 86/609/EEC of 24 November 1986 and with the approval of an ethics committee.

In vivo EdU treatment and multicolor flow cytometry analysis

The experimental protocol is summarized in Fig 1 . Mice received three intraperitoneal injections of 1 mg of EdU (5-ethynyl-2'-deoxyuridine, Life technologies) at 0, 1 and 16 hours and were killed half an hour after the last injection to remove thymus and spleen. Cell suspensions were obtained by mechanical disruption of organs in PBS + 3% newborn calf serum at 4°C, were then washed and submitted to CD16/CD32 blockage with 2.4G2 hybridoma cell supernatant for 10 minutes. Cell membrane labeling was then done for 20 minutes with antibodies coupled with fluorochromes Percp, APC-H7, and APC. After washing, the cells were fixed and permeabilized with the ClickIT fixative, and EdU was revealed by Click-IT containing Alexa488-azide. The last step consisted of staining with antibodies labeled with PE and PE-cy7 for 30 minutes and washing with Perm/Wash buffer (Becton-Dickinson). 1 µL FxCycle TM violet stain (Life Technologies) was added per tube, 30 minutes before acquisition. Cell parameters were acquired on LSR2 equipped with 405, 488, and 633 nm lasers (Becton-Dickinson) at a rate of 2000 cells/sec with DIVA software. Cells were analysed with Flowjo (www.flowjo.com) on the basis of structural cell markers FSC-A (Forward Scatter-Area), SSC-A (Side Scatter-Area), FSC-W (Forward Scatter-Width) to remove debris, apoptotic cells, and doublets.

Multicolour expression of specific markers CD3, CD4, CD8, CD25, CD44 allowed identification of cell sub-population phenotype by hierarchical gating. Then, DNA content and EdU label enabled construction of bi-dimensional dot plots to determine cell kinetics in each population. Dead cells were quantified in subG0/G1 apoptotic cells and cell that display a DNA content >4N, corresponding to macrophages engulfing dying cells before apoptosis.

Mathematical model

We designed an ODE (Ordinary Differential Equation) model describing cell dynamics during an EdU labelling experiment. Cells cycle through G0/G1, S, G2/M and are either unlabelled or labelled by EdU, giving six populations of cells represented by the six following ODEs (system of equations ( 1) ):

!" !" = 2𝑎 ! 𝑀 + 2𝛼𝑎 !! 𝑀 ! -𝑎 ! + 𝑑 ! 𝐺 !" !" = 𝑎 ! 𝐺 -𝑎 ! + 𝑑 ! 𝑆 -0, 𝛽 𝑆 !" !" = 𝑎 ! 𝑆 -𝑎 ! + 𝑑 ! 𝑀 !"! !" = 2(1 -𝛼)𝑎 ! ! 𝑀′ -𝑎 ! ! + 𝑑 ! ! 𝐺′ !"! !" = 𝑎 !! 𝐺 ! -𝑎 ! ! + 𝑑 ! ! 𝑆 ! + 0, 𝛽 𝑆 !"! !" = 𝑎 !! 𝑆′ -𝑎 ! ! + 𝑑 ! ! 𝑀′ (1) (1)
Populations are described in Table 1, parameters in Table 2, the equivalent state transition

diagram is in Fig 2.

Table 1. Cell populations described in the model and initial values used in the fit according to our hypotheses.

Population name

Population description Initial value used for the fit

G

Cells in G0 or G1 phase, unlabelled (G 0 ,S 0 ,M 0 ) and (G' 0 ,S' 0 ,M' 0 ) ; and the labelling protocol, i.e., the beginning and end of pulse and chase phases.

(𝐺 + 𝐺′) !! = !""! ! ! ! !! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! S Cells in S phase, unlabelled (𝑆 + 𝑆′) !! = !""! ! ! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! M Cells in G2 or M phase, unlabelled (𝑀 + 𝑀′) !! = !"" ! ! !! ! ! ! !! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! G' Cells in G0

Fit of simulations to experimental results

The results of simulations were fitted to experimental results represented by EdU/DNA flow cytometry dot plots. On the dot plot, three groups of cells can be individualized and their relative percentages determined according to EdU labelling and DNA content: G0/G1 EdU -, G2/M EdU -, and EdU + . We therefore defined an experimental result by three percentages: (G exp ,M exp ,(G'+S'+M') exp ). To fit our model to experimental results and to restrict the search for parameters, we made hypotheses listed in Table 3:

- -Hyp.3: There is no shedding of EdU label during the time of the experiment, i.e., α=0.

-Hyp.4: Cells are in a steady-state, and so is the distribution of cells in each cycle phase ((G+G') SS ,(S+S') SS ,(M+M') SS ) as described by system of equations (2):

!(!!! ! ) !" = 0 !(!!! ! ) !" = 0 !(!!! ! ) !" = 0
(𝐺 + 𝐺′) !! = !""! ! ! ! !! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! (𝑆 + 𝑆′) !! = !""! ! ! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! (𝑀 + 𝑀′) !! = !"" ! ! !! ! ! ! !! ! !! ! ! ! !! ! !!! ! ! ! ! ! ! !! ! ! ! !! ! (2)
with

(𝐺 + 𝐺′) !! + (𝑆 + 𝑆′) !! + (𝑀 + 𝑀′) !! = 100(%).
As a consequence, the total number of cells during the experiment is constant and percentages can be used. This also gives a conservation equation as follows:

a M M+a M' M'=d G G+ d S S+ d M M+ d G' G'+ d S' S'+ d M' M'.
-Hyp.5: There is no exit during S or G2/M phase (15), i.e., d S =d M =d S' =d M' =0.

-Hyp.6: One injection of EdU is modelled by a one-hour pulse [START_REF] Cheraghali | Pharmacokinetics and bioavailability of 5-ethyl-2'-deoxyuridine and its novel (5R,6R)-5-bromo-6-ethoxy-5,6-dihydro prodrugs in mice[END_REF].

-Hyp.7: S phase lasts 6.5 hours (15), i.e., a S =1/6.5.

G 0 =(G+G') SS ,S 0 =(S+S') SS ,M 0 =(M+M') SS ,G' 0 =0,S' 0 =0,M' 0 =0) .
To fit our model to experimental results, we run simulations over a twodimensional range of parameters (a G ,a M ). The time step chosen with the integration tool was 0.01 hour (36 seconds). The best fit is the couple of parameters which minimizes the Euclidean distance between (G exp ,M exp ,(G'+S'+M') exp ) and (G sim ,M sim ,(G'+S'+M') sim ).

All results obtained have a Euclidean distance to experiment of less than 0.05. We deduce from the best fit the mean duration of G0/G1 phase 1/a G , the mean duration of G2/M phase 1/a M , and the proliferation rate p, which with our hypotheses is equal to 200a M a G /(2a M +13a M a G +2a G ) (in percentage of cells dividing per hour).

The identifiability of parameters obtained during the fitting procedure and some illustrative results are given in S1 protocol and S6 Table . We also give standard deviations for our estimates of proliferation rates in total thymus and spleen, calculated with the Hessian matrix (variance-covariance matrix method).

Statistics

Statistics were performed using the R software (www.r-project.org). Statistical significance level was fixed at α=0.05 (Type 1 error). Before applying comparison tests, all data were assessed for normality by Shapiro-Wilk's test and for equality of variances by Levene's test. These results led us to the following choice of statistical tests: (i) if the normality hypothesis did not hold, we used the non-parametric Mann-Whitney test (ii) otherwise, we used the Student's two-tailed unpaired T-test with a hypothesis of either equal variances or unequal variances, depending on the result of Levene's test. We computed arithmetical means and standard deviations (variance square root) for each group (n=4).

Software

The mathematical model of cell proliferation and its fitting to experimental data were developed and simulated using Python (http://www.python.org). To allow highthroughput analyses, we developed a program in Python to implement the fitting procedure on tables of values extracted from the dot plot manual gates with FlowJo software or from published results recovered from 2-dimentional descriptions of cell kinetics. Graphs were produced with JMPpro11.

S5 Table. Evolution of thymocyte numbers according to strains, ages, and differentiation stages.

Numbers are given as millions of cells in the thymus (see Fig 7 and8). The ratio between EdU + and dead cells gives a performance of cell expansion. DP total represents the sum of the DPe (gated on CD4 hi CD8 hi ) and DP late (gated on CD4 med CD8 med ). . Proliferation rates with confidence intervals and standard deviations of single mice, for total thymus and total spleen. Proliferation rates (%/day) with confidence intervals (CI= IC_min -IC_max) and standard deviations (sd) calculated with the use of the Hessian matrix of all sixteen mice in whole thymus and whole spleen. Confidence intervals and standard deviations are calculated as described in the S1 Protocol.

S6 Table

S1 Protocol: Identifiability of parameters and calculation of confidence intervals and standard deviations in individual mice allowing for mathematical proliferation rate and cell cycle phase duration estimation

To illustrate the fitting procedure which allows estimating proliferation rates and cycle phase durations from experimental data, we choose the free parameters to be (S+S') SS and (M+M') SS , which we will call S SS and M SS for convenience, i.e. proportions of cells in S phase and G2/M phase at steady-state (Hyp.4), since S SS is directly proportional to proliferation rate. In this case, proliferation rate is equal to p=a M M SS =a S S SS =S SS /6.5 (Hyp.7). This choice of parameters is arbitrary and has no incidence on results, since all results are fixed by two free parameters. Mean duration of G0/G1 phase is 1/a G with a G =a S S SS /G SS =a S S SS /(100-S SS -M SS ), and mean duration of G2/M phase is 1/a M with a M =a S S SS /M SS . After having run 160,000 simulations (S SS and M SS each ranging from 0.1% to 40% every 0.1%), the best fit is the pair of parameters, which minimizes the Euclidean distance between experiment and simulation.

To take into account the variability of parameters giving a good fit, we define "confidence intervals" for our fitting estimates of proliferation rates, corresponding to the variability of proliferation rates with a distance to experiment of less than 1. This limit of 1 is based on the fact that our flow cytometry data (and also our fitting grid) is precise to the hundredth of a percent (0.01%), and we therefore judge it appropriate to define 1 as a limit distance between experiment and simulation for calculating our confidence intervals. It is however of course arbitrary, and the size of confidence intervals are indeed directly related to the value of this limit.

An illustrative example of the confidence interval is show below. Confidence intervals for proliferation rates of all sixteen mice in whole thymus and whole spleen are shown in

S6 Table.
Identifiability of parameters obtained during the fitting procedure, with a graphical representation. The example is for total thymocytes of a B6 18-month old mouse. After having run 160,000 simulations (S SS and M SS each ranging from 0.1% to 40% every 0.1%), the best fit is the pair of parameters which minimizes the Euclidean distance between experiment and simulation. Euclidean distance (ms) between experiment and simulation is plotted here. Right panel (B) corresponds to zooming in on the red square on the left (A). Point B is the best fit (Euclidean distance = 0.11), it corresponds to the pair S SS =5.7% and M SS =8.6%, i.e. proliferation rate is 21.0 % (/day), G0/G1 phase is 97.7 hours and G2/M phase is 9.8 hours. (In our fitting procedure, we use a resolution of 0.01% so that the best fit is actually S SS =5.74% and M SS =8.55% with Euclidean distance of 0.01, but for graphical purposes we do not use this resolution in the example). In order to quantify the variability of the fit according to parameters, we take all fits with a maximum Euclidean distance with experience of 1. These are contained inside a zone centered by point B and bordered by points A and C. In this case, S SS varies between 5.4% and 6.1% (red double arrow), i.e. between proliferation rates of 19.9% and 22.5%.

On the other hand, M SS varies between 7.3 and 9.8% (yellow double arrow). Thus G0/G1 phase varies between 93.9 and 100.1 hours, and G2/M phase between 7.9 and 11.6 hours, leading to an inter-mitotic time varying from 4.5 to 4.9 days in this mouse. "Confidence intervals" for our fitting estimates of proliferation rates correspond to the variability of proliferation rates with a distance to experiment of less than 1 (red double arrow).

Fig. 2 .

 2 Fig. 2. Thymocyte proliferation in mice according to cell differentiation stage, age and genetic origin. Cell proliferation was assessed in C57BL/6 and FVB mice of 2 or 18 months as in Fig 1. (A) Multi-colour flow cytometry and hierarchical gating in

Fig. 3 .

 3 Fig. 3. State transition model of cell cycle and EdU labeling. According to the FACS dot plot, the evolution of cells through the cell cycle and EdU labelling can be represented by a transition diagram. Unlabelled cells can transit successively in the three phases of the cell cycle: G0/G1 (G), S (S), and G2/M (M). Terms in a X correspond to rates of entry into the next phase of cells in X phase (X being either G, S, or M). The exit of cells, either due to death, differentiation, or migration is represented by terms in d X . During the pulse-phase the cells in S phase incorporate EdU (red arrow) and enter into S' with rate β and continue the cycle to M' and G'. During the chase phase, in the absence of EdU, the labelled cells lose labelling upon several cell divisions (blue arrow) (delabelling) with rate α. Unlabelled cells can enter into S phase but remain unlabelled.

Fig. 6 .

 6 Fig. 6. Signatures of lymphocyte proliferation according to strains and ages. The upper panel shows the correlation between the proliferation rate and percentages of EdU + cells observed 16 h after the pulse/chase/pulse experiment in thymus of B6 and FVB mice (80 values) and according to the age of mice (2 months: blue circle, 18 months: red circles). Hierarchical unsupervised clustering showing the clusters of mice according to the rates of proliferation. Values and statistics are given in Tables S1, S2, S3, S4.

Fig. 8 .

 8 Fig. 8. Quantification of cell proliferation and death in thymus. For each stage of thymocyte differentiation, the total number of live cells (green), of accumulated EdU + (red) cells, and of dead cells (blue) after the pulse/chase/pulse are shown (upper panel).The percentage of EdU + cells (red) and dead cells (blue) among each differentiation stage identified is given (lower panel). Quantification is done by FACS and given as the median (line) of 4 mice per group (points). In DN2 and DN3 cells in some B6 mice death is so rare that quantification is below 10 cells/thymus and sometimes up to 0 and cannot be represented on the log scale. The ratio of EdU + /dead cells gives a performance of cell expansion.

β

  a G' Rate of entry into S phase of G1 phase unlabelled, labelled cells a G = a G ' = parameter to fit a S , a S' Rate of entry into G2 phase of S phase unlabelled, labelled cells a S = a S ' = 1/6.5 a M , a M' Rate of entry into cell division of M phase unlabelled, labelled cells a M = a M ' = parameter to fit d G , d G' Rate of cell death, differentiation or migration of G0/G1 phase unlabelled, labelled cells d G = d G' = a G d S , d S' Rate of cell death, differentiation or migration of S phase unlabelled, labelled cells d S = d S' = 0 d M , d M' Rate of cell death, differentiation or migration of G2/M phase unlabelled, labelled cells d M = d M' = 0 α Proportion (0≤α≤1) of labelled M phase cells shedding label during cell division α = 0 (no cell loses its label with cell division) Rate of labelling of S phase cells during the pulse phase β = ∞ (instant labelling by EdU) p Proliferation rate Result of the fit G, S, M are numbers of cells in G0/G1, S, G2/M phase of the cell cycle, respectively. G', S', M' are the corresponding EdU labelled cell numbers. a X terms correspond to rates of entry into the next phase of cells in X phase (X being either G, S, M or G', S', M'). The exit of cells, either due to death, differentiation, or migration is represented by d X terms. EdU labelling intensity is modelled by two terms: β is the rate of labelling of S phase cells by EdU, while α is the rate of loss of labelling by labelled cells after division. To distinguish pulse and chase phases of the experiment, we have introduced the symbol [0,β] , meaning that during the chase phase (no EdU) we use 0, while during the pulse phase we use β in the equation (equivalent to a Dirac function multiplied by β , with the Dirac function equal to 1 during pulse and 0 during chase). The proliferation rate is defined as p = a M M+a M' M' , i.e., it is the number of cells in G2/M phase, unlabelled (M) and labelled (M'), which divide at rate a M (unlabelled) or a M' (labelled). A simulation is run by assigning: the twelve parameters of the cell cycle (a G ,a S ,a M ,d G ,d S ,d M ) and (a G' ,a S' ,a M' ,d G' ,d S' ,d M' ); the two parameters corresponding to delabelling and labelling, respectively (α,β) ; the initial state of the cell populations

Hyp. 1 :

 1 Parameters during pulse and chase are equal, reducing to six the number of parameters of the cell cycle, i.e., (a G ,a S ,a M ,d G ,d S ,d M )=(a G' ,a S' ,a M' ,d G' ,d S' ,d M').-Hyp.2: EdU cell labelling is instantaneous, i.e., β=∞. While mathematically not integrable as such, this is very simply implemented in the computer model by transferring at the start of the pulse phase all unlabelled S cells to the labelled S compartment, and by allowing the unmarked G1 cells to only go into the labelled S compartment during pulse phase.

S2Fig.

  Heterogeneity of duration of cycle phase and intermitotic time in thymic and splenic sub-populations, and variability with strains of mice and aging. The G0/G1 and G2/M estimated durations are in hours. The inter-mitotic time durations are given in days. The lines represent median values of 4 mice per group.S1 Table. Parameter values obtained by fitting the mathematical model to whole spleen and whole thymus.Mean and standard deviation with n=4 mice per groupS4 Table. Comparison of proliferation rates according to strain, age, organ, and T cell populations.Total represents the whole organ. CD4 and CD8 mature T cells are observed in the thymus and spleen. Mean proliferation rates per day are given during the transition of cells from DN1 to DP_CD3 + then CD4 or CD8 in thymus and from naïve (CD44 lo ) to effector/memory (CD44 hi ) differentiation in spleen; Foxp3 cells are a subpopulation of CD4 that are CD44 hi . These values are only indicative, since the transition of cells from one stage to another and to death are not modelled.

  Proliferation rates/day are significantly different (p<0.05) comparing :Fig 6). The heterogeneity of cell dynamics is revealed at the various steps of differentiation. If decomposition of the population is performed so as to identify rare populations as DN1 to DN4, then to DP that increase CD3 expression through maturation while modulating CD4 and CD8 co-expression, then one can observe oscillations in proliferation (S4 Table and S1 Fig median panel), as previously observed in young B6 mice using other protocols

	In thymus : B6 2M vs B6 18M: 2M>18M in DP, CD8 ; FVB 2M vs FVB 18M:
	2M>18M in DN and CD4 ; B6 2M vs FVB 2M:B6>FVB in DP CD3 + , CD8 ; B6 18M vs
	FVB 18M:B6>FVB in CD8
	In spleen : B6 2M vs B6 18M: 2M>18M in CD4 CD44 hi ; FVB 2M vs FVB 18M:
	2M>18M in CD4 CD44 hi ;B6 2M vs FVB 2M:B6>FVB in CD4, CD4 Foxp3 + ; B6 18M
	vs FVB 18M:B6>FVB in CD4 Foxp3 + , CD8 CD44 hi
	Because the intragroup variability of proliferation rates is low while inter-group
	heterogeneity is high (S1 Fig), unsupervised hierarchical clustering allows identification
	of clusters of mice according to age, strain, and differentiation stage giving particular
	signatures (

Table 3 . Biological and experimental hypotheses, allowing for constraint of parameter values.

 3

	Hypothesis	Hypothesis description	Consequence for parameters
	number
	1	EdU labelling does not affect cell-	(a G ,a S ,a M ,d G ,d S ,d M )=(a G' ,a S' ,a M' ,d G' ,d S' ,d M' )
		cycle kinetics
	2	Labelling of S phase cells is	β=∞
		instantaneous during the pulse
		phase
	3	Labelled cells do not divide	α=0
		enough times during the
		experiment to shed label
	4	No dynamics perturbation under	(G 0 =(G+G') SS ,S 0 =(S+S') SS ,M 0 =(M+M') SS
		physiological conditions: steady-	)
		state cell number in each phase of cell cycle and in each cell population	a M M+a M' M'=d G G+ d S S+ d M M+ d G' G'+ d S' S'+ d M' M'
	5	No cell death or differentiation	d S =d M =d S' =d M' =0
		during S phase and G2/M phase
	6	Bio-disponibility of EdU is one	None
		hour after an injection
	7	S phase lasts 6.5 hours		a S =1/6.5
	With these hypotheses, ten out of the twelve parameters of the cell cycle are fixed,
	except for a M and a G , since:
	(a G ,a S =1/6.5,a M ,d G =a G ,d S =0,d M =0)=(a G' ,a S' ,a M' ,d G' ,d S' ,d M' ). The two parameters of EdU
	labelling	intensity	are:	(α=0,β=∞).	The	initial	state	is:
	(

Sample Age (months) Proliferation (%/day) IC_min IC_max sd Proliferation (%/day) IC_min IC_max sd

				Thymus			Spleen
	B6_43	2	29.4	28.1	30.6 0.18	3.6	1.8	5.2 0.08
	B6_44	2	36.3	35.1	37.7 0.14	3.1	1.1	4.8 0.19
	B6_45	2	33.2	31.8	34.7 0.09	3	1.1	4.4 0.05
	B6_46	2	31	29.5	32.5 0.07	3.5	1.8	5.2 0.18
	B6_09	18	21.2	19.9	22.5 0.16	10.8	9.6	11.8 0.08
	B6_10	18	14.1	12.9	15.1 0.06	20.6	19.2	21.8 0.12
	B6_11	18	7.8	6.3	8.9 0.05	12.8	11.8	14 0.09
	B6_12	18	22	20.7	23.3 0.12	9.8	8.5	11.1 0.19
	FVB_31	2	18.1	16.6	19.2 0.06	3.5	1.8	5.2 0.18
	FVB_32	2	18.9	17.7	20.3 0.11	2.8	0.7	4.8 0.10
	FVB_33	2	14.2	13.3	15.5 0.19	2.3	0.4	4.1 0.05
	FVB_34	2	18.9	17.7	20.3 0.12	3.7	1.8	5.2 0.07
	FVB_14	18	10.6	9.6	11.8 0.13	5.4	3.7	7	0.16
	FVB_15	18	8.8	7.4	10 0.12	5.1	3.3	6.6 0.04
	FVB_16	18	7.3	5.9	8.5 0.06	7.2	5.5	8.9 0.16
	FVB_17	18	11.1	10	12.2 0.02	11.3	10	12.6 0.14

S3 Table. Splenocyte dynamics parameter values in B6 and FVB mice aged 2 and 18 months.

Estimated percentage of proliferation/day from the model, allowing estimation of duration of G0/G1, G2/M phase, and the potential inter-mitotic time for various cell populations. Whole CD4-CD3 hi and CD8-CD3 hi cells are decomposed into CD44 hi (effector/memory) and CD44 lo (naïve) cells showing the heterogeneity of dynamics according to the granularity of populations. CD4-Foxp3 are regulatory T cells Foxp3 hi . These values are only indicative, since transition of cells from one stage to another and to death are not modelled. Statistical analysis obtained by fitting the procedure for the proliferation rate (%/day), in spleen. Statistics f are given for the populations of total spleen; CD4, CD4 CD44 lo , CD4 CD44 hi , CD4 FoxP3+; CD8, CD8 CD44 lo and CD8 CD44 hi splenocytes. B6_2M: 2 month-old B6 mice, B6_18M: 18 month-old B6 mice, FVB_2M: 2 month-old FVB mice, FVB_18M: 18 month-old FVB mice. Level of significance of statistical tests: > (resp. <) indicates that in population p, group a (on the left) has a mean which is superior (resp. inferior) to group b (on the right) with a level of significance of p<0.05; >> (resp. <<) is for p<0.01; >>> (resp. <<<) is for p<0.001.

Calculation of standard deviations for estimates of proliferation rates.

To compute standard deviations for our estimated parameters, we use the Hessian matrix defined as ∇ ! 𝑓(𝑥 * ) where 𝑥 * is the point in parameter space which corresponds to the best fit, i.e. (𝑆 !! * , 𝑀 !! * ).

We compute the Hessian matrix in the following way:

where we use the Jacobian matrix :

with 𝑟 ! (𝑥 * ) = 𝐺 !"# (𝑥 * ) -𝐺 !"# and 𝑟 ! 𝑥 * = 𝑀 !"# 𝑥 * -𝑀 !"# (and we neglect the term 𝑄(𝑥 * ) with second-order derivatives which approximates to zero close to 𝑥 * ).

We compute numerically the Jacobian in the following way:


Véronique Maume-Deschamps

Didier Rullière

Antoine Usseglio-Carleve

Expectile prediction through asymmetric kriging

come

Introduction

Kriging (see [START_REF] Krige | A statistical approach to some basic mine valuation problems on the witwatersrand[END_REF]) aims at predicting the conditional mean of a random field (Z t ) t∈T given the values Z t 1 , ..., Z t N of the field at some points t 1 , ..., t N ∈ T , where typically T ⊂ R d . It seems natural to predict, in the same spirit as Kriging, other functionals. In our study, we focus on expectiles for elliptical random fields. We also did a similar work for quantiles (see [START_REF] Maume-Deschamps | Spatial Quantile Predictions for Elliptical Random Fields[END_REF]). [START_REF] Cambanis | On the theory of elliptically contoured distributions[END_REF] give the representation : the random vector X ∈ R d is elliptical with parameters µ ∈ R d and Σ ∈ R d×d , if and only if

Elliptical Distributions

X = µ + RΛU (d) , (1)
where ΛΛ T = Σ, U (d) is a d-dimensional random vector uniformly distributed on S d-1 (the unit disk of dimension d), and R is a non-negative random variable independent of U (d) . Furthermore, X is said consistent if :

R d = χ d (2) Distribution Gaussian 1 Student, ν > 0 χ ν √ ν Unimodal Gaussian Mixture n k=1 π k δ θ k Laplace, λ > 0 1 √ E( 1 λ ) Slash β (a, 1)
Table 1: Some consistent distributions

We can now define the notion of elliptical random fields. Indeed, a random field (X(t)) t∈T is R-elliptical if ∀n ∈ N, ∀t 1 , ..., t n ∈ T , the vector (X(t 1 ), ..., X(t n )) is (R, n)-elliptical.

, if X = (X 1 , X 2 ) is a (R, d 1 + d 2 )-elliptical random vector with parameters µ = (µ 1 , µ 2 ), we know that X 2 |X 1 is (R * , d 2 )
-elliptical with parameters :

µ 2|1 = µ 2 + Σ 21 Σ -1 11 (x 1 -µ 1 ) Σ 2|1 = Σ 22 -Σ 21 Σ -1 11 Σ 12 (3)

Elliptical expectiles

If X is a random variable, the α-expectile of X is given by :

e α (X) = arg min x∈R E[φ α (X -x)] where the scoring function φ α is (1 -α)x 2 1 {x<0} + αx 2 1 {x>0} = αx 2 + (1 -2α)x 2 1 {x<0} .
(4) Let X be a (R, 1)-elliptical random variable, with parameters µ and σ 2 . We define the function Ψ R as follows :

Ψ R (x) = x -∞ c 1 g 1 (y 2 )dy + 1 x +∞ x yc 1 g 1 (y 2 )dy where c 1 g 1 (t) = 1 2 f R ( √ t), f R (t)
being the density of R. Then, the α-expectile of X is given by :

e α (X) = µ + σΨ -1 R       α 2α -1       (5)
We propose, in Maume-Deschamps et al. (2016a), several ways to compute

Ψ -1 R    α 2α-1    .
Figure 2: Computation of the 0.95 gaussian expectile Now, we consider X = (X 1 , X 2 ), X 1 ∈ R N , X 2 ∈ R a (R, N + 1)-elliptical random vector. According to (5), the α-expectile of X 2 |X 1 is given by :

e α (X) = µ 2|1 + Σ 2|1 Ψ -1 R *       α 2α -1       (6)
In the general case, the term

Ψ -1 R *    α 2α-1    is difficult to compute.
This is why we propose some other predictors.

Expectile Regression

From now one, we consider the following context: (X(t)) t∈T is an R-elliptical random field. We consider N observations at locations t 1 , ..., t N ∈ T , called (X(t 1 ), ..., X(t N )). In order to predict the value of X 2 = X(t) given X 1 = X(t 1 ), ..., X(t N ), we approximate X(t) by : êα

(X 2 |X 1 = x 1 ) = β * T x 1 + β * 0 , (7)
where β * and β * 0 are solutions of the following minimization problem

(β * , β * 0 ) = arg min β∈R N ,β 0 ∈R E[φ α (X 2 -β T X 1 -β 0 )]. (8)
In the case α = 1 2 , the solution of ( 8) is exactly the kriging vector. Otherwise, this problem leads to an expectile regression, introduced by [START_REF] Newey | Asymmetric Least Squares Estimation and Testing[END_REF]. In our context of elliptical random fields, we are able to solve this minimization problem, and then define the Expectile Regression Predictor :

êα (X 2 |X 1 = x 1 ) = µ 2|1 + Σ 2|1 Ψ -1 R    α 2α-1    (9) Furthermore, its distribution is êα (X 2 |X 1 ) ∼ E 1       µ 2 + Σ 2|1 Ψ -1 R       α 2α -1       , Σ 21 Σ -1 11 Σ 12 , R       (10)

Extremal expectiles

In this section, the aim is to establish a relation between Ψ -1 R and Ψ -1 R * for extremal values of α. For that, we do an assumption : their exist 0 < < +∞ and γ ∈ R such as :

lim x→+∞ 1 -Ψ R * (x) 1 -Ψ R (x γ ) = (11)
Under this assumption, we can define Extreme Conditional Expectiles Predictors :

     êα↑ (X 2 |X 1 = x 1 ) = µ 2|1 + σ 2|1     Ψ -1 R     1 -α-1 (2α-1)         1 γ êα↓ (X 2 |X 1 = x 1 ) = µ 2|1 -σ 2|1     Ψ -1 R     1 -α (2α-1)         1 γ (12) Distribution γ Gaussian 1 1 Student, ν > 0 N ν + 1 Γ ( ν+N+1 2 ) Γ ( ν 2 ) Γ ( ν+N 2 ) Γ ( ν+1 2 )   1 + q 1 ν   N+ν 2 ν N 2 +1 ν+N ν-1 ν+N-1 Unimodal GM 1 min(θ 1 ,...,θ n ) N exp    - min(θ 1 ,...,θ n ) 2 2 q 1    n k=1 π k θ N k exp    - θ 2 k 2 q 1     Slash, a > 0 N a + 1 2 1-a 2 (a-1)Γ ( N+1+a 2 ) q N+a 2 1 (N+a)(N+a-1)Γ ( N+a 2 ) Γ ( 1+a 2 ) χ 2 N+a (q 1 )
Table 2: Some examples Thanks to the paper of [START_REF] Djurčić | Strong Asymptotic Equivalence and Inversion of Functions in the Class Kc[END_REF], we are able to prove that these predictors êα↑ and êα↓ are asymptoticaly equivalent to the theoretical expectiles respectively when α → 1 and α → 0.

Figure 1 :

 1 Figure 1: Slash random field with a Matérn 3 2 kernel

Figure 3 :

 3 Figure 3: Levels of quantile α = 0.9995 and α = 0.0005