Modeling and analysis of the DeltaNotch dependent boundary formation in the Drosophila large intestine
 Fei Liu†^{1, 2}Email author,
 Deshun Sun†^{1},
 Ryutaro Murakami^{3} and
 Hiroshi Matsuno^{3}
https://doi.org/10.1186/s1291801704558
© The Author(s) 2017
Published: 21 September 2017
Abstract
Background
The boundary formation in the Drosophila large intestine is widely studied as an important biological problem. It has been shown that the DeltaNotch signaling pathway plays an essential role in the formation of boundary cells.
Results
In this paper, we propose a mathematical model for the DeltaNotch dependent boundary formation in the Drosophila large intestine in order to better interpret related experimental findings of this biological phenomenon. To achieve this, we not only perform stability analysis on the model from a theoretical point of view, but also perform numerical simulations to analyze the model with and without noises, the phenotype change with the change of Delta or Notch expression, and the perturbation influences of binding and inhibition parameters on the boundary formation.
Conclusions
By doing all these work, we can assure that our model can better interpret the biological findings related to the boundary formation in the Drosophila large intestine.
Keywords
Background
The large intestine of Drosophila embryo is a middle and large region of the Drosophila hindgut and is subdivided into dorsal and ventral domains, between which a onecellwide strand of boundary cells forms bilaterally for wildtype embryos [1–4]. The large intestine is a multicellular system that involves a number of cells composed of three cell types, dorsal, ventral and boundary cells, organized in a singlelayered epithelial tube [5]. For such developmental patterning problems, different kinds of computational strategies [6] have been proposed, e.g., signaling gradients [7] and activatorinhibitor systems [8], and many computational techniques are adopted, e.g., ordinary/partial differential equation [9] and colored Petri nets [10]. See [11] for a review. However, for the boundary formation in the Drosophila large intestine, the mechanism has been widely explored in vivo, e.g., [1–4, 12], but rarely from the computational point of view, e.g., [13].
It has been shown that the DeltaNotch signaling pathway plays an essential role in the boundary formation of the Drosophila large intestine [1, 12]. In fact, the DeltaNotch pathway is considered as one of the six major signaling pathways in cells, which is active in developing embryos at different phases [14]. Both Notch and Delta proteins are transmembrane proteins, where Notch proteins act as receptors and Delta proteins as ligands. When Delta ligands in a cell bind to Notch receptors in neighboring cells, all the cells in a system may evolve and finally form different types of patterns [15–18].
In order to understand the mechanism of the DeltaNotch pathway, there have been some mathematical models proposed. For example, Collier et al. proposed a simple ordinary differential equation (ODE) model for the DeltaNotch signaling pathway, and discussed numerical simulation of multiplecell systems [19]. Boareto et al. devised a theoretical framework that includes a couple of ODEs to explore the effects of Jagged in cellfate determination [20]. Sprinzak et al. gave a model of a set of ODEs for describing mutual inactivation of Notch and Delta and used this model to illustrate how cisinteractions between Notch and Delta generate mutually exclusive signaling states [21]. Specifically, Matsuno et al. analyzed the mechanism of the Notchdependent boundary formation in the Drosophila large intestine and built a hybrid Petri net model to numerically explore how the DeltaNotch pathway affects the boundary formation in twodimensional space [13].
In this paper, we propose a mathematical model for the DeltaNotch dependent boundary formation in the Drosophila large intestine based on the work of [13], aiming at better interpreting related biological findings and further making predictions. Compared with the existing work in this area, our work has the following main contributions.
(1) We give a mathematical model for the DeltaNotch dependent boundary formation in the Drosophila large intestine, which can better interpret relevant biological findings, so far obtained in the lab. Moreover, this model can be numerically simulated efficiently. To make the model both mathematically and biologically sound, we do the following analysis work.
(2) We perform local stability analysis on the twocell model to make the model mathematically sound. The analysis confirms that the model would reach an equilibrium, which corresponds to that the system would result in a stable pattern in either normal or mutant conditions. The model is also stable when some conditions are satisfied, which means even if there are small disturbances of parameters (e.g., small environmental noises), the system would still converge to a stable state after a period of time. We use the local stability analysis result to determine appropriate parameter values that make the system stable.
(3) We perform numerical simulations with and without noises of Notch expression, which shows that both deterministic and random simulation results are consistent with experimental observations. We further analyze how phenotypes change with the change of Delta or Notch expression.
(4) We analyze the perturbation influences of binding and inhibition parameters on the boundary formation. As the boundary formation is affected by many environmental factors or noises, to make the model more realistic, we need to consider noises into the model by finding those parameters that have significant influences on the dynamics of the model when they vary. According to the analysis result, we add appropriate random noise items to key parameters of the model.
This paper is structured as follows. In the section of methods, we introduce relevant biological background of the DeltaNotch dependent boundary formation in the Drosophila large intestine, and present a mathematical model for this biological phenomenon. In the section of results and discussions, we give stability analysis, simulation analysis, and parameter perturbation analysis of the model. Finally, the conclusion is given.
Methods
In this section, we introduce relevant biological background and give a description of the mathematical model we develop, as well as the simulation and analysis methods which we use.
Boundary formation in the Drosophila large intestine

Scenario 1. In a wildtype embryo, a onecellwide strand of boundary cells forms at the interface of dorsal and ventral domains. See Fig. 3a for an illustration.

Scenario 2. In an embryo where overexpression of Notch or Delta proteins happens, the number of boundary cells would change. If we fix overexpression of Notch proteins, the number of boundary cells would increase with the decrease of Delta background. See Fig. 3 b1 to b2 for an illustration.

Scenario 3. However, if we fix the background of Delta proteins, the number of boundary cells would increase with the increase of overexpression of Notch proteins. See Fig. 3 c1 to c2 for an illustration.

Scenario 4. During the boundary formation in the Drosophila large intestine, environmental factors such as temperature may heavily affect the boundary formation, thus resulting in different (random) phenotypes. Therefore, the model to be built should incorporate such random noises into some parameters of interest.
In what follows, we will construct and analyze our model step by step.
Mathematical model
Here, 1≤i≤N C. D _{ i }, N _{ i } and A _{ i } denote the concentration of Delta proteins, inactive and active Notch proteins, respectively, in the i ^{ t h } cell. In \(\frac {\lambda }{1+\Delta \cdot {{A}_{i}}}\), λ represents the production rate of Delta proteins, and Δ the inhibition coefficient of active Notch proteins A _{ i }. d _{1} represents the degradation rate of Delta proteins. \(\sum \limits _{NG(i)}{{{D}_{i}}}\) represents the concentration of the Delta proteins (in the i ^{ t h } cell) that bind to the Notch proteins of the contacting cells of the i ^{ t h } cell, and N G(i) denotes all the neighbors of the i ^{ t h } cell. f _{1} represents the binding rate of Delta ligands to Notch receptors. λ _{ N } and d _{2} represent the production rate and the degradation rate of inactive Notch proteins, respectively. f _{2} represents the binding rate of Delta ligands to Notch receptors. In \(\frac {a{{N}_{i}}}{b{{D}_{i}}+{{N}_{i}}}\), a represents the conversion rate of Notch proteins from the inactive state to the active state, while b describes the inhibition effect of Delta on Notch in the same cell i. d _{3} is the degradation rate of active Notch proteins.
Simulation and analysis methods
We encode our mathematical model and perform all simulations with Matlab 2016. Besides, we employ the following analysis methods for our model.
Local stability analysis
We conduct local stability analysis to explore stability conditions of our model, with which we further determine values of parameters.
Assume x ^{∗} is an equilibrium point, i.e., f(x ^{∗})=0. J ^{∗} is the Jacobian matrix evaluated at x ^{∗}. The equilibrium point x ^{∗} is stable if all the eigenvalues of the characteristic equation of J ^{∗} have negative real parts [22].
Sensitivity analysis
Results and discussion
In what follows, we present analysis and simulation results for the model.
Local stability analysis
In the following, we will take λ (the production rate of Delta proteins) as an example to study the equilibria and stability of System (2).
Please note that the equilibria of System (2) we obtain above mean that the system will converge to a stable state after a period of time.
It is obvious that all of the six eigenvalues of Eq. (4) are negative, which gives the following conclusion.
Lemma 1
The equilibrium E _{0} is locally asymptotically stable.
Namely, M _{1} M _{2}−M _{3}>0, and M _{1} M _{2}−M ^{′} _{3}>0. Because M _{3}>M ^{′} _{3}, we only need the condition M _{1} M _{2}−M _{3}>0 to be satisfied. That is, if \(\frac {{{M}_{3}}}{{{M}_{1}}{{M}_{2}}}<1\), the eigenvalues of F _{1}+F _{2}=0 and F _{1}−F _{2}=0 are negative; thus, E _{1} is locally stable.
Lemma 2
If \(\frac {{{M}_{3}}}{{{M}_{1}}{{M}_{2}}}<1\), the equilibrium E _{1} is locally asymptotically stable; if\(\frac {{{M}_{3}}}{{{M}_{1}}{{M}_{2}}}>1\), it is unstable.
Please note that the equilibria and stability analysis above show that the model would reach an equilibrium in either case of λ, which corresponds to that the system would result in a stable pattern in either normal or mutant conditions. If there are small disturbances of parameters (e.g., small environmental noises), the system would always converge to a stable state after a period of time when λ=0 according to Lemma 1, but the system would converge to a stable state for λ>0 when the stability conditions in Lemma 2 are satisfied. Besides, we will also use Lemma 2 to carefully choose appropriate parameter values for preserving the stability of the system. This is expected to be shown in the the boundary formation model of the Drosophila large intestine.
Deterministic simulation results
However, it is impossible to analyze the model above with a large number of cells in theory, e.g., a model with 60 cells would result in 180 dimensional ordinary differential equations. In this section, we will explore the model using numerical simulation.
Parameter values of the model
Parameter  Value 

d _{1}  0.01 
d _{2}  0.01 
d _{3}  0.01 
f _{1}  0.665 
f _{2}  0.665 
a  0.012 
b  69 
Δ  10^{6} 
λ  [0, 12] 
λ _{ N }  [0. 0.1] 
Besides, we determine a boundary cell using the following rule. If the concentration of the active Notch proteins (A) in a cell is equal to or greater than 1.0, we regard the cell as a boundary cell.
Phenotype change due to the change of Delta or Notch expression
During the boundary formation in the Drosophila large intestine, different Delta or Notch expression (normal, misexpression or overexpression) may result in different boundary cell distributions (phenotypes) [1]. In our model, we try to incorporate all these phenotypes reported so far. As shown in System (1), Delta or Notch expression is controlled by the rates λ and λ _{ N }, respectively. Therefore, we need to map different Delta and Notch expression observed in experiments to different values of λ and λ _{ N } in the model.
To do this, we explore the phenotype change due to the change of Delta or Notch expression in their parameter space and then determine what parameter value results in what kind of phenotype. As well, we cannot adopt a theoretical analysis due to the large number of equations involved, and instead, we still use numerical simulation. Besides, each time we tune only one parameter, either λ or λ _{ N }, by fixing the other.
Phenotype change due to the change of Delta expression
Phenotype change due to the change of Delta expression
Expression  Phenotype 

12≥λ>0.2  At the first three rows are boundary cells 
0≤λ≤0.2  At all five rows are boundary cells 
From the table above, we can see when λ>0.2, the first three rows are stable boundary cells, while when λ≤0.2, all the five rows are boundary cells. This means there are only two phenotypes appearing due to the change of Delta expression.
Phenotype change due to the change of Notch expression
Phenotype change due to the change of Notch expression
Expression  Phenotype 

0≤λ _{ N }<0.0001  Only at the third row are boundary cells 
0.0001≤λ _{ N }<0.0098  At the first three rows are boundary cells 
0.0098≤λ _{ N }<0.01  At the first four rows are boundary cells 
0.01≤λ _{ N }≤0.1  At all five rows are boundary cells 
Table 3 shows with the increase of λ _{ N }, the number of boundary cells increases, which results in a couple of different phenotypes. Besides, we should notice that a tiny parameter change would result in a big change of the number of boundary cells, so we can deduce that parameter λ _{ N } is very sensitive.
According to the analysis result above, we finally map different Delta and Notch expression observed in experiments to different values of λ and λ _{ N } in the model. Besides, we can reproduce all known biological phenotypes we have found by combining different values of λ and λ _{ N }.
Simulation results with noises of Notch expression
As described in [1], the boundary formation of the Drosophila large intestine is greatly affected by temperature; the random effect of temperature may cause fluctuations of gene expression levels, resulting in different phenotypes, in which boundary cells may randomly be distributed in the dorsal and ventral domains. In this section, we will discuss how to use our model to obtain these random phenotypes.
where e is a small noise.
Perturbation influences of binding and inhibition parameters on boundary formation
During the development of the Drosophila large intestine, there are many environmental factors or noises, which may vary and even heavily affect the boundary formation. To make the model more realistic, we need to consider noises into the model, e.g., adding a random term to each parameter of interest, and analyze the model in noisy conditions. In the section above, we have discussed the random effects of Notch expression on the boundary formation. But in this section, we will further analyze the perturbation influences of the binding and inhibition parameters of the DeltaNotch pathway on the boundary formation. That is, we will explore the inhibition parameters, a and b, and binding parameters, f _{1} and f _{2}.
For these parameters, which ones are of interest? To answer this, we may alternatively find those parameters which have significant influences on the dynamics of the model when they vary. For this, we use sensitivity analysis on the mathematical model to study the influence of parameter perturbation [28]. After that, we can design experiments by considering noises in vivo, and carefully add appropriate noise items to key parameters. We will explore this in the case of overexpression of Notch proteins. For the other cases, we can do the same exploration.
Next, with the parameter setting given in Fig. 8b (of course we can use any other setting), we will explore the influence of parameter perturbation on the boundary formation by changing parameters, a, b, f _{1} and f _{2}, to different perturbation intensities.
Perturbation influence of parameter a
From Fig. 10, we can see that a small perturbation of parameter a may cause a big change of the number of boundary cells (mean), i.e., the bigger the perturbation intensity is, the less the number of boundary cells is (from 36 boundary cells without perturbation to about 15 boundary cells with a perturbation intensity ε=10^{−4}). Furthermore, from the error bars, we know that the bigger the perturbation intensity is, the bigger the standard error is, although the smaller the mean is.
Perturbation influence of parameter b
Figure 11 shows that a big perturbation causes a small change (e.g., there are 36 boundary cells without perturbation and around 34 boundary cells with the perturbation intensity ε=0.2), which is quite different from the effect of parameter a. Therefore, the perturbation influence of parameter b is much smaller than that of a.
Perturbation influence of parameters f _{1} and f _{2}
From Fig. 12, we can see, similar to parameter b, the perturbation of parameter f does not have an obvious influence on the number of boundary cells in the given perturbation intensities.
Conclusion
In this paper, we give a mathematical model for the DeltaNotch dependent boundary formation in the Drosophila large intestine. We aim to use this model to better interpret related biological phenomena and therefore we perform not only theoretical but also simulation analysis to achieve the goal. By combining different analysis techniques, we finally confirm that the model is both mathematically and biologically sound and is sufficient to interpret related experimental findings of the boundary formation in the Drosophila large intestine. Moreover, by modulating parameters, our simulation can generate various abberant patterns of boundary cells that have not been reported in biological studies so far. Biologically, these varieties of phenotypes are assumed to be results of variation of gene products in the DeltaNotch pathway. To further verify the validity of our simulation method, we are now trying to develop a biological experimental system to manipulate levels of forcedgene expression. In a next step, we will shift the purpose of the model from interpretation to prediction. That is, we will use this model to make different predications in different conditions and then help biologists to have an idea of detailed mechanisms of complicated biological systems, and to design experiments to validate the predictions obtained by our simulation.
Declarations
Acknowledgment
We would like to thank Prof. Monika Heiner for fruitful discussions and good comments. We also would like to thank the anonymous reviewers for their constructive comments.
Funding
This work was mainly supported by National Natural Science Foundation of China (61273226), which also financed the articleprocessing charge. FL was supported by JSPS Invitation Fellowships for Research in Japan (L15542) for his stay in Japan during this research. None of the funding agencies had any role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.
Availability of data and materials
The data supporting the results of this article are included within the article.
About this supplement
This article has been published as part of BMC Systems Biology Volume 11 Supplement 4, 2017: Selected papers from the 10th International Conference on Systems Biology (ISB 2016). The full contents of the supplement are available online at https://bmcsystbiol.biomedcentral.com/articles/supplements/ volume11supplement4.
Authors’ contributions
FL conceived the study and developed the model. DS performed the simulations and analyzed the model. HM and RM provided the biological knowledge. All authors read and approved the final manuscript.
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Authors’ Affiliations
References
 Takashima S, Yoshimori H, Yamasaki N, Matsuno K, Murakami R. Cellfate choice and boundary formation by combined action of Notch and engrailed in the Drosophila hindgut. Dev Genes Evol. 2002; 212:534–41.View ArticlePubMedGoogle Scholar
 Murakami R, Shiotsuki Y. Ultrastructure of the hindgut of Drosophila larva, with special reference to the domains identified by specific gene expression patterns. J Morphol. 2001; 248:144–50.View ArticlePubMedGoogle Scholar
 Takashima S, Murakami R. Regulation of pattern formation in the Drosophila hindgut by wg, hh, dpp, and en. Mech Dev. 2001; 101:79–90.View ArticlePubMedGoogle Scholar
 Hamaguchi T, Takashima S, Okamoto A, Imaoka M, Okumura T, Murakami R. Dorsoventral patterning of the Drosophila hindgut is determined by interaction of genes under the control of two independent gene regulatory systems, the dorsal and terminal systems. Mechanisms of Development. 2012; 29:236–43.View ArticleGoogle Scholar
 Kumichel A, Knust E. Apical localisation of crumbs in the boundary cells of the Drosophila hindgut is independent of its canonical interaction partner stardust. PLoS ONE. 2014; 9(5):94038.View ArticleGoogle Scholar
 Wang Y, Zhang SX, Chen L. Modelling biological systems from molecules to dynamical networks. BMC Systems Biology. 2012; 6(Suppl 1):1.View ArticleGoogle Scholar
 Grimm O, Coppey M, Wieschaus E. Modelling the bicoid gradient. Development. 2010; 137:2253–64.View ArticlePubMedPubMed CentralGoogle Scholar
 Murray JD. Mathematical Biology. Berlin: Springer; 2003.Google Scholar
 Nakamasu A, Takahashi G, Kanbe A, Kondo S. Interactions between zebrafish pigment cells responsible for the generation of Turing patterns. PNAS. 2009; 106(21):8029–34.View ArticleGoogle Scholar
 Liu F, Blätke M, Heiner M, Yang M. Modelling and simulating reactiondiffusion systems using coloured Petri nets. Comput Biol Med. 2014; 53:297–308.View ArticlePubMedGoogle Scholar
 Morelli LG, Uriu K, Ares S, Oates AC. Computational approaches to developmental patterning. Science. 2012; 36:187–91.View ArticleGoogle Scholar
 Iwaki D, Lengyel JA. A DeltaNotch signaling border regulated by engrailed/invected repression specifies boundary cells in the Drosophila hindgut. Mech Dev. 2002; 114:71–84.View ArticlePubMedGoogle Scholar
 Matsuno H, Murakami R, Yamane R, Yamasaki N, Fujita S, Yoshimori H, Miyano S. Boundary formation by notch signaling in Drosophila multicellular systems: experimental observations and gene network modeling by genomic object net. In: Proceedings of Pac Symp Biocomput: 37 January 2003; Kauai, Hawaii. World Scientific Press: 2003. p. 152–63.Google Scholar
 Hayward P, Kalmar T, Arias AM. Wnt/notch signalling and information processing during development. Development. 2008; 135:411.View ArticlePubMedGoogle Scholar
 Muskavitch MA. DeltaNotch signaling and Drosophila cell fate choice. Dev Biol. 1994; 166(2):415–30.View ArticlePubMedGoogle Scholar
 Bray SJ. Notch signalling: a simple pathway becomes complex. Nat Rev Mol Cell Biol. 2006; 7:678–89.View ArticlePubMedGoogle Scholar
 Kopan R, Ilagan MXG. The canonical Notch signaling pathway: unfolding the activation mechanism. Cell. 2009; 137:216–33.View ArticlePubMedPubMed CentralGoogle Scholar
 Baker NE. Notch signaling in the nervous system: pieces still missing from the puzzle. BioEssays. 2000; 22:264–73.View ArticlePubMedGoogle Scholar
 Collier JE, Monk NA, Maini PK, Lewis JH. Pattern formation by lateral inhibition with feedback: a mathematical model of DeltaNotch intercellular signaling. J Theor Biol. 1996; 183(4):429–46.View ArticlePubMedGoogle Scholar
 Boareto M, Jolly MK, Lu M, Onuchic JN, Clementi C, BenJacob E. JaggedDelta asymmetry in Notch signaling can give rise to a sender/receiver hybrid phenotype. PNAS. 2015; 22(5):402–9.View ArticleGoogle Scholar
 Sprinzak D, Lakhanpal A, LeBon L, Santat LA, Fontes M, Anderson GA, GarciaOjalvo J, Elowitz MB. Cisinteractions between Notch and Delta generate mutually exclusive signalling states. Nature. 2010; 465:86–90.View ArticlePubMedPubMed CentralGoogle Scholar
 Boyce WE. Elementary Differential Equations and Boundary Value Problems, Binder Ready Version, (10th Edition). USA: Wiley; 2012.Google Scholar
 Hamby DM. A review of techniques for parameter sensitivity analysis of environmental models. Environ Monit Assess. 1994; 32(2):135–54.View ArticlePubMedGoogle Scholar
 Turányi T. Sensitivity analysis of complex kinetic systems: Tools and applications. J Math Chem. 1990; 5(3):203–48.View ArticleGoogle Scholar
 Seydel R. Practical Bifurcation and Stability Analysis. New York: Springer; 2010.View ArticleGoogle Scholar
 Lakshmikantham V, Leela S, Martynyuk AA. Stability Analysis of Nonlinear Systems. Switzerland: Springer; 2015.View ArticleGoogle Scholar
 Khanh NH, Huy NB. Stability analysis of a computer virus propagation model with antidote in vulnerable system. Acta Math Sci. 2016; 36B(1):49–61.View ArticleGoogle Scholar
 Shivamoggi BK. Perturbation Methods for Differential Equations. New York: Springer; 2003.View ArticleGoogle Scholar