Abstract
Many genes and their regulatory relationships are involved in developmental phenomena. However, by chemical information alone, we cannot fully understand changing organ morphologies through tissue growth because deformation and growth of the organ are essentially mechanical processes. Here, we develop a mathematical model to describe the change of organ morphologies through cell proliferation. Our basic idea is that the proper specification of localized volume source (e.g., cell proliferation) is able to guide organ morphogenesis, and that the specification is given by chemical gradients. We call this idea “growth-based morphogenesis.” We find that this morphogenetic mechanism works if the tissue is elastic for small deformation and plastic for large deformation. To illustrate our concept, we study the development of vertebrate limb buds, in which a limb bud protrudes from a flat lateral plate and extends distally in a self-organized manner. We show how the proportion of limb bud shape depends on different parameters and also show the conditions needed for normal morphogenesis, which can explain abnormal morphology of some mutants. We believe that the ideas shown in the present paper are useful for the morphogenesis of other organs.
Similar content being viewed by others
Avoid common mistakes on your manuscript.
1 1. Introduction
Recent advances in molecular biology have provided much information on developmental phenomena. For each developmental event, many responsible genes have been identified, and regulatory networks of these genes and spatio-temporal patterns of their expressions have been clarified. On the other hand, developmental dynamics include spatial elements as well. Each cell has to recognize its own position within a given organ and behave appropriately. Concepts such as positional information and morphogen have been proposed as a solution to this issue (Wolpert, 1969; Wolpert et al., 1998; Jaeger and Reinitz, 2006; Morishita and Iwasa, 2008). Typically, the concentration gradients of morphogens (diffusive chemicals) synthesized at the boundary of an organ provides positional information to the cells spatially distributed. Each cell decodes its own positional value by sensing and processing the morphogen concentrations with intra- and intercellular chemical reactions including complex feedback loops and lateral inhibition (Gierer and Meinhardt, 1972; Meinhardt and Gierer, 2000). Based on its positional value, each cell can respond in diverse manners, such as by differentiation, division, apoptosis, and movement (see Fig. 1).
However, developmental process does not end with spatially heterogeneous cellular responses. These lead to the deformation and change in the boundary of an organ, which may also alter the position of the morphogen source. Thus, the positional information may change in accordance with the deformation and the growth of the organ. Hence, development is a dynamical cycle of interaction among morphogen synthesis, decoding of positional value, and mechanical deformation of organ. With chemical information alone, we cannot fully understand this cycle.
Previous theoretical studies of development with mechanical processes include those with discrete descriptions of tissues and those with continuous description. In the cellular Potts models, each cell is represented as a cluster of grid points with a volume constraint. Cell surface energy depends on cell types, and cell aggregation and sorting can be realized (Graner and Glazier, 1992; Glazier and Graner, 1993; Mombach et al., 1995, 2001; Savill and Hogeweg, 1997; Jiang et al., 1998; Mombach, 1999; Scott et al., 1999; Hogeweg, 2000; Maree and Hogeweg, 2001; Kesmir and De Boer, 2003; Savill and Sherratt, 2003; Zajac et al., 2003; Turner et al., 2004; Zeng et al., 2004; Merks and Glazier, 2005; Merks et al., 2006). In a vertex dynamics model, a polygon formed by linking several vertices represents a cell or a cluster of cells. Each vertex is influenced by forces acting on it and constraints such as minimization of cell boundaries. The vertex dynamics model has been adopted for morphogenesis in Fundulus epiboly and notochord development in Xenopus laevis (Stein and Gordon, 1982; Weliky and Oster, 1990; Weliky et al., 1991; Brodland and Chen, 2000a, 2000b; Chen and Brodland, 2000; Nagai and Honda, 2001; Brodland, 2002; Honda et al., 2004), as well as for plant morphogenesis (Rudge and Haseloff, 2005). In the center dynamics model, a node represents a cluster of cells and receives forces from its neighboring nodes, which are linked based on actual spatial configuration. The center dynamics model has been adopted to describe cell aggregation, locomotion, and rearrangement (Honda, 1978; Honda et al., 1979, 1984; Graner and Sawada, 1993; Palsson and Othmer, 2000; Palsson, 2001; Meineke et al., 2001). Drasdo et al. dealt with one-layer growing tissues and explained spatial patterns such as cleavage, blastulation, and gastrulation (Drasdo et al., 1995; Drasdo and Forgacs, 2000; Drasdo, 2000; Drasdo and Loeffler, 2001). The immersed Boundary Method (IBM) was introduced by Peskin (1972) to describe flow patterns around heart valves. It consists of the interaction between incompressible viscous fluid and elastic boundary structures immersed in the fluid. IBM has been adopted to flagella motion and swimming of microorganisms (Peskin, 1972; Fauci and Peskin, 1988; Fauci, 1990; Dillon et al., 1995, 1996; Bottino, 1998; Bottino and Fauci, 1998; Dillon and Fauci, 2000; Lubkin and Li, 2002), and to the growth of trophoblast bilayer and the outgrowth of the vertebrate limb bud (Dillon and Othmer, 1999; Rejniak et al., 2004).
In this paper, we attempt to answer a fundamental question: what determines organ morphologies that change through tissue growth. Our basic idea is that the proper specification of localized volume source (realized by cell proliferation) is able to guide organ morphogenesis, and that the specification is given by inhomogeneous chemical distributions (e.g., morphogen gradients). We call the idea “growth-based morphogenesis.” We discovered that this morphogenetic mechanism works if the tissue is elastic for small deformation, but plastic for large deformation. That is, when the shape of the tissue is deformed with small magnitude, it tends to go back to the original shape (elasticity), which stabilizes the present organ shapes against perturbations. However, if the deformation is large, the tissue would not go back to the original one and stay in a new shape (plasticity), which allows the organ shapes to change.
We adopt a modeling framework based on the center dynamics model to describe the change of organ morphologies. Our model consists of the mechanical interaction between epithelial and mesenchymal tissues with different physical properties and topologies, which are commonly observed in many developmental processes.
We illustrate the above concept by using an example of vertebrate limb bud formation and elongation processes, the system with a long history of research (Wolpert et al., 1998). In the limb development, the apical ectodermal ridge (AER), a thick epithelium located at the tip is closely related to the elongation. Its removal by surgery leads to apoptosis at the tip of the limb bud, and consequently incomplete outgrowth (Dudley et al., 2002). FGF family is considered to be a candidate that promotes active cell proliferation (or specifies the location of volume source) at the tip. Lu et al. observed in their experiments that the limb bud shape becomes wider when the Fgf expression level is increased by genetic operation and also showed that the mutant has abnormality in the digit number. Verheyden et al. (2005) showed that the mutant without the expression of FGF receptors has a wider limb bud in the antero-posterior direction and loses distal structures such as fingers in later stages. In addition, Niederreither et al. observed a branched limb bud in the mutant that lacks the middle part of the AER. These suggest that abnormalities of limb morphology might be caused by dysfunction of the mechanism specifying the volume source.
Previous theoretical studies on the limb development have focused on skeletal patterns (Newman and Frisch, 1979; Oster et al., 1983, 1985; Murray et al., 1988; Hentschel et al., 2004; Izaguirre et al., 2004; Chaturvedi et al., 2005; Cickovski et al., 2005; Miura et al., 2006). They do not deal with the growth of the limb shape itself. An exception is the classical work by Ede and Law (1969), who modeled limb bud shape by a set of points on the lattice space. The outgrowth was based on cell division that places a daughter cell on the nearest empty lattice point. Oriented division is required to achieve elongated shape of limb bud. However, their model considers no mechanical forces such as the elasticity between cells. Another exception is the work by Dillon and Othmer (1999). They modeled limb bud elongation processes by using IBM. In their model, the width of the limb is constrained by elastic fibers connecting between the anterior and posterior boundaries to maintain the realistic limb shape.
In our model, a limb bud protrudes and elongates from a flat lateral plate with a constant limb width in a self-organized manner. The morphology of the bud is controlled by the spatio-temporal pattern of the volume source, unlike existing models for limb bud formation. The model℉s results can explain abnormal morphology of some mutants (Lu et al., 2006; Verheyden et al., 2005; Niederreither et al., 2002). We investigate the relation between the proportion of limb bud shape and various system parameters. We conclude that the balance of elasticity between different tissues is important to achieve normal morphogenesis with smooth boundaries.
2 2. Model
A vertebrate limb bud consists of two major tissues, mesenchyme, and epithelium (Wolpert et al., 1998). In this study, these two tissues are described by a network of finite number of nodes. Each node corresponds to the location of the center of a mesenchymal (M) cell or an epithelial (E) cell (see Fig. 2). We can also regard a node as a cluster of M (or E) cells instead of a single cell. In the following, we refer to a node representing an M cell as an M-node and an E cell as an E-node.
We start with a set of E and M nodes lined in parallel, corresponding to a flat lateral plate. By cell division, the region filled with M-nodes expands, and the chain of E-nodes surrounding it becomes longer.
Our modeling is two-dimensional whilst the real limb-bud is three dimensional. The main purpose of this study is to illustrate a basic idea that specifying spatio-temporal patterns of active cell proliferation is able to guide organ morphologies, rather than showing a detailed realistic description of limb-bud formation.
2.1 2.1. Generation of tissue network
In this model, the epithelium is regarded as a one-layer sheet. Adjacent epithelial cells form with junctions. Occasionally an E-node divides, and the links between E-nodes are relinked as explained later. This modeling of epithelium is based on that it is often regarded as a one-layer sheet, and that adjacent epithelial cells have tight junctions with each other (Alberts et al., 2002).
In the region enclosed by a chain of E-nodes, are M-nodes. To consider the linking between M-nodes, we first separate the region into sub-triangles each having a node at the vertex by applying Delaunay triangulation to the set of all nodes, and then each edge of triangles including M-nodes are lined (see Fig. 2(1)). We denote the Delaunay triangulation for a set S of nodes in the plane by DT(S). In triangulation DT(S) no node in S is inside the circumcircle of any triangle (Sloan, 1987). The Delaunay triangulation is done at each time step. Since it is mathematically equivalent to the Voronoi partition, we can also describe the tissue as a set of polygonal cells (Fig. 2(2)). For the boundary of the organ, we adopted a specific partition method by Bottino et al. (2002).
Each link has an optimal distance where its potential energy becomes minimal. The nodes have an elastic property—the length of the links between them tends to return to the optimal one when deviated. Thus, when the shape of the organ is deformed by a small magnitude, the shape tends to return to its original one (see Fig. 3).
The shape also has plasticity—if the organ is largely deformed due to cell movement, division, and so on, the shape would not return to the original one. The plastic behavior occurs because internodal links are re-linked based on Delaunay triangulation (see Fig. 3). In the region with high cell proliferation activity, cells move frequently, and the neighboring relation among nodes change. In contrast, in the region where cell proliferation activity is low and cellular movement is slow, internodal links tend to remain unchanged. Plasticity is required for proper growth of an organ’s shape.
2.2 2.2. Internodal potential and cellular movement
According to previous studies (Nagai and Honda, 2001; Meineke et al., 2001; Drasdo, 2000), we assume that inertial forces acting on the cellular network can be neglected and that the displacement of each node within a small time interval Δt occurs in the direction to reduce the potential energy:
, where x i and Φ i are the positional vector and the total potential energy for the node i. μ i is a friction (viscous) coefficient of the node. For simplicity, we assume μ i = μ for all i.
The total potential energy for an M-node i,Φ i , is given as follows:
where Φ MM ij is an internodal potential energy between M-nodes and the first summation in the right-hand side is for all M-M links including the M-node i. Φ ME ij is an energy between M and E nodes and the second summation is for all E-M links including the M-node. Similarly, the total potential energy acting on an E-node i, Φ E i , is given as follows:
, where Φ EE ij is an internodal energy between E-nodes. The first summation in the righthand side is for all E-E links including the E-node i. Φ EM ij (= Φ ME ji is energy between E and M nodes and the second summation is for all E-M links including the E-node. In the following, the energy Φ k ij (k, l ∈ {M,E}) is a function of the distance between two linked nodes. To be specific, it is given in a manner similar to intermolecular potential as follows:
, where r ij = |r ij | ≡ |x j − x i | and Φ kl ij (k, l} ∈ {M,E}) are the distance and a potential energy between two nodes, i and j whose positional vectors are x i and x j , respectively. ɛ kl (k, l ∈ {M,E} determines the gradient of the energy, which is related to the elasticity of tissues. The equilibrium distances between linked nodes are given by \(r_{kl}^0 = \sqrt 3 {\sigma _{kl}}\) (see Fig. 4(a)). Around an equilibrium distance, internodal links work as a linear spring, while a steep gradient for r ij < r 0 kk represents a strong repulsion, and hence limited compressibility of each cell. The elasticity vanishes for very large distance (r ij ≫ r 0 kk ).
The elastic properties differ between epithelium and mesenchyme. Epithelial cells form tight junction with neighboring cells and form an epithelial layer to prevent mesenchymal cells from leaking out. In our model, this can be realized as ɛ EE > ɛ MM . As shown in the results section, the elastic balance between epithelium and mesenchyme is an important factor for normal morphogenesis.
2.3 2.3. Organ growth through cell proliferation
Organ growth is controlled by “net” proliferation rate at each location, defined as the number of cell divisions minus that of cell death “per unit time.” For example, no apparent change in the number of cells is observed when the rates of cell division and apoptosis are the same.
In the vertebrate limb development, the AER (apical ectodermal ridge) is closely related to its elongation. For example, when the AER is removed by surgery, apoptosis occurs at the tip of limb bud (up to about 200 µm from the tip), which leads to incomplete elongation (Dudley et al., 2002). Further, the rate of cell division at the tip is more frequent than that in the proximal region in middle stages or later (e.g., from st. 23 in chick) although the difference is hardly observed at initial stages of limb bud formation (Hornbruch and Wolpert, 1970). In addition, according to experiments about the fate map of cells (Vargesson et al., 1997; Sato et al., 2007), the cell fate at the tip (200∼300 µm from the tip) has a large variation for the future position along the proximo-distal (P-D) axis. These experimental observations imply that the net proliferation rate at the tip of a limb bud (e.g., 200∼300 µm from the tip) is higher than that in the proximal region (at least in middle stages or later).
In this study, we conjecture that the net proliferation rate is controlled by a gradient signal from the tip of a limb bud, which we call an “AER-signal.” FGF family is a candidate of this AER-signal because FGFs are diffusive proteins synthesized at the AER (Martin, 1998; Sun et al., 2000, 2002) and the expression patterns of down stream genes of FGFs are reported to have a gradient along the P-D axis (Pascoal et al., 2007).
It should be noted that, as shown later, we observe that the early phase of the limb bud formation can be performed without a clear gradient of the net proliferation rate along the P-D axis. Instead, the bud can be formed under the condition of a uniform rate of net proliferation. However, a difference of net proliferation between the distal and proximal regions is required for the late phase of the limb bud elongation (see below).
Many chemicals are involved in the limb morphogenesis such as SHH from the Zone of Polarizing Activity (ZPA) (Johnson and Tabin, 1997; Niswander, 2002; Pagan et al., 1996; Niswander et al., 1994; Laufer et al., 1994; McGlinn and Tabin, 2006). However, there are experiments showing that the limb continues to elongate after the removal of the ZPA (Pagan et al., 1996). Since the purpose of this study is to develop a minimum model to describe limb bud formation and elongation processes, we do not include the ZPA into the model.
2.3.1 2.3.1. AER and diffusion of AER-signal
Complex genetic interactions are known to regulate the position and the size of AER, including well-studied AER-ZPA feedback, but a detailed mechanism is not identified completely (Johnson and Tabin, 1997; Niswander, 2002). In the present model, we assume that the position and the size of AER are given in advance: NAER consecutive E-nodes are AER. The AER-signal is synthesized in the AER and transported into the mesenchymal tissue by diffusion.
Let c m i be the AER-signal concentration at M-node i, and c AER k } that at E-node k of AER. The diffusion of the AER-signal occurs only between linked nodes. The dynamics of its concentration is given by the following reaction-diffusion equation on the nodal network:
, where D is a diffusion constant of the AER-signal. For each node, the flux of a chemical is proportional to the difference between the focal node and its neighboring node. The summation in the second term of the right-hand side is for all M-nodes j linked with the node i, and that in the third term is for all E-nodes k linked with the node i. c AER k is assumed to be constant for all time t. The last term indicates the degradation of the AER-signal at each M-node, where γ is a degradation rate constant.
In the following numerical studies, we assume that the AER-signal diffusion occurs much faster than cell proliferation and movement. Thus, the spatial pattern of the AER-signal is always at the quasi-equilibrium state when we discuss the dynamics of limb bud morphogenesis.
2.3.2 2.3.2. Cell division
M-nodes divide stochastically with the probability P div dt within a time interval dt. P div is assumed to be proportional to the concentration of the AER-signal, that is, P div = f div c M i , where f div is the cell division frequency per unit AER-signal concentration.
When an M-node divides, a dividing direction is chosen randomly, and two daughter nodes are placed as follows (see Fig. 4(b)):
, where x, x daughter1, and x daughter2 are the positional vectors of the mother node and the two daughter nodes. e is a unit vector randomly chosen and d div is the distance between the two daughter nodes just after the division. Since we chose the value of d div as 40–50% of internodal equilibrium length, daughter cells do not receive very strong forces. Further, we confirmed that a small change of d div did not lead to a large change of the limb shape.
Organ growth is observed throughout the limb mesenchyme after the AER removal (Dudley et al., 2002). Thus, the mesenchyme may have a low, constant proliferation rate (basal rate) uniformly over the whole region in addition to the proximally-biased rate enhancement due to the AER-signal. We examined the case where a uniform basal rate is added to the AER-dependent proliferation and confirmed that results obtained below are hardly affected if the added basal value is small. In the following analysis, for simplicity, the basal net proliferation rate is not included.
We assume that the cell density in epithelial layer stays near a standard value because the detailed mechanism of the proliferation of epithelial cells is still unknown. Due to the proliferation of M-nodes, epithelial membrane swells and each E-E segment is stretched. A new E-node is added in the middle of an E-E segment when the length of the E-E segment exceeds twice its equilibrium distance.
The procedure for numerical calculations is described in the Appendix.
3 3. Results
3.1 3.1. Growth of the limb bud
Figure 5 shows a typical temporal change in the shape of a growing limb bud produced by the model. The process of limb bud growth consists of two phases. In the early phase, a small bud protrudes from a flat lateral plate and swells until it has a certain width W. Here, the limb bud width is defined as the size measured in the anterior-posterior direction. In the late phase, the limb continues to elongate perpendicularly to the lateral plate with a constant width W. Since we stopped the calculation when the total node number reached a certain value N e (N e = 800 in Fig. 5), the final length of the limb, L, is almost proportional to 1/W.
We can explain the developmental process more intuitively by showing the spatial distribution of AER-signal concentration (c M i ), the magnitude of energy gradient (∥∇Φ i ), and the energy level (ΔΦ i )} in the mesenchyme. Here ΔΦ i is defined as Φ i − Φ 0 i , where Φ 0 i is the basal energy where all internodal links have their equilibrium length. At the tip of the limb bud, the AER-signal concentration is high (see Fig. 6(a2)), which leads to frequent cell divisions. Bright color at the tip in Fig. 6(a3) indicates the area with a high energy level, that is, compressed region. The pressure induced by the net cell proliferation is the driving force for the elongation of limb bud. As shown in Fig. 6(a4), cellular movement is active at the tip of the limb bud (bright color in the figure shows steep energy gradient).
3.2 3.2. What controls the limb bud shape?
There are many parameters for mechano-chemical properties of the tissues. These parameters can be classified into two categories: One category includes factors determining the spatial pattern of the net cell proliferation in the mesenchyme, such as AER size, AER-signal diffusivity, AER-signal expression level, and division frequency of M-nodes. These factors control the aspect ratio of the limb bud shape at the final pattern, defined as, L/W ∝ L 2. In contrast, the second category of parameters includes the internodal potentials determining the elasticity of tissues and cell mobility. Changes in the parameters in this category can cause morphological anomalies, as shown later.
To study parameter dependences of the limb shape, we chose a standard set of parameters giving the shape similar to observed one, and then examined the model with modified parameters. Parameters are changed one by one around the standard value with all the other parameters fixed. The results are summarized as follows:
3.2.1 3.2.1. Aspect ratios of the limb bud
(1) AER size A fixed number of consecutive E-nodes act as AER, the source of the AER-signal. We call the number of these cells the AER size, and denote it by N AER (see Fig. 2). As illustrated in Fig. 7, the limb bud becomes wider as N AER increases. For a larger AER size, the spatial distribution of the AER-signal extends and the area of the net proliferation becomes broader in the anterior-posterior direction (see Fig. 6(b)). When N AER changes, only the limb width changes without morphological anomalies during the transient state (data not shown).
Figure 8 illustrates that the aspect ratio L/W in the final pattern decreases as AER size N AER increases. The effect of the AER size on the aspect ratio was much stronger than the effects of other factors. This suggests that the positioning and maintenance of the AER and the control of its size are very important to realize the accurate proportion of the limb bud.
(2) Diffusivity and expression level of AER-signal The diffusivity of the AER-signal is determined by its diffusion constant D and its degradation rate γ (see Eq. (4)). In a one-dimensional diffusion with linear degradation, with the boundary concentration given by the source strength, the concentration of the chemical at equilibrium decreases with the distance x from the source exponentially: in proportion to \(\exp \left[ { - x/\sqrt {D/\gamma } } \right]\). The quantity \(\sqrt {D/\gamma } \) is an indicator of the diffusivity of the AER-signal. An increase in \(\sqrt {D/\gamma } \) makes the limb bud broader and shorter, which is similar to the change caused by an increase in the AER size (see Figs. 7 and 8).
In the apical region with strong volume source, the deformation of the nearby tissues is large and the neighborhood relationship among cells is easy to change (i.e., the plasticity of tissue). In contrast, in the proximal region with low activity, the configuration remains unchanged because the deformation is small (the elasticity of tissue). As a consequence, the volume expansion occurring on the tip of the limb bud leads to the limb bud elongation, but does not affect the proximal part at all. When the diffusivity \(\sqrt {D/\gamma } \) is small, the volume source is restricted to the limb bud tip. Thus, the limb extends in the distal direction with a slim width. In contrast, when the diffusivity is large, the volume source expands to more proximal region, which makes the limb bud expand more to sideways, i.e., the antero-posterior direction as well.
The shape of the limb bud formed also depends on the expression level of the AERsignal at its source, c AER (see Figs. 7 and 8). For a larger c AER, active proliferation area extends proximally, and the shape becomes shorter and wider.
(3) Division frequency When the frequency of cell division per unit concentration of the AER-signal, f div , becomes higher, the limb bud shape becomes broader and shorter (see Figs. 7 and 8). As shown in Fig. 6(d), the proliferation area with compression expands proximally as f div increases, which leads to a more sideward growth.
3.2.2 3.2.2. Morphological anomalies
(4) Internodal potential The change in the three factors discussed above modified the aspect ratio of the limb bud, but no morphological anomalies were caused. In contrast, if the balance among internodal potential changes, many different abnormal shapes were obtained. For example, Fig. 9 shows the transient and final patterns of limb bud shapes with different ratios of ɛ EE to ɛ MM. Shapes with smaller (larger) ɛ EE are similar to those with larger (smaller) ɛ MM.
When the ratio of ɛ EE to ɛ MM is small, the transient shape becomes a rugged hemicircle and the organ surface of the final pattern becomes also rugged. This may be because the internal pressure due to cellular proliferation freely disperses since the elasticity between E-E links is relatively weak. These bumps do not seem to originate from the internodal equilibrium distance or from the finiteness of the node number, because we observed bumps of a similar size when the internodal equilibrium distance halved and the node number increased.
In contrast, when the ratio of ɛ EE to ɛ MM is large, the transient shape becomes Gaussian-like and the final one is balloon-like, because elasticity in the epithelium is strong.
The change in ɛ EM somewhat affects the smoothness of organ surface, but hardly affects the proportion of the limb bud.
3.3 3.3. The region for cell proliferation
In the above analysis, the division of M-nodes was regulated by the concentration of the AER-signal, and the position of the AER moves as the organ grows. If instead the area of active cell proliferation is fixed, rather than specified by the concentration of the AER-signal, the model can no longer generate normal elongation of limb bud. In Fig. 10(a), mesenchymal cells can proliferate only in a rectangular window beneath the initial epithelial layer, and the size and location of the window is fixed during development. Although a normal-shaped limb bud can be formed from the flat lateral plate in the early phase, the organ shape gradually becomes balloon-like in the late phase of the simulation. This remains valid for different sizes of the window of active cell division (data not shown).
The results suggest that the proliferation area needs to be dynamically specified during organ growth to have an elongated shape of the limb bud. Moreover, the specification is possibly realized by the gradient of chemicals (morphogen) released from the organ boundary that moves with the growth and deformation. In contrast, the spatial gradient of proliferation rate is not necessarily required for the initiation of bud formation.
Figures 10(b) and (c) clearly demonstrate the importance of the position of the source of the AER-signal. We assumed that the proliferation depends on the concentration of the AER-signal, but the position of the AER moves with time (the trajectories of AER-center are indicated by bold white lines in both figures). In Fig. 10(b), the position of AER moves from the tip of limb bud to the posterior side. The tissue growth tracks the trajectory of the AER, which leads to a limb bud that bends toward the posterior direction. Such bending is observed in the actual development of vertebrate limb bud. It may be realized when cell proliferation area is affected by posterior tissue through chemicals such as SHH secreted from ZPA that is another morphogen source locating at posterior mesenchyme. We also performed a numerical experiment in which the AER bifurcates with time. Similarly, the tissue growth tracked the trajectories of the divided AERs, and ended up with a branched limb bud (see Fig. 10(c)). In the experiment by Niederreither et al. (2002), a branched limb bud was observed in the mutant that lacks the middle part of the AER.
4 4. Discussion
We have developed a mechano-chemical model to describe the change of organ morphology. We assumed that the prime force of development is the volume source realized by cell proliferation. Our basic idea is that the proper specification of localized volume source is able to guide organ morphogenesis, where the specification is achieved by chemical gradients (e.g., FGF in the case of limb bud outgrowth). We may call this concept as the principle of “growth-based morphogenesis.” We have found that this morphogenetic mechanism works if the tissue is elastic for small deformation and plastic for large deformation.
To illustrate our idea, we studied the developmental process of vertebrate limb bud. The model is composed of the mechanical interaction between epithelium and mesenchyme with different physical properties, which is realized by adopting different intercellular potentials. The location of active cell proliferation is specified by the gradient of FGF from the AER. We examined the relationship between system parameters and the limb morphology.
The results can be summarized as the following three points:
First, the aspect ratio of the limb morphology, i.e., the ratio of limb bud length to limb bud width, is determined by the spatial pattern of volume source realized by cell proliferation in the mesenchyme. If the proliferation area is limited to the tip of the limb bud, the limb elongates in the distal direction and the aspect ratio (length/width) becomes large. In contrast, if the proliferation area spreads to more proximal part, the shape becomes broader and the aspect ratio becomes smaller. The mechanism can be explained as follows: In the area with active cell proliferation, the deformation of the nearby tissue is large and the configuration among cells is easy to change due to the plasticity of tissue. In contrast, in the area with low activity, the configuration remains unchanged due to the elasticity of tissue. Thus, the spatial pattern of cell proliferation determines the proportion of the limb bud. Factors that decide the spatial pattern include the AER size, the diffusivity of the AER-signal, the expression level of the AER-signal at the source, and the frequency of cell division.
Second, the balance of elasticity between mesenchymal and epithelial tissues is very important for the normal growth of limb bud. When the elasticity of epithelial tissues is strong, it is difficult for mesenchymal cells to push the epithelial layer distally, which makes the limb bud elongation difficult. In contrast, a strong elasticity of mesenchymal tissues leads to bumpy blocks because of the easy deformation of the epithelial layer. We believe that the elastic balance between mesenchyme and epithelium is required for normal morphogenesis with smooth organ boundary not only for the vertebrate limb bud formation, but also for many other developmental processes.
Third, normal elongation of limb bud is not observed any longer if the area with high proliferation activity does not dynamically change with the growth and deformation, although a small bud was formed in all cases. In addition, the shape of limb bud is strongly affected by the spatial pattern of the area for active cell proliferation, which in turn is specified by the trajectory of the morphogen source.
Description of a tissue shape that is both elastic and plastic is a key to the modeling of development with growth. In the model for intracellular processes, Bottino (1998) examined different linking rules to reproduce viscoelastic properties of the actin network within a cell.
In the fluid-dynamics model by Dillon and Othmer (1999), the anterior and posterior boundaries of the limb bud are connected by elastic links. This constraint prevents the area from being balloon-like, and it makes possible the extension of the area in the proximodistal direction. In contrast, in our model, the limb bud shape develops in a manner guided by the location of volume sources.
Genes and their products related to developmental processes affect the mechanochemical properties of cells. If a mutant has abnormalities in the morphology compared to wild types, the result of this paper suggests that the mutant is likely to have abnormalities in genes or their products whose functions are related to the process of specifying the area of cell proliferation and/or to the elasticity of tissues. For example, Lu et al. observed in their experiments that the limb bud shape becomes wider when the Fgf expression level is increased by genetic operations (see Fig. 1 in Lu et al., 2006). This observation is consistent with the prediction of the present paper (see Section 3.2). In addition, the experiment also showed that the mutant has abnormality in the digit number, which clearly suggests that the abnormalities in the position of morphogen source and in the morphogen amount that determine active proliferation area cause critical morphological deficiencies. Another example is the observation by Verheyden et al. (2005), who showed that the mutant without the expression of FGF receptors has a wider limb bud in the antero-posterior direction and loses distal structures such as fingers in later stages. Since FGF receptors are required to specify active proliferation area by FGF signal, their absence may make the organ difficult to elongate as expected by the results in Section 3.3 of the present paper.
In this study, the size and the position of the AER are treated as given. It is known that these are regulated by the interaction between the AER and the ZPA (Niswander et al., 1994; Laufer et al., 1994; McGlinn and Tabin, 2006; Hirashima et al., 2008). The ZPA is another morphogen source located at a posterior mesenchyme in the vertebrate limb development, and is required to maintain the AER activity (Wolpert et al., 1998). The ZPA also plays an important role to give the positional information along anterior-posterior axis. An extended model including the interaction between AER and ZPA will be an important target of future theoretical study. Further, in the present study, the driving force for cell migration is assumed to be the pressure caused by cell proliferation. However, there is an ongoing discussion about what is the driving force to cellular migration. Li and Muneoka (1999), for example, showed that the AER has a chemoattractive function and distal mesenchymal cells move toward the AER. We think that the investigations on the relation between the organ shape and the factors other than pressure due to proliferation are also important future works. In addition, our modeling deals with a two-dimensional projection. The fact that the dorso-ventral distance in a real limb bud is smaller than the antero-posterior distance may present a problem for a three-dimensional extension of this model. We believe that most of qualitative properties observed in this paper can be seen in the three-dimensional model. However, clarifying the difference between 2D and 3D systems is an important open issue.
References
Alberts, B., Johnson, A., Lewis, J., Raff, M., Roberts, K., Walter, P., 2002. Molecular Biology of the Cell, 4th edn. Garland, New York.
Bottino, D.C., 1998. Modeling viscoelastic networks and cell deformation in the context of the immersed boundary method. J. Comput. Phys. 147, 86–13.
Bottino, D.C., Fauci, L.J., 1998. A computational model of ameboid deformation and locomotion. Eur. Biophys. J. 27, 532–39.
Bottino, D., Mogilner, A., Roberts, T., Stewart, M., Oster, G., 2002. How nematode sperm crawl. J. Cell. Sci. 115, 367–84.
Brodland, G.W., 2002. The differential interfacial tension hypothesis (DITH): a comprehensive theory for the self-rearrangement of embryonic cells and tissues. J. Biomech. Eng. 124, 188–97.
Brodland, G.W., Chen, H.H., 2000a. The mechanics of heterotypic cell aggregates: insights from computer simulations. J. Biomech. Eng. 122, 402–07.
Brodland, G.W., Chen, H.H., 2000b. The mechanics of cell sorting and envelopment. J. Biomech. 33, 845–51.
Chaturvedi, R., Huang, C., Kazmierczak, B., Schneider, T., Izaguirre, J.A., Glimm, T., Hentschel, H.G.E., Glazier, J.A., Newman, S.A., Alber, M., 2005. On multiscale approaches to 3-dimensional modeling of morphogenesis. J. R. Soc. Interface 2, 237–53.
Chen, H.H., Brodland, G.W., 2000. Cell-level finite element studies of viscous cells in planar aggregates. J. Biomech. Eng. 122, 394–01.
Cickovski, T., Huang, C., Chaturvedi, R., Glimm, T., Hentschel, H.G.E., Alber, M., Glazier, J.A., Newman, S.A., Izaguirre, J.A., 2005. A framework for three-dimensional simulation of morphogenesis. IEEE/ACM Trans. Comput. Biol. Bioinf. 2, 1545–963.
Dillon, R.H., Fauci, L.J., 2000. An integrative model of internal axoneme mechanics and external fluid dynamics in ciliary beating. J. Theor. Biol. 207, 415–30.
Dillon, R., Othmer, H.G., 1999. A mathematical model for outgrowth and spatial patterning of the vertebrate limb bud. J. Theor. Biol. 197, 295–30.
Dillon, R., Fauci, L., Gaver, D. III, 1995. A microscale model of bacterial swimming, chemotaxis and substrate transport. J. Theor. Biol. 77, 325–40.
Dillon, R., Fauci, L.J., Fogelson, A., Gaver, D. III, 1996. Modeling biofilm processes using the immersed boundary method. J. Comput. Phys. 129, 57–3.
Drasdo, D., 2000. Buckling instabilities of one-layered growing tissues. Phys. Rev. Lett. 84, 4244–247.
Drasdo, D., Forgacs, G., 2000. Modeling the interplay of generic and genetic mechanisms in cleavage, blastulation, and gastrulation. Dev. Dyn. 219, 182–91.
Drasdo, D., Loeffler, M., 2001. Individual-based models to growth and folding in one-layered tissues: intestinal crypts and early development. Nonlinear Anal. 47, 245–56.
Drasdo, D., Kree, R., McCaskill, J.S., 1995. Monte Carlo approach to tissue-cell populations. Phys. Rev. E 52, 6635–657.
Dudley, A.T., Ros, M.A., Tabin, C.J., 2002. A re-examination of proximodistal patterning during vertebrate limb development. Nature 418, 539–44.
Ede, D.A., Law, J.T., 1969. Computer simulation of vertebrate limb morphogenesis. Nature 221, 244–48.
Fauci, L.J., 1990. Interaction of oscillating filaments: a computational study. J. Comput. Phys. 86, 294–13.
Fauci, L.J., Peskin, C.S., 1988. A computational model of aquatic animal locomotion. J. Comput. Phys. 77, 85–08.
Gierer, A., Meinhardt, H., 1972. A theory of biological pattern formation. Kybernetik 12, 30–9.
Glazier, J.A., Graner, F., 1993. Simulation of the differential driven rearrangement of biological cells. Phys. Rev. E 47, 2128–154.
Graner, F., Glazier, J.A., 1992. A two-dimensional extended Potts model for cell-sorting. Phys. Rev. Lett. 69, 2013–016.
Graner, F., Sawada, Y., 1993. Can surface adhesion drive cell rearrangement? Part II: a geometrical model. J. Theor. Biol. 164, 477–06.
Hentschel, H.G.E., Glimm, T., Glazier, J.A., Newman, S.A., 2004. Dynamical mechanisms for skeletal pattern formation in the vertebrate limb. Proc. R. Soc. Lond. B 271, 1713–722.
Hirashima, T., Iwasa, Y., Morishita, Y., 2008. Distance between AER and ZPA is defined by feed-forward loop and is stabilized by their feedback loop in vertebrate limb bud. Bull. Math. Biol. 70, 438–59.
Hogeweg, P., 2000. Evolving mechanism of morphogenesis: on the interplay between differential adhesion and cell differentiation. J. Theor. Biol. 203, 317–33.
Honda, H., 1978. Description of cellular patterns by Dirichlet domains: the two-dimensional case. J. Theor. Biol. 72, 523–43.
Honda, H., Morita, T., Tanabe, A., 1979. Establishment of epidermal cell columns in mammalian skin: computer simulation. J. Theor. Biol. 81, 745–59.
Honda, H., Yamanaka, H., Dan-Sohkawa, M., 1984. A computer simulation of geometrical configurations during cell division. J. Theor. Biol. 106, 423–35.
Honda, H., Tanemura, M., Nagai, T., 2004. A three-dimensional vertex dynamics cell model of space-filling polyhedra simulating cell behavior in a cell aggregate. J. Theor. Biol. 226, 439–53.
Hornbruch, A., Wolpert, L., 1970. Cell division in the early growth and morphogenesis of the chick limb. Nature 226, 764–66.
Izaguirre, J.A., Chaturvedi, R., Huang, C., Cickovski, T., Coffland, J., Thomas, G., Forgacs, G., Alber, M., Hentschel, G., Newman, S.A., Glazier, J.A., 2004. COMPUCELL, a multi-model framework for simulation of morphogenesis. Bioinfomatics 20, 1129–137.
Jaeger, J., Reinitz, J., 2006. On the dynamic nature of positional information. BioEssays 28, 1102–111.
Jiang, Y., Levine, H., Glazier, J.A., 1998. Possible cooperation of differential adhesion and chemotaxis in mound formation of Dictyostelium. Biophys. J. 75, 2615–625.
Johnson, R.L., Tabin, C.J., 1997. Molecular models for vertebrate limb development. Cell 90, 979–90.
Kesmir, C., De Boer, R.J., 2003. A spatial model of germinal center reactions: cellular adhesion based sorting of B cells results in efficient affinity maturation. J. Theor. Biol. 222, 9–2.
Laufer, E., Nelson, C.E., Johnson, R.L., Morgan, B.A., Tabin, C.J., 1994. Sonic hedgehog and Fgf-4 act through a signaling cascade and feedback loop to integrate growth and patterning of the developing limb bud. Cell 79, 993–003.
Li, S., Muneoka, K., 1999. Cell migration and chick limb development: chemotactic action of FGF-4 and the AER. Dev. Biol. 211, 335–47.
Lu, P., Minowada, G., Martin, G.R., 2006. Increasing Fgf4 expression in the mouse limb bud causes polysyndactyly and rescues the skeletal defects that result from loss of Fgf8 function. Development 133, 33–2.
Lubkin, S.R., Li, Z., 2002. Force and deformation on branching rudiments: cleaving between hypotheses. Biomech. Model. Mechanobiol. 1, 5–6.
Maree, A.F.M., Hogeweg, P., 2001. How amoeboids self-organize into a fruiting body: multicellular coordination in Dictyostelium discoideum. Proc. Natl. Acad. Sci. USA 98, 3879–883.
Martin, G.R., 1998. The roles of FGFs in the early development of vertebrate limbs. Genes Dev. 12, 1571–586.
McGlinn, E., Tabin, C.J., 2006. Mechanistic insight into how Shh patterns the vertebrate limb. Curr. Opin. Genet. Dev. 16, 426–32.
Meineke, F.A., Potten, C.S., Loeffler, M., 2001. Cell migration and organization in the intestinal crypt using a lattice-free model. Cell Prolif. 34, 253–66.
Meinhardt, H., Gierer, A., 2000. Pattern formation by local self-activation and lateral inhibition. BioeEssays 22, 753–60.
Merks, R.M.A., Glazier, J.A., 2005. A cell-centered approach to developmental biology. Physica A 352, 113–30.
Merks, R.M.A., Brodsky, S.V., Goligorksy, M.S., Newman, S.A., Glazier, J.A., 2006. Cell elongation is key to in silico replication of in vitro vasculogenesis and subsequent remodeling. Dev. Biol. 289, 44–4.
Miura, T., Shiota, K., Morriss-Kay, G., Maini, P.K., 2006. Mixed-mode pattern in Doublefoot mutant mouse limb—Turing reaction–diffusion model on a growing domain during limb development. J. Theor. Biol. 240, 562–73.
Mombach, J.C.M., 1999. Simulation of embryonic cell self-organization: a study of aggregates with different concentrations of cell types. Phys. Rev. E 59, R3827–R3830.
Mombach, J.C.M., Glazier, J.A., Raphael, R.C., Zajac, M., 1995. Quantitative comparison between differential adhesion models and cell sorting in the presence and absence of fluctuations. Phys. Rev. Lett. 75, 2244–247.
Mombach, J.C.M., de Almeida, R.M.C., Thomas, G.L., Upadhyaya, A., Glazier, J.A., 2001. Bursts and cavity formation in Hydra cells aggregates: experiments and simulations. Physica A 297, 495–08.
Morishita, Y., Iwasa, Y., 2008. Optimal placement of multiple morphogen sources. Phys. Rev. E 77, 041909.
Murray, J.D., Maini, P.K., Tranquillo, R.T., 1988. Mechanochemical models for generating biological pattern and form in development. Phys. Rep. 171, 59–4.
Nagai, T., Honda, H., 2001. A dynamic cell model for the formation of epithelial tissues. Philos. Mag. B 81, 699–19.
Newman, S.A., Frisch, H.L., 1979. Dynamics of skeletal pattern formation in developing chick limb. Science 205, 662–68.
Niederreither, K., Vermot, J., Schuhbaur, B., Chambon, P., Dolle, P., 2002. Embryonic retinoic acid synthesis is required for fore limb growth and anteroposterior patterning in the mouse. Development 129, 3563–574.
Niswander, L., 2002. Interplay between the molecular signals that control vertebrate limb development. Int. J. Dev. Biol. 46, 877–81.
Niswander, L., Jeffrey, S., Martin, G.R., Tickle, C., 1994. A positive feedback loop coordinates growth and patterning in the vertebrate limb. Nature 371, 609–12.
Oster, G.F., Murray, J.D., Harris, A.K., 1983. Mechanical aspects of mesenchymal morphogenesis. J. Embryol. Exp. Morphol. 78, 83–25.
Oster, G.F., Murray, J.D., Maini, P.K., 1985. A model for chondrogenic condensations in the developing limb: the role of extracellular matrix and cell tractions. J. Embryol. Exp. Morphol. 89, 93–12.
Pagan, S.M., Ros, M.A., Tabin, C., Fallon, J.F., 1996. Surgical removal of limb bud sonic hedgehog results in posterior skeletal defects. Dev. Biol. 180, 35–0.
Palsson, E., 2001. A three-dimensional model of cell movement in multicellular systems. Future Gener. Comput. Syst. 17, 835–52.
Palsson, E., Othmer, H.G., 2000. A model for individual and collective cell movement in Dictyostelium discoideum. Proc. Natl. Acad. Sci. 97, 10448–0453.
Pascoal, S., Andrade, R.P., Bajanca, F., Palmeirim, I., 2007. Progressive mRNA decay establishes an mkp3 expression gradient in the chick limb bud. Biochem. Biophys. Res. Commun. 352, 153–57.
Peskin, C.S., 1972. Flow patterns around heart valves: a numerical method. J. Comput. Phys. 10, 252–71.
Rejniak, K.A., Kliman, H.J., Fauci, L.J., 2004. A computational model of the mechanics of growth of the villous trophoblast bilayer. Bull. Math. Biol. 66, 199–32.
Rudge, T., Haseloff, J., 2005. A computational model of cellular morphogenesis in plants. In: ECAL 2005, LNAI, vol. 3630, pp. 78–7. Springer, Berlin.
Sato, K., Koizumi, Y., Takahashi, M., Kuroiwa, A., Tamura, K., 2007. Specification of cell fate along the proximal-distal axis in the developing chick limb bud. Development 134, 1397–406.
Savill, N.J., Hogeweg, P., 1997. Modeling morphogenesis: from single cells to crawling slugs. J. Theor. Biol. 184, 229–35.
Savill, N.J., Sherratt, J.A., 2003. Control of epidermal stem cell clusters by Notch-mediated lateral induction. Dev. Biol. 258, 141–53.
Scott, E.L., Britton, N.F., Glazier, J.A., Zajac, M., 1999. Stochastic simulation of benign avascular tumour growth using the Potts model. Math. Comput. Model. 30, 183–98.
Sloan, S.W., 1987. A fast algorithm for constructing Delaunay triangulations in the plane. Adv. Eng. Softw. 9, 34–5.
Stein, M.B., Gordon, R., 1982. Epithelia as bubble rafts: a new method for analysis of cell shape and intercellular adhesion in embryonic and other epithelia. J. Theor. Biol. 97, 625–39.
Sun, X., Lewandoski, M., Meyers, E.N., Liu, Y.H., Maxson, R.E., Martin, G.R. Jr., 2000. Conditional inactivation of Fgf4 reveals complexity of signalling during limb bud development. Nat. Genet. 25, 83–6.
Sun, X., Mariani, F.V., Martin, G.R., 2002. Functions of FGF signalling from the apical ectodermal ridge in limb development. Nature 418, 501–08.
Turner, S., Sherratt, J.A., Cameron, D., 2004. Tamoxifen treatment failure in cancer and the nonlinear dynamics of TGF β. J. Theor. Biol. 229, 101–11.
Vargesson, N., Clarke, J.D.W., Vincent, K., Coles, C., Wolpert, L., Tickle, C., 1997. Cell fate in the chick limb bud and relationship to gene expression. Development 124, 1909–918.
Verheyden, J.M., Lewandoski, M., Deng, C., Harfe, B.D., Sun, X., 2005. Conditional inactivation of Fgfr1 in mouse defines its role in limb bud establishment, outgrowth and digit patterning. Development 132, 4235–245.
Weliky, M., Oster, G., 1990. The mechanical basis of cell rearrangement I. Epithelial morphogenesis during Fundulus epiboly. Development 109, 373–86.
Weliky, M., Minsuk, S., Keller, R., Oster, G., 1991. Notochord morphogenesis in Xenopus laevis: simulation of cell behavior underlying tissue convergence and extension. Development 113, 1231–244.
Wolpert, L., 1969. Positional information and the spatial pattern of cellular differentiation. J. Theor. Biol. 25, 1–7.
Wolpert, L., Brockes, J., Jessell, T., Lawrence, P., Meyerowitz, E., 1998. Principles of Development. Oxford Univ. Press, Oxford.
Zajac, M., Jones, G.L., Glazier, J.A., 2003. Simulating convergent extension by way of anisotropic differential adhesion. J. Theor. Biol. 222, 247–59.
Zeng, W., Thomas, G.L., Glazier, J.A., 2004. Non-Turing stripes and spots: a novel mechanism for biological cell clustering. Physica A 341, 482–94.
Acknowledgement
This work was done by the support of Grant-in-Aids for Scientific Research by JSPS. We thank the following people for their very helpful comments: D. Drasdo, F.G. Feugier, P.K. Maini, T. Miura, J. Nakabayashi, K. Sato, and K. Tamura.
Open Access
This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
Author information
Authors and Affiliations
Corresponding author
Appendix: Procedure of numerical simulations
Appendix: Procedure of numerical simulations
1.1 Initial conditions
We start with a flat network composed of 1 layer of E-nodes and 3 layers of M-nodes (Fig. 5(a)). This initial configuration corresponds to the lateral plate without the protrusion of a limb bud yet. Several central E-nodes are defined as AER with the constant concentration of the AER-signal, c AER. The concentration of the AER-signal in mesenchyme is set to be zero, that is, c M j (0) = 0 for all M-nodes j.
1.2 Updating states
In each time step, the following processes are done.
-
(i)
Specifying AER. An E-node i satisfying {i|min i |y i − y o |} is defined as the center of AER, where y i and y o are y-coordinates of each E-node and a constant.
-
(ii)
Diffusion. According to Eq. (4), the concentrations of the AER-signal c M j are calculated for all j until their values become close to the equilibrium ones.
-
(iii)
Cell movement. One node i is randomly chosen. Depending on the node type of i, the gradient of potential energy ∇Φ M i or βΦ E i is calculated. The position of the nodes x i is updated to x i − ∇Φ i dt/μ i (Φ i } ∈ {Φ M i ,Φ E i ,Φ E i ). This process is repeated for N times, where N is the total number of E and M-nodes. The results were the same between if the update of all the cells occurs simultaneously and if the cells are updated one by one according to random order.
-
(iv)
Division of E-nodes. If the length of an E-E segment becomes the twice of its equilibrium distance due to cell movements, a new E-node is added in the middle of the E-E segment. For each division, links for all E and M-nodes are relinked by the Delaunay triangulation.
-
(v)
Division of M-nodes. One M-node is randomly chosen. According to the probability P div = f div c M i , whether the M-node divides or not is determined. If the node divides, the positions of daughter nodes are calculated based on Eqs. (5), and links for all E and M-nodes are also relinked by the Delaunay triangulation. This process is repeated N times.
-
(vi)
Relinking. For every time step, the link for all M-M and E-M links are relinked by using the Delaunay triangulation to determine the neighboring relation for each node.
For almost all simulation runs, we start with the total number of nodes N = 205 (41 E-nodes and 164 M-nodes).When N becomes 800 by the proliferation of nodes, M-nodes stop their divisions. Each run is stopped when the change in a network configuration becomes sufficiently small, which we call a “final” state.
Rights and permissions
Open Access This is an open access article distributed under the terms of the Creative Commons Attribution Noncommercial License (https://creativecommons.org/licenses/by-nc/2.0), which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
About this article
Cite this article
Morishita, Y., Iwasa, Y. Growth Based Morphogenesis of Vertebrate Limb Bud. Bull. Math. Biol. 70, 1957–1978 (2008). https://doi.org/10.1007/s11538-008-9334-1
Received:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1007/s11538-008-9334-1