A mathematical model of cell fate selection on a dynamic tissue
AA mathematical model of cell fate selection on adynamic tissue
Domenic P.J. Germano a , James M. Osborne a, ∗ a School of Mathematics and Statistics, The University of Melbourne, Parkville, Victoria3010, Australia
Abstract
Multicellular tissues are the building blocks of many biological systems and or-gans. These tissues are not static, but dynamically change over time. Even if theoverall structure remains the same there is a turnover of cells within the tissue.This dynamic homeostasis is maintained by numerous governing mechanismswhich are finely tuned in such a way that the tissue remains in a homeostaticstate, even across large timescales. Some of these governing mechanisms includecell motion, and cell fate selection through inter cellular signalling. However,it is not yet clear how to link these two processes, or how they may affect oneanother across the tissue. In this paper, we present a multicellular, multiscalemodel, which brings together the two phenomena of cell motility, and inter cel-lular signalling, to describe cell fate selection on a dynamic tissue. We find thatthe affinity for cellular signalling to occur greatly influences a cells ability to dif-ferentiate. We also find that our results support claims that cell differentiation isa finely tuned process within dynamic tissues at homeostasis, with excessive cellturnover rates leading to unhealthy (undifferentiated and unpatterned) tissues.
Keywords:
Multicellular modelling, Multiscale modelling, Cell fate selection
1. Introduction
Organs are comprised of dynamic multicellular tissues, with many healthytissues existing in a state of dynamics homeostasis. This balance is achievedthrough multiple interacting governing mechanisms. One tissue where this bal-ance is crucial is the intestinal epithelium, which protects the intestine and colon ∗ Corresponding author
Email addresses: [email protected] (Domenic P.J. Germano), [email protected] (James M. Osborne) a r X i v : . [ q - b i o . T O ] J un gainst digestive by-products, and also aids in the absorption of water and nu-trients [1]. The intestinal epithelium is one of the fastest self renewing tissuesin the human body, and this renewal is controlled within millions of test-tube-like structures which line the intestinal walls, known as the crypts of Lieberkhn[16]. At the base of each crypt, resides a stem cell population, which is directlyresponsible for the replenishment of the epithelial cells through mitosis. Afterdividing, cells migrate up the walls of the crypt, differentiating as they go, whereupon reaching the top of the crypt they form the epithelial lining of the intestineand are removed by sloughing to prevent overcrowding. In humans, this renewalprocess usually takes 6–7 days [16]. Therefore, healthy crypt homeostasis is con-trolled by the fine balance of cell proliferation, migration, differentiation, andapoptosis [16].As cells migrate up the crypt, cell differentiation and cell fate selection iscontrolled by another finely tuned mechanism, known as Delta-Notch signalling,which occurs between neighbouring cells [4]. The key mechanisms controllingDelta-Notch signalling are the Notch receptors and Delta ligands, both of whichare transmembrane proteins. The pathway is activated when Delta ligandson one cell bind to Notch receptors of a neighbouring cell. As a response,intracellular reactions are triggered on the said neighbouring cell, which in turnallows target gene expression, and thus leads to cell differentiation [2].Mutations of cells within the crypts (particularly the stem cells) may triggerthe development of colorectal cancer [10, 15]. For this reason, an understandingas to how crypt dynamics are maintained at homeostasis, and what may leadto malignancies is highly desirable. Mathematical and computational modelsprovide an ideal framework to study these dynamical tissues. Previous mod-elling work has provided insight into how Delta-Notch signalling occurs withindynamic tissues. One early model of Delta-Notch signalling was proposed byCollier et al. [4], where a tissues of static cells interact with their nearest neigh-bours through lateral inhibition. The authors found that their model reproducedNotch patterning, analogous to those found in living systems. In 2011 Buskeet al. [3] coupled a logic based Delta-Notch signalling with an over-lappingspheres model of cell dynamics with position dependent proliferation. Usingtheir model, the authors were able to reproduce the correct ratios of cell typespresent at steady state. The model can also describe and predict dynamicbehaviour of the epithelium at both steady state, and also following the intro-duction of mutant cells. The authors also showed that the intestinal epitheliumis capable of complete recovery following eliminations of each subpopulations of2ell within the crypt epithelium. The Collier et al. model for Delta Notch hasbeen coupled with numerous multicellular models. Osborne et al. [14] present acomparison of coupling the signaling model with five biomechanical models andshow that as long as the use of biomechanical model is appropriate (i.e there areno model artefacts), then the biomechanical models are equivalent. The Collieret al. model was extended to include other aspects of signalling in the Crypt inKay et al. [9]. In the paper the authors analyse this mode on pairs of connectedcells and show that the Delta Notch patterning can influence cell differentiation.Furthermore, it is know that during development, cell morphology changes,and therefore the contact geometry between neighbouring cells changes too.Since Notch signalling is mediated via transmembrane proteins, these morpho-logical changes could influence cellular communication. To investigate, Shayaet al. developed a model of Notch signalling which is dependent upon the con-tact area between neighbouring cells [17]. They found that contact area biasescellular differentiation, where smaller cells are more likely to differentiate intothe primary cell fate.It is clear that understanding the interplay between cell fate signalling andcell dynamics is crucial in furthering our knowledge of tissue and organ develop-ment and function. To the authors knowledge no previous study has investigatedhow properties of the inter-cellular signalling, and cell turnover, influences celldifferentiation and fate selection. In this paper we present a multiscale multi-cellular model which couples cell dynamics with Delta-Notch signalling and cellfate selection.The remainder of this paper is structured as follows, we begin by presentingour multicellular multiscale model of cell fate selection in a dynamic tissue inSection 2. In Section 3, we first present an investigation into how subcellulardynamics and tissue geometry influences pattern formation on static tissues,before demonstrating the influence of cell turnover on patterning. Finally inSection 4, we discuss our results and relate them back to biological dynamicaltissues like the colorectal crypt.
2. Model
Here, we model the tissue as a collection of discrete, interacting, individualcells, more commonly referred to as a multicellular model. In a multicellu-lar model, cells are represented as a single point (or a collection of points) inspace, and allows details at the cellular level and tissue level to be included [14].3pecifically, they allow the inclusion of cell population turnover and cellular sig-nalling, which occur at differing timescales.Below, we discuss the multicellular,multiscale model we use to couple these processes, in order to analyse cell fateselection within a dynamic tissue.
To describe cell dynamics, we use a lattice-free, cell-centred model [12]. Thenet force on a given cell i , F Net i , is found by balancing the force due to neigh-bouring cell interactions, F Interactions i , and the viscous forces on the cell, F Viscous i ,as proposed by Meineke et al. [12]: F Net i = F Interactions i + F Viscous i , ∀ i. (1)Due to the highly viscous environment that cells occupy, we assume that cellmotion is over-damped [5], and therefore viscous forces dominate allowing us toneglect inertial terms and all motion is determined by a force balance, F Net i = . We follow Meineke et al. [12] and model the force due to neighbouring cellinteractions using a linear Hooke’s law acting at cell-centres to describe attrac-tion and repulsion between neighbouring cells. We also assume viscous forcesacting on the cell oppose the direction of motion. This leads to the equationsof motions, for cell i , at position r i : ν d r i dt = (cid:88) j ∈ M i k spij ( | r ij | − s ij ) ˆ r ij , ∀ i, (2)where r ij = r j − r i is the displacement between cells i and j , ν is the dragcoefficient applied to all cell centres, M i the neighbouring cells of cell i , k spij thespring constant between cells i and j (here taken to be constant so k spij = k sp ),and s ij = s ij ( t ) is the length separation between cells i and j , which has theform: s ij ( t ) = ε + τ i (1 − ε ) , τ i ≤ , , otherwise , (3)with ε > τ i ≥ i at time t . All parameters are given in Table 1. Note all distancesare measured in cell diameters (cd) which we take to be the average size of acrypt epithelial cell 10 µm [5]. 4 .1.2. Cell Population Turnover Our aim is to get a cell turnover rate of γ divisions (and deaths) per hour percell. In order to do this we model each cell to have its own cell cycle duration, T , which is sampled from a Uniform distribution : T ∼ U (cid:18) γ , γ (cid:19) , (4)where γ is the target cell turnover rate. When a cell reaches its cell cycle dura-tion, it is labelled the parent cell, with position r p , and subsequently proliferatesinto two daughter cells, with positions r i and r j , displaced at a separation ε > n : r i = r p + ε n , r j = r p − ε n . (5)To represent a turnover of cells, and to maintain a fixed number of cells in thetissue, whenever a cell divides we randomly select a cell to be removed in thesame timestep. The model of intracellular signalling we employ was initially described byCollier et al. [4]. Collier et al. consider a simplified model of a cell, focusing onthe Delta-Notch pathway, where inhibited cells have a reduced ability to inhibitother cells. The key mechanisms of the Delta-Notch pathway are the Notchreceptors and Delta ligands, both of which are transmembrane proteins. Thepathway is activated when Delta ligands on one cell bind to Notch receptors ofneighbouring cells. As a response, intracellular reactions are triggered, whichin turn allows target gene expression, and thus leads to cell differentiation [2].The biochemical model assumptions are summarised by the following:1. Only cells in direct contact may interact via Delta-Notch signalling.2. The rate of Notch production is an increasing function of the amount ofDelta present in neighbouring cells.3. The rate of Delta production is a decreasing function of the amount ofNotch within the same cell.4. The rate of both Notch and Delta decay obey exponential laws, with rates µ and ρ respectively. similar results can be obtained by sampling from a Gamma distribution of T ∼ Γ (48 , γ ).
5. A cell’s fate is determined by the amount of Notch within the cell.These model assumptions give rise to the following (non-dimensional) mathe-matical model, for cell i , with Notch level N i and Delta level D i : dN i dt = µ ¯ D ki a + ¯ D ki − µN i , ∀ i, (6) dD i dt = ρ
11 + bN hi − ρD i , ∀ i. (7)where a is the Notch affinity constant (which dictates how readily Delta ligandsbind to Notch receptors), and b the Delta affinity constant (which dictates howreadily Delta ligands are produced), and k and h the exponent for Notch andDelta synthesis respectfully. Lastly, ¯ D i is the mean level of Delta within neigh-bouring cells (of cell i ). On cell division the level of Delta and Notch in thedaughter cells are assigned to be the same as the parent cell. We constrain the tissue to lie in a domain of size L x × L y , where L x , L y ∈ R ,with horizontally and vertically periodic boundaries (a toroidal domain). Dueto the honeycomb structure of centre based cells in equilibrium it is more in-formative to describe the tissue size by referring to the number of horizontallyand vertically stacked cells, C x × C y , where C x , C y ∈ N . The relation betweenthe number of cells horizontally stacked ( C x ), and the horizontal length ( L x ) is L x = C x . However, because we initialise the tissue on a honeycomb, hexagonallattice, which is the equilibrium state for our biomechanical model [18], the rela-tion between the number of cells vertically stacked ( C y ), and the vertical length( L y ) is L y = √ C y . We treat each cell as a distinct agent, which interacts withit’s neighbouring cells. These neighbouring cells are determined by a DelaunayTriangulation between cell centres. The shape of each cell is given by a Voronoitessellation, which is the natural dual of the Delaunay Triangulation. Figure 1shows a typical tissue geometry with L x = 6cd and L y = 3 √ C x = C y = 6).6 Figure 1: A typical tissue geometry. The black dots show cell centres. The dashed linesshow the neighbouring cells defined via a Delaunay Triangulation, with black dashes showingimmediate neighbouring cells and red dashes describing neighbours due to periodicity. Thesolid blue lines describe the cell shapes, given by a Voronoi tesselation.
Unless stated otherwise, the initial conditions we use for the biochemicalmodel are uniform homogeneous. That is, for cell i , the initial levels of Notchand Delta within the cell are: N i (0) = N , D i (0) = D , ∀ i, (8)where 0 ≤ N ≤ ≤ D ≤
3. Results
We first consider the behaviour of a static tissue (i.e. cells have no turnoverrate). Therefore, cells remain stationary and only interact with a fixed set ofsix neighbouring cells. In Section 3.3 we relax this assumption and allow cellsto proliferate, undergo apoptosis and move.
In their 1996 paper [4], Collier et al. explained that a default patterning ofone primary cell (low Notch level) for every two secondary cells (high Notchlevel) is the dominant ordering of cells at steady state, as shown in Figure 2c.This ordering is such that the neighbouring set of each primary cell is exactlysix secondary cells (see Figure 2k), while the neighbouring set of each secondarycell consists of three primary and three secondary cells in an alternating fashion(see Figure 2l). 7 o m og e n e o u s Steady State Notch (a)
Intracellular Notch Solution (b) S i n g l e S ee d (c) (d) R a nd o m (e) (f) R a nd o m (g) (h) U nfi t D o m a i n (i) (j)(k) (l) Figure 2:
In silico experiments with parameter values a = 0 . b = 100, µ = ρ = 10, k = h = 2, γ → ∞ , ν = 1 and k sp = 50. a, c, e, g and i show the steady state Notchpatterning for differing initial conditions of Notch. b, d, f, h and j show the associated timeseries solutions of intracellular Notch. Note d includes inset of early time solutions to showinitial heterogeneity. Videos of these simulations can be found in SI Movie 1 – SI Movie 5. kshows a primary cell, and l shows a secondary cell in a default patterning. − ) then the default Notch patterning was obtained, see Figure2c. Moreover, looking at the time series solutions for the default patternedtissue (Figure 2d), we can see that the patterning propagates from the initiallyseeded cell (see SI Movie 2 for a video of the simulation).It should be noted, that the system is quite sensitive to initial conditions ofDelta and Notch. This can be seen by running two simulations from differentinitial random initial conditions, Figures 2e (2f, see SI Movie 3 for a video of thesimulation)) and 2g (2h, (see SI Movie 4 for a video of the simulation)). Theparticular random conditions used in Figure 2e are able to achieve the defaultNotch patterning, and we see from the time series solution Figure 2f that thetissue is able to reorient itself and establish distinct patterning. However, thoserandom initial conditions used within Figure 2g are suitable, resulting insteadin a partly patterned tissue, containing unresolved errors, as we see from thetime series solution Figure 2h, the tissue does attempt to achieve the defaultpatterned state, but distinct (stable) errors occur (see SI Movie 5 for a video ofthe simulation).Patterning initiates from irregularities within the Delta and Notch levels ofcells, this irregularity then propagates radially. If we were to consider infinitetissue structures, then the default patterning shown in Figure 2c would emerge,referred to as a period 3 pattern [4]. However, when considering periodic bound-ary conditions, which is often the case in in silico experiments [13, 7], combinedwith the understanding of how irregularities propagate, we find that the pat-terning present in Figure 2c is not always possible. We see from simulationson different tissue geometries starting with a single seeded cell (the most stableinitial condition for patterning) that patterning is not possible on domains ofarbitrary size. For example, Figure 2i shows the pattern generated on a tissueof 10 × C x = 10 and C y = 8), which show errors in the patterning. Thetime series solutions, Figure 2j shows that initially, the patterning propagatesfrom the initial irregularity. However, as this propagation meets up with itself,a discrepancy occurs which results in the error shown down the middle of the9issue.To determine the suitable geometry size, C x and C y , such that the tissuesupports a default pattern, we first observe that the patterning is a period 3pattern, meaning that if we start at a given cell, and we move away from it ina straight line, then as we move, we would notice that the pattern repeats itselfevery 3 cells. Due to the packing of hexagonal lattices with periodic boundaryconditions, we require that the number of vertical cells be C y = 2 n , where n ∈ N . We note that there are three axes of symmetry on hexagonal lattices:the horizontal axis, another at 60 o to the horizontal, and another at 120 o to thehorizontal. Due to the periodicity of the tissue, we now define what is called atorus knot [11], but on a lattice: if we start at a given location, and move alongone of the axes of symmetry, we will always end up back where we started. Forthe tissue to support patterning, we require that the length of all torus knotsalong each of the 3 axes of symmetry be divisible by three. The easiest and mostobvious torus knot is horizontal and is shown in Figure 3a, which for a generaltissue has size C x × n . Requiring that the length of this knot be divisible by 3,means we must have C x = 3 m , where m ∈ N , resulting in tissue sizes of 3 m × n .However, it is not immediately obvious that the length of the non-trivial knots,shown by Figures 3b and 3c, are also divisible by 3. The length of these knotscan be written as: knot length = size of tissuenumber of disjoint knots . (9)For a tissue of 3 m × n cells, the size of the tissue is simply 6 mn . To find thenumber of disjoint knots, we first observe that for every 2 n steps we traversediagonally, we advance n rows horizontally, due to periodicity. Therefore, thenumber of disjoint knots which can be supported on a tissue of size 3 m × n isgoing to be given by:number of disjoint knots = gcd (3 m, n ) , (10)which gives the length of the knots as:knot length = 6 mn gcd (3 m, n ) , (11)which is always divisible by 3. 10 a) (b) (c) Figure 3: Examples of knots on a 6 × Therefore, we can say that only tissues of sizes 3 m × n cells, where m, n ∈ N will support a default Notch pattern shown in Figure 2c. Tissues of other sizescannot support a default patterning, but rather a partly patterned tissue witherrors, analogues to that shown by Figure 2g. The simulations presented inFigures 2c and 2i (along with other geometries, not shown for brevity) confirmthis. a, b Control Cell Patterning
In previous work, Collier et al. stated that the affinity constants of bothNotch and Delta control the existence of cell fate patterning, but did not statehow. By further analysing the multicellular system, we are able to describe theeffects that the affinity constants have on the existence of cell fate patterning,specifically we can say which parameter values allow patterning. By consideringa static tissue which exhibits a default patterned state, such as that of Figure2c, we observe that there are at most two different cell configurations, up torotation. The first is that of the primary fate cell, shown in Figure 2k, and thesecond is the secondary fate cell, shown in Figure 2l. We therefore can writedown the governing Notch and Delta equations of a primary fate cell in termsof N p and D p , and similarly those for a secondary fate cell in terms of N s and D s as follows: dN p dt = µ (cid:32) ¯ D kp a + ¯ D kp − N p (cid:33) , dN s dt = µ (cid:18) ¯ D ks a + ¯ D ks − N s (cid:19) ,dD p dt = ρ (cid:18)
11 + bN hp − D p (cid:19) , dD s dt = ρ (cid:18)
11 + bN hs − D s (cid:19) . (12)11e then make the observation from Figure 2k that the neighbours of all primarycells, consist of six secondary cells, to then write ¯ D p as the following:¯ D p = 16 (cid:88) j ∈ M p D j = 16 (6 D s ) = D s . (13)Similarly, the neighbours of all secondary cells consists of three primary andthree secondary cells, arranged in an alternating fashion, to then write ¯ D s asthe following: ¯ D s = 16 (cid:88) j ∈ M s D j = 16 (3 D p + 3 D s ) = 12 ( D p + D s ) . (14)Substituting Equations (13) and (14), into Equation (12), leads to the simplifiedsystem: dN p dt = µ (cid:32) D sk a + D sk − N p (cid:33) , dN s dt = µ (cid:32) (cid:2) ( D p + D s ) (cid:3) k a + (cid:2) ( D p + D s ) (cid:3) k − N s (cid:33) ,dD p dt = ρ (cid:32)
11 + bN ph − D p (cid:33) , dD s dt = ρ (cid:18)
11 + bN sh − D s (cid:19) . (15)Using the above set of equations, we were able to efficiently numerically deter-mine how the affinity values affect the existence of cell fate patterning. Steadystates for Equation (15) are shown, for varying affinity parameters, in Figures4a and 4b . 12 teady State Notch with a , b = 100 N o t c h Steady State of Notch with a, b = (a) Steady State Notch with b , a = 0 . N o t c h Steady State of Notch with b, a = (b) Relationship between a crit and b crit b crit a c r i t Relationship between a crit and b crit PatterningNo Patterning
MulticellularAlgebraic (c)
Figure 4: a and b shows the steady state Notch levels of the primary cells (blue), secondarycells (red) and the homogeneous level (black). Solid lines show the stable equilibrium anddashed lines show the unstable equilibrium, not numerically observed. a shows the Notchlevels with varying affinity rate a , and b shows the Notch level with varying affinity rate b . c shows the relationship between the critical values a crit and b crit and the region withinparameter space where patterning is permitted. The red solid line is results obtained from thesimplified equations, and the blue markers are results obtained from multicellular simulations.Parameter values of k = h = 2 and µ = ρ = 1. We define the critical values a crit and b crit to be the bifurcation points in a and b respectfully where Notch patterning ceases to exist. For example, consid-ering Figure 4a, for a fixed Delta affinity of b = 100, the critical Notch affinity is a crit ≈ . a = 0 . b crit ≈ . a crit and13 crit can be found by solving the simplified system, given in Equation (15), in-voking the steady state condition, and equating N p = N s and D s = D p . Theresult is a polynomial in both the critical values a crit and b crit , which can thenbe solved numerically. The relationship is shown in Figure 4c by the red solidline. The same relationship may also be found by solving the full multicellularsystem to steady state, and performing a parameter sweep to identify both a crit and b crit , also shown in Figure 4c with blue markers.From the results shown in Figure 4 and the above discussion, we can thereforesay that the affinity constants, a , which controls how readily Delta ligandsbind to Notch receptors, and b , which controls how readily Delta is produced,completely govern the existence of Notch patterning on appropriate domains (i.e.domains where patterning is possible). When Delta-Notch binding events do notreadily occur, or Delta is not easily expressed, then cells cannot differentiate. We now consider how a finite cell turnover rate, γ , affects the patterning(and therefore differentiation) of the tissue in a dynamic steady state. We fixthe biomechanical and biochemical parameter values of the model to those ofTable 1. Note that the affinity parameters are chosen to allow patterning.Parameter Description Value µ Rate of notch degradation 10 ρ Rate of delta degradation 10 a Notch affinity 0.01 b Delta affinity 100 k Notch synthesis exponent 2 h Delta synthesis exponent 2 ν Drag coefficient 1 k sp ij Spring strength 50 ε initial cell separation 0.001 I p Notch primary cell threshold 0.1 I s Notch secondary cell threshold 0.6
Table 1: Biomechanical, biochemical and patterning proportion parameter values.
To quantify the patterning throughout the tissue, we specify a Notch thresh-old for differentiated cells, I p , for primary cells and I s for secondary cells, with I p < I s . We then categorise each of the cells within the tissue as either primary14 N i ≤ I p ), secondary ( N i ≥ I s ), or undifferentiated ( I p < N i < I s ). Thesethreshold values can be determined by solving Equation (15) to steady state fora specified set of parameters (see Table 1). One then chooses I p and I s suitably.For example, considering the biochemical parameter values of a = 0 . b = 100, k = h = 2, then the steady state Notch values for primary cells is N p ≈ . N s ≈ . I p = 0 . I s = 0 . t ), is then given by the total number of differentiated cells, divided by thetotal number of cells within the tissue at time t . However, due to the dynamicnature of the system, we further consider a moving average of the instantaneouspattering, Ω( t ), as: ¯Ω ( t ; δ ) = 1 δ (cid:90) tt − δ Ω( τ ) dτ. (16)Where 0 < δ < t is the length of the interval we sample. Figure 5 showsan example of two tissues evolving over time with different cell turnover ratesof γ = 0 . γ = 1 respectively. Figure 5e shows the instantaneous (blue)and average (red) proportion of patterning throughout the tissue. We observethat the cell turnover rate, γ influences the proportion of patterning withinthe tissue. Specifically larger levels of cell turnover result in less patterning(differentiation). 15 = . (a) t = 0hrs (b) t = 40hrs (c) t = 80hrs (d) t = 120hrs Ω Patterning of Tissue with Time time (hrs) (e) γ = (f) t = 0hrs (g) t = 40hrs (h) t = 80hrs (i) t = 120hrs Figure 5: Figures a - d show a typical in silico experiment of a toroidal tissue, with acell turnover rate of γ = 0 . in silico experiment of a toroidal tissue, with a cellturnover rate of γ = 1 division per hour per cell (see SI Movie 7 for a video of the simulation).Figure e shows the corresponding evolution of the patterning for both γ = 0 . γ = 1 (bottom), on a tissue size of L x = 12cd and L y = 6 √ To Further investigate the influence of the cell turnover rate on tissue pat-terning, we performed multiple in silico experiments, with varying cell turnover16ates of γ ∈ (cid:2) − , (cid:3) , results are shown in Figure 6, on a tissue size of L x = 12cd and L y = 6 √ γ , is increased theDelta-Notch patterning is inhibited and cell differentiation suffers as a result.This is because as γ increases, the neighbouring cells any given cell interacts withchanges more rapidly, preventing the cell from achieving chemical equilibrium,which results in undifferentiated cells and hence inhibited Notch patterning.Specifically, as we increase the cell turnover rate from γ = 10 − , to γ = 0 . γ = 0 . γ = 2, the tissue transitions from being in a patterned state to a ho-mogeneous state, exhibiting no distinct cell fates. Further increases in γ simplyresults in the homogeneous state. Effects of Cell Turnover Rate on Patterning -3 -2 -1 Figure 6: Effects of cell turnover rate on patterning on a toroidal tissue, on a tissue size L x = 12cd and L y = 6 √ in silico experiments. Blue error bars represent the 95% confidence interval.
4. Discussion
Throughout this paper, we have presented a mathematical model for cell fateselection on a dynamic tissue. We couple the model of cell signalling via contact17nhibition of Collier et al. [4], with the model of cell dynamics of Meineke et al.[12].On static tissues, we extended the work of Collier et al. [4] to observe thatboth the tissue geometry, and initial Delta-Notch conditions govern the steadystate distribution within the tissue. Specifically, we see that both can disruptregular patterning, and even inhibit patterning altogether. We further observedthat the affinity constants, a (affinity for Delta ligands to bind to Notch recep-tors) and b (affinity for Delta expression) completely govern Notch patterningon appropriate domains. Through an exploration of the parameter space, wedetermined the critical thresholds, a crit and b crit which permit Notch patterning.These results suggest that when the Delta-Notch binding event does not occurreadily, then cells do not become differentiated and cell fate selection does notoccur. Similar results occur when Delta ligands are not readily produced.Moreover, by observing a rotational symmetry of the steady state Notchpatterning, we were able to deduce a reduced system of ordinary differentialequations which is capable of capturing the steady state behaviour of the sys-tem, while significantly reducing the computational complexity in comparisonto the full multicellular system. We were then able to find an algebraic expres-sion relating the bifurcation parameters a crit and b crit .Lastly, we looked at how cell turnover rate affects Notch patterning withindynamic tissues. We found that a fast cell turnover rate inhibits Notch pat-terning, resulting in a homogeneous tissue. For cell turnover rates ranging from γ = 0 . γ = 2 divisions per hour, the tissue transitionsfrom being almost fully patterned to a homogeneous tissue, with no distinctNotch patterning. These results suggest that the cell turnover rate plays a cru-cial role in maintaining healthy tissues, with higher uncontrolled cell turnoverrates leading to a malignant, unhealthy epithelium.Future avenues of study include an analysis of how realistic tissue geometries,such as those of the colonic crypt, behave. Specifically, within the crypts, it isknown that the level of Wnt signalling a cell receives governs cell proliferationwithin the crypt [8]. It is further known that Wnt signalling is known to begreatest at the base of the crypt, and decreases up the crypt axis [8]. Onepossible way to model this phenomenon in 2D would be to consider a cylindricalgeometry, with a section near the base which is allowed to proliferate and cellsloughing near the top, similar to those presented in [13, 19]. This geometry18ould then determine how cell differentiation evolves as cells migrate towardsthe top of the crypt There is also a known cross-talk between the Notch and theWnt pathways within the cell [9]. Including these dynamics would further allowus to obtain full, realistic model of cellular signalling. Finally, the above modelnaturally resides on a 2D tissue. However, the development and dynamics ofreal biological systems is rarely captured with 2D projections, but rather resideas 2D surfaces in 3D space. Thus the development of a 3D model which iscapable of describing realistic tissue deformations is needed [6]. Acknowledgements
This research was supported by an Australian Government Research Train-ing Program (RTP) Scholarship (awarded to DPJG).
Supplementary InformationSI Movie 1.
Video of simulation from Figures 2a and 2b, which shows howhomogeneous Notch initial conditions lead to a homogeneous tissue.
SI Movie 1.mp4
SI Movie 2.
Video of simulation from Figures 2c and 2d, which shows how aperturbation away from homogeneous Notch initial conditions lead to apatterned tissue. This also shows how the patterning radially propagatesfrom the perturbed cell.
SI Movie 2.mp4
SI Movie 3.
Video of simulation from Figures 2e and 2f, which shows randominitial Notch conditions may lead to a patterned tissue.
SI Movie 3.mp4
SI Movie 4.
Video of simulation from Figures 2e and 2f, which shows randominitial Notch conditions may lead to only a partially a patterned tissue.This shows how Notch patterning is sensitive to initial conditions.
SI Movie 4.mp4
SI Movie 5.
Video of simulation from Figures 2e and 2f, which shows howthe tissue size effects the final patterned tissue. Here, the tissue sizeis not suitable, resulting in a mismatch in patterning as the patterningpropagation meets up with itself, due to periodicity.
SI Movie 5.mp4 I Movie 6.
Video of simulation from Figures 5a – 5d, which shows a dynamictissue evolving over time, with a cell turnover rate of γ = 0 .
1, from times t = 0hrs to t = 20hrs. Here, we can see that the tissue reaches mostlypatterned state. SI Movie 6.mp4
SI Movie 7.
Video of simulation from Figures 5f – 5i, which shows a dynamictissue evolving over time, with a cell turnover rate of γ = 1, from times t = 0hrs to t = 20hrs. Here, we can see that the tissue reaches a partiallypatterned state, as the cell do not sufficiently inhibit their neighbours. SI Movie 7.mp4
References [1] D. H. Alpers, A. N. Kalloo, N. Kaplowitz, C. Owyang, and D. W. Powell.
Textbook of gastroenterology . John Wiley & Sons, 2011.[2] S. J. Bray. Notch signalling: a simple pathway becomes complex.
Naturereviews Molecular cell biology , 7(9):678, 2006.[3] P. Buske, J. Galle, N. Barker, G. Aust, H. Clevers, and M. Loeffler. A com-prehensive model of the spatio-temporal stem cell and tissue organisationin the intestinal crypt.
PLoS computational biology , 7(1):e1001045, 2011.[4] J. R. Collier, N. A. Monk, P. K. Maini, and J. H. Lewis. Pattern formationby lateral inhibition with feedback: a mathematical model of delta-notch in-tercellular signalling.
Journal of theoretical Biology , 183(4):429–446, 1996.[5] J. C. Dallon and H. G. Othmer. How cellular movement determines thecollective force generated by the dictyostelium discoideum slug.
Journal oftheoretical biology , 231(2):203–222, 2004.[6] S.-J. Dunn, I. S. N¨athke, and J. M. Osborne. Computational models reveala passive mechanism for cell migration in the crypt.
PLoS One , 8(11), 2013.[7] A. G. Fletcher, J. M. Osborne, P. K. Maini, and D. J. Gavaghan. Im-plementing vertex dynamics models of cell populations in biology within aconsistent computational framework.
Progress in biophysics and molecularbiology , 113(2):299–326, 2013. 208] A. Gregorieff and H. Clevers. Wnt signaling in the intestinal epithelium:from endoderm to cancer.
Genes & development , 19(8):877–890, 2005.[9] S. K. Kay, H. A. Harrington, S. Shepherd, K. Brennan, T. Dale, J. M.Osborne, D. J. Gavaghan, and H. M. Byrne. The role of the hes1 crosstalkhub in notch-wnt interactions of the intestinal crypt.
PLoS computationalbiology , 13(2):e1005400, 2017.[10] F. J. A. Khalek, G. I. Gallicano, and L. Mishra. Colon cancer stem cells.
Gastrointestinal cancer research: GCR , (Suppl 1):S16, 2010.[11] C. Livingston and M. A. of America.
Knot Theory . Number v. 24 in Carusmathematical monographs. Mathematical Association of America, 1993.ISBN 9780883850275.[12] F. A. Meineke, C. S. Potten, and M. Loeffler. Cell migration and organi-zation in the intestinal crypt using a lattice-free model.
Cell proliferation ,34(4):253–266, 2001.[13] J. M. Osborne, A. Walter, S. Kershaw, G. Mirams, A. Fletcher, P. Path-manathan, D. Gavaghan, O. Jensen, P. Maini, and H. Byrne. A hybridapproach to multi-scale modelling of cancer.
Philosophical Transactions ofthe Royal Society of London A: Mathematical, Physical and EngineeringSciences , 368(1930):5013–5028, 2010.[14] J. M. Osborne, A. G. Fletcher, J. M. Pitt-Francis, P. K. Maini, and D. J.Gavaghan. Comparing individual-based approaches to modelling the self-organization of multicellular tissues.
PLoS computational biology , 13(2):e1005387, 2017.[15] P. Salama and C. Platell. Colorectal cancer stem cells.
ANZ journal ofsurgery , 79(10):697–702, 2009.[16] M. Shanmugathasan and S. Jothy. Apoptosis, anoikis and their relevance tothe pathobiology of colon cancer.
Pathology international , 50(4):273–279,2000.[17] O. Shaya, U. Binshtok, M. Hersch, D. Rivkin, S. Weinreb, L. Amir-Zilberstein, B. Khamaisi, O. Oppenheim, R. A. Desai, R. J. Goodyear,et al. Cell-cell contact area affects notch signaling and notch-dependentpatterning.
Developmental cell , 40(5):505–511, 2017.2118] D. W. Thompson. On growth and form.
On growth and form. , 1942.[19] I. M. M. Van Leeuwen, G. R. Mirams, A. Walter, A. Fletcher, P. Mur-ray, J. Osborne, S. Varma, S. J. Young, J. Cooper, B. Doyle, J. Pitt-Francis, L. Momtahan, P. Pathmanathan, J. P. Whiteley, S. J. Chap-man, D. J. Gavaghan, O. E. Jensen, J. R. King, P. K. Maini, S. L.Waters, and H. M. Byrne. An integrative computational model for in-testinal tissue renewal.
Cell Proliferation , 42(5):617–636, 2009. doi:10.1111/j.1365-2184.2009.00627.x. URL https://onlinelibrary.wiley.com/doi/abs/10.1111/j.1365-2184.2009.00627.xhttps://onlinelibrary.wiley.com/doi/abs/10.1111/j.1365-2184.2009.00627.x