Elevated apoptosis impairs epithelial cell turnover and shortens villi in TNF-driven intestinal inflammation

The intestinal epithelial monolayer, at the boundary between microbes and the host immune system, plays an important role in the development of inflammatory bowel disease (IBD), particularly as a target and producer of pro-inflammatory TNF. Chronic overexpression of TNF leads to IBD-like pathology over time, but the mechanisms driving early pathogenesis events are not clear. We studied the epithelial response to inflammation by combining mathematical models with in vivo experimental models resembling acute and chronic TNF-mediated injury. We found significant villus atrophy with increased epithelial cell death along the crypt-villus axis, most dramatically at the villus tips, in both acute and chronic inflammation. In the acute model, we observed overexpression of TNF receptor I in the villus tip rapidly after TNF injection and concurrent with elevated levels of intracellular TNF and rapid shedding at the tip. In the chronic model, sustained villus atrophy was accompanied by a reduction in absolute epithelial cell turnover. Mathematical modelling demonstrated that increased cell apoptosis on the villus body explains the reduction in epithelial cell turnover along the crypt-villus axis observed in chronic inflammation. Cell destruction in the villus was not accompanied by changes in proliferative cell number or division rate within the crypt. Epithelial morphology and immunological changes in the chronic setting suggest a repair response to cell damage although the villus length is not recovered. A better understanding of how this state is further destabilised and results in clinical pathology resembling IBD will help identify suitable pathways for therapeutic intervention.

where Ω is the set of n (=3) experimental conditions or mouse models, which include healthy, acute and chronic inflammation, ai is the position of the labelled front at time t=0 in condition i, and Fi is an indicator (dummy) variable with value equal to 1 for condition i and 0 otherwise.
The likelihood function for O, the set of observed positions of the labelled front in our data, was therefore given by where 'Gamma' denotes the probability density function of a gamma distribution with mean value equal to 1 and variance 10 3 . Table S1 and Fig S2 show the posterior estimates and the fitting diagnosis plots, respectively, for each experimental dataset.

Analysis of the number of dead cells in the CVEU using TUNEL labelling experiments
We sought to compare the number of dead cells quantified in control conditions with those estimated in the chronic inflammation mouse model and in the acute inflammation mouse model at the peak of the process (1-1.5h post TNF injection) and during recovery (6 h post-TNF injection) in different regions of the CVEU. We assumed that the number of TUNELpositive cells in the crypt, villus body or villus tip followed a Poisson distribution, whose parameter θ was expressed as a function of the experimental conditions by where Ω and Fi are as defined above, and bi denotes the number of dead cells per CVEU in experimental condition i. The likelihood function for O, the set of observed number of TUNELpositive cells in our experimental conditions, is therefore given by 1 1 where 'Poisson' denotes the probability density function of the Poisson distribution with parameter wθ and d is the observed number of TUNEL-positive cells in the w sampled CVEUs (See Table S2). We selected the conjugate prior distribution for this likelihood: π (bi) ~ Gamma (0.001, 0.001), for i =1..n, where "Gamma" is as defined above. Therefore, the posterior distribution for b1..bn was also a gamma distribution with parameters 0.001  (Table S2).

Estimation of cell death rates in the CVEU using TUNEL labelling experiments
To quantify the temporal dynamics of cell death in our experimental models, we developed a mathematical model that describes TUNEL labelling dynamics in epithelial cells.
We assumed that epithelial cells that lose viability become TUNEL-positive and detach from the epithelium. We defined H as the number of TUNEL-negative (healthy) cells and D as the number of TUNEL-positive (dead) cells in the CVEU. Thus, TUNEL labelling dynamics can be described by the following system of ordinary differential equations: where δ is the cell death rate (h -1 ) and γD (9) Parameter estimation was performed using Bayesian inference. We assumed that both the number of TUNEL positive and negative cells in CVEU followed a Poisson distribution with parameters equal to H (t) and D (t) described in equations (8) and (9), respectively, which were dependant on experimental conditions and collected in an unique expression, θ, as follows: where Ω and Fi are as defined above, Hi and Di correspond to the expressions in equations (8)- In addition, we assumed that the rate of detachment of TUNEL-positive cells, γD, was not affected by the experimental conditions used in this work. This implied that the detachment rate estimated in acute inflammation conditions was considered a good estimate for all our experimental conditions. This assumption was partially supported by the estimation of similar values for the detachment rate in ileum and duodenum which, in contrast, exhibited different death rates during acute inflammation (Table S3). The steady-state number of TUNEL-positive and negative cells are reported in Fig 3B and Table S4 for each experimental condition.
Estimation of cell division rate by modelling BrdU cell labelling dynamics across a CVEU described by three-compartments.
We used the three-compartment model developed in our previous work 4  We have modelled the propagation of BrdU labelling across these three compartments after a pulse of BrdU (see Fig 4A). Following BrdU injection, proliferative and nonproliferative BrdU-labelled cells are generated within the crypt and transferred onto the villus once labelled cells reach the crypt-villus boundary, which is modelled with a threshold for the number of labelled cells within the crypt. We assumed a common value for the rate of cell death in the two crypt compartments, proliferative and non-proliferative, but this rate might differ for villus cells. Cell shedding from the villus is initiated when labelled cells reach the villus tip or, equivalently, the number of labelled cells in the villus reaches a threshold value.
As in our previous work 1,4 , we restrict our attention to the period before the cell BrdU content has been diluted below the detection limit, which is observed after 4-5 generations 5 .
With these assumptions, the temporal dynamics of our model satisfy the following system of ordinary differential equations: where LCP (t) and LCQ (t) denote the numbers of proliferative and non-proliferative BrdU-       for the intervals t0 ≤ t < t1 and t1 ≤ t < t2, i.e. before labelled cells reach the villus tip. We define . These solutions are: The values of the parameters in equations (16)-(23) were identified by using BrdU datasets as well as the results on the size of the proliferative compartment and on cell death along the CVEU reported in this work.
Our experimental period was restricted to the interval t0 ≤ t < t2, before labelled cells reach the villus tip so that the shedding rate of labelled cells from the tip during this period was descendants. However, the estimation of these two parameters is feasible under the assumption that the ratio between the number of proliferative, CP, and non-proliferative, CQ, cells in the crypt is constant in our experimental conditions, which is supported by our previous results ( Fig 2E) and that BrdU labelling does not affect the rate at which cells enter the quiescent state.
Under these assumptions, we have that As described above, we assumed that the ratio between the number of proliferative, CP, and non-proliferative, CQ, cells in the crypt is assumed to be constant under our experimental conditions and not affected by BrdU uptake and given the structure of our model we have that . Under these assumptions, equation (24) can be expressed as: Estimates for CP and CQ were obtained as described above. The replacement of γC in equations (16-23) by the expression in equation (25)  To compare statistically the value of the parameter λ, cell division rate, between our mouse models, we conducted Bayesian inference to estimate the posterior distribution of λ using MCMC methods. As in previous sections, we assumed that the number of labelled cells in both the crypt and the villus had a Poisson distribution with parameter θ, which was dependent on time and experimental conditions as follows: , , where Ω and Fi are as defined above; LC,i = LCP,i + LCQ,i and LV,i correspond to the number of BrdU-labelled cells in the crypt and villus, respectively, described in equations (16)-(23) for experimental condition (mouse model) i; ξC is an indicator variable with value equal to 1 for crypt cells, or, otherwise 0.
The distribution function for the likelihood of the dataset was: where lC,h and lV,h are the observed number of labeled cells in the crypt and villus, respectively, in the wh samples collected at the th sampling time and m is the number of time points. We selected non-informative prior distributions for the parameters: π (λi) ~ Gamma (0.001, 0.001), which is defined above.
Fitting diagnosis plots and posterior estimates can be found in Fig S5 and Table S5, respectively.