Inverse problem for cell division rate in population dynamics

. In population dynamics, more precisely in cell structured population, experimental measures of the temporal dynamics of some parameters are often difﬁcult or even impossible because of the lack of appropriate technology. It is the case of the division rate parameter for a structured population of cells. So theoretical and numerical approaches are needed to recover such parameters by indirect observation of the solution of the modeling system. In this work we present how to estimate the division rate from a growth-fragmentation type equation as done in the work of Doumic and Tine (see [1]) and following the work done in [2] and [3]. The estimation strategy is based on the observation of the long-term solution of the model for a size structured population when the growth and division rates are fully general. The issue of the theoretical and numerical results are discussed here.


Introduction
The observation of a population of cells (proteins, bacteria, etc.) allows to understand the dynamics of the considered population such as its distribution, the behavior of its parameters (growth rate, death rate, division rate).This understanding depends on technologies that biologists have for experiences.For instance, one knows that from cytometric flow technology it is possible to measure at each time the distribution of the population what is not the case for instance for the measurement of the division rate because of the lack of technology able to follow a single cell from its birth to its division.So, to estimate this division rate which is a crucial parameter (in the case of pathogenic cells as E. Coli bacteria), mathematical and numerical approaches are needed.Then without lost of generality let consider the following growth-fragmentation type model: where n(t, x) is the density repartition function of the cells structured by their size variable x at time t.We assume the cell population to evolve only by growth and division process without death and birth events.The growth rate g is a given or measured function which depends on the size x.The total rate a e-mail: leon-matar.tine@univ-lyon1.frITM Web of Conferences of division for cell of size x is denoted by B(x), so the product B(y) (x, y) is the rate at which a single cell of size y divides in to cells of size respectively x and y − x.The kernel division (x, y) is required to satisfy: With this model (1) we can formally deduce by simple integration that, in one hand the time variation of the total number of cells depends on the division rate B what means In the other hand, one deduce that the time variation of the total mass of the population depends on the growth rate g what means The aim of this paper is to show how to estimate the division rate B from the knowledge of measurements on the density repartition n(t, x) and on the growth rate g(x).Let us recall that these assumed known quantities can be experimentally measured by cytometric flows technology, while the division rate need theoretical and numerical approaches for its estimation.Estimating this division rate is crucial if one is interested of pathologic cells (as bacteria or misfolding proteins) which cause disease.Work in this theoretical estimation exist for the simple case g(x) ≡ 1 and where the individual cell divides with equal mitosis (see [2,3]), that means, (x, y) = x=y/2 .This equal mitosis change the integral term in (1) into 4B(2x)n(t, 2x).Here we investigate the case where the kernel division is quite general and the growth rate g not necessarily constant.Of course, for the model (1), there are many results on the well-posedness and the essential properties.For more detail one can refer to [5][6][7][8][9].The fact that, at any time we can have observation of the size distribution of the cell population allows us to focus the estimation of the division rate on the time asymptotic behavior of the model (1).So the paper is organized as follow: in Sect. 1, we present the asymptotic behavior of the model and related results proved by generalized relative entropy principle (GRE).In Sect. 2 we present the inverse problem of estimation the division rate B thanks to two regulation methods.The Sect. 3, is devoted to the numerical methods and results.The last section is for the conclusion.

Asymptotic behavior of the model
The asymptotic behavior of the model ( 1) is based on the GRE principle (see [8]) which need basic hypothesis on the parameters.
With the notation we assume: the kernel division fulfill in addition to (2), the relations The growth rate and the division rate satisfy WMLS 2014 and for all compact set The cells of size zero forming is avoided thanks to the assumption For the cells of small size, the growth is assumed to be more important than de division, that means For the cells of big size, in order to avoid the formation of infinite size cells at finite time (gelation phenomenon), the growth and division rates satisfy Thanks to the general assumptions ( 3)-( 8), the existence, uniqueness and properties of eigenelements associated to (1) is proved in [4].In particular it is shown by the GRE principle that where ( , N, ) is the unique eigentriplet solution of the following stationary equation and its adjoint equation: The estimation approaches of the division rate we develop here is based on the problem ( 9) and the assumption on the knowledge of g, N or their experimental measurements.Of course there is two fundamental estimates for the eigenvalue from simple integration of (9) which are written as follow

Estimation of the division rate
In this part, we investigate the question of recovering the division rate B when the growth rate g is given and a sequence of measurements N such that, with an adequate norm they fulfill, N − N < and | − | < for and N satisfying Eq. ( 9).

ITM Web of Conferences
The possible noise on the measures can make the problem ill posed what is inherent to almost inverse problem.What's happen in the problem (9) when N is replaced by N and one wants to quantify the term * x (g(x)N (x)).To be more clear, lets consider the following sample.
Example 3.1: Let consider f ∈ C ∞ (R) and for arbitrary n ∈ N * .We define for all > 0, the approximation f (x) = f (x) + sin nx .It is obvious to see that f − f ∞ ≤ means that f as an good approximation of f .But when one looks to the derivative f (x) = f (x) + n cos nx , it fulfill f − f ∞ = n and we imagine when n is chosen as big as possible.
In order to come up this difficulty of controlling the term * x (gN ), we propose two methods of regularization named Filtering and Quasi-reversible methods that we present here below by introducing the new unknown H = BN, solution of the approximated problem

Filtering method
Here, we introduce the mollifier and define the following regularized problem For this regularized problem we prove the existence and uniqueness (see [1]) by using the Lax-Miligram theorem which is based on the bilinear form x x p u(x)u(y) (x, y) dy dx, that we show the continuity and coercivity.In addition we have by priori estimates and Minkowski inequality, the following approximations where C(p, , N) is a constant depending on p, and N .

Quasi-reversible method
The key point of this approach is to add another derivative term in the problem in order to balance the term source of ill posed term.So the regularized problem is casted as follow for k ∈ R + , k ≥ 2 and > 0.
Here, we prove the existence and uniqueness of the solution (see [1]) by the Banach fixed point theorem.We begin by writing the regularized problem with a general unknown u as: WMLS 2014 where F (x) = −c , * x gN − N .So we define the application in L 1 (R + , x dx) such that We prove that T is a contracting Lipschitz function.
In addition we have by priori estimates and Minkowski inequality, the following approximations where C(p, N) is a constant depending on p and N .

Numerical simulation for the inverse problem of estimating B
This part is devoted to the numerical test of the estimated division rate B from the two proposed regularization methods.Experimental data on the long time behavior of the cell repartition function and on the growth rate are needed.In [1], we have no data available, so we numerical generate synthetic data on N by solving first the direct problem (1) with given g, and B; what allows to obtain the asymptotic profile of the size distribution function of cells.Once this asymptotic profile is obtained, we keep the same parameters g, and we forget the division rate B that we want to recover.By adding noise on N , we numerically recover the division rate B thanks to our two regularized methods and compare B and B which was taken into account in the resolution of the direct problem (1).

Results on the direct problem
The simulation of the direct problem follows a finite volume method applied to (1).For the length domain L and grid number points n max we have . So, for a time step t we define t k = k t, k ∈ N and take the following approximations We use time splitting and solve first the growth part of (1), * t n + * x (g n) = 0, as follows where (gn) k i±1/2 are incoming and outgoing flows.In second, we solve the gain and lost part by cell division as follows In order to get the asymptotic profile, we renormalize the discrete solution: For different initial conditions and parameter values on g, and B satisfying the general assumption of existence and uniqueness of solution our scheme converge to the same asymptotic behavior as one can see in figure (Fig. 1) below.Let recall that the time step is computed with respect the CFL (Courant-Friedrichs-Lewi) stability condition, t < 1/ max 0≤i≤n max B i + g i x , with the notation B i = B(i x) and g i = g(i x).Once results on the direct problem are obtained, means the asymptotic behavior, we numerically simulate the regularized problems in order to estimate the division rate.

Results on the regularized problems
For the quasi-reversibility method, the associated regularized problem ( 12) is discretized as follows WMLS 2014 ).Down center: g(x) = x 1/3 , B is here a tray-function (see 13).
We then rewrite this discrete equation under matrix shape A × H = L where A appears as a triangular matrix, so we deduce the solution H by solving this simple iteration The non-oscillatory of this scheme is ensured by the constraint A i,i = 1 + i − 2 (i, i) x > 0. For the filtering method, we apply the same discretization as in quasi-reversibility case to (11  We then write the discrete equation of (11) under matrix shape A × H = L with A an upper triangular matrix and deduce the approximation of H = BN as follows The numerical results show a good accuracy of our regularization methods as depicted in the pictures below in the case where we assume that we have exact measures ( = 0) that means N = N .These pictures compare the exact value on B taken in the direct problem and its estimated values obtained by 01003-p.8 WMLS 2014 solving the quasi-reversibility or filtering methods.We define a Tray-function for B as follows In order to take into account the uncertainty in measurement, we add some multiplicative noise which is uniformly distributed in N = max(N(1 + l ), 0) with l ∈ −1 2 , 1 2 and ∈ [0, 1].
So, with N considered as the measured data, we obtain the following results with various choice of and the fixed parameters g(x) = x, (x, y) = 1 y I {x<y} , B(x) = exp(−0.08(x− 13) 2 .In the simulations we investigate for any , the optimal regularization parameter which gives the best estimation B , of the division rate B by minimizing the L 2 -norm of the following relative error

Conclusion
This study presents how to recover the division rate for a general growth-fragmentation model for a size structured population where measurements on the time asymptotic profile of the density.Two method of regularization are highlighting here in order to overcome the difficulty of the ill-posed inverse problem.Theoretical estimates are built with strong analysis on the well-posedness of the regularized problem and numerical methods are implemented.As show on the above numerical illustrations, the estimation of the division rate is quite good.However, some issues remain to be addressed, specially for small size x, our algorithm need to be performed such as the appearance of oscillations.