# Modeling mammary organogenesis from biological first principles: Cells and their physical constraints

Progress in Biophysics and Molecular Biology

We developed a mathematical model of mammary gland based on proper biological principles: the default state of cells and the principle of organization.

Abstract

Abstract In multicellular organisms, relations among parts and between parts and the whole are contextual and interdependent. These organisms and their cells are ontogenetically linked: an organism starts as a cell that divides producing non-identical cells, which organize in tri-dimensional patterns. These association patterns and cells types change as tissues and organs are formed. This contextuality and circularity makes it difficult to establish detailed cause and effect relationships. Here we propose an approach to overcome these intrinsic difficulties by combining the use of two models; 1) an experimental one that employs 3D culture technology to obtain the structures of the mammary gland, namely, ducts and acini, and 2) a mathematical model based on biological principles. The typical approach for mathematical modeling in biology is to apply mathematical tools and concepts developed originally in physics or computer sciences. Instead, we propose to construct a mathematical model based on proper biological principles. Specifically, we use principles identified as fundamental for the elaboration of a theory of organisms, namely i) the default state of cell proliferation with variation and motility and ii) the principle of organization by closure of constraints. This model has a biological component, the cells, and a physical component, a matrix which contains collagen fibers. Cells display agency and move and proliferate unless constrained; they exert mechanical forces that i) act on collagen fibers and ii) on other cells. As fibers organize, they constrain the cells on their ability to move and to proliferate. The model exhibits a circularity that can be interpreted in terms of closure of constraints. Implementing the mathematical model shows that constraints to the default state are sufficient to explain ductal and acinar formation, and points to a target of future research, namely, to inhibitors of cell proliferation and motility generated by the epithelial cells. The success of this model suggests a step-wise approach whereby additional constraints imposed by the tissue and the organism could be examined in silico and rigorously tested by in vitro and in vivo experiments, in accordance with the organicist perspective we embrace.

Keywords: Ductal morphogenesis, Mathematical models, Organicism, Organizational closure, Acinar morphogenesis, Mammary gland morphogenesis

Modeling mammary organogenesis from biological first principles: cells and their physical constraints

Maël Montévila,b, Lucia Speronic, Carlos Sonnenscheind,e,c, Ana M. Sotod,c

aLaboratoire "Matière et Systèmes Complexes" (MSC), UMR 7057 CNRS, Université Paris 7 Diderot, 75205 Paris Cedex 13, France

bInstitut d’Histoire et de Philosophie des Sciences et des Techniques (IHPST) - UMR 8590, 13, rue du Four, 75006 Paris, France

cDepartment of Integrative Physiology and Pathobiology, Tufts University School of Medicine, Boston, MA USA

dCentre Cavaillès, République des Savoirs, CNRS USR3608, Collège de France et École Normale Supérieure, Paris, France

eInstitut d’Etudes Avancées de Nantes, France.

### Abstract

In multicellular organisms, relations among parts and between parts and the whole are contextual and interdependent. These organisms and their cells are ontogenetically linked: an organism starts as a cell that divides producing non-identical cells, which organize in tri-dimensional patterns. These association patterns and cells types change as tissues and organs are formed. This contextuality and circularity makes it difficult to establish detailed cause and effect relationships. Here we propose an approach to overcome these intrinsic difficulties by combining the use of two models; 1) an experimental one that employs 3D culture technology to obtain the structures of the mammary gland, namely, ducts and acini, and 2) a mathematical model based on biological principles.

The typical approach for mathematical modeling in biology is to apply mathematical tools and concepts developed originally in physics or computer sciences. Instead, we propose to construct a mathematical model based on proper biological principles. Specifically, we use principles identified as fundamental for the elaboration of a theory of organisms, namely i) the default state of cell proliferation with variation and motility and ii) the principle of organization by closure of constraints.

This model has a biological component, the cells, and a physical component, a matrix which contains collagen fibers. Cells display agency and move and proliferate unless constrained; they exert mechanical forces that i) act on collagen fibers and ii) on other cells. As fibers organize, they constrain the cells on their ability to move and to proliferate. The model exhibits a circularity that can be interpreted in terms of closure of constraints.

Implementing the mathematical model shows that constraints to the default state are sufficient to explain ductal and acinar formation, and points to a target of future research, namely, to inhibitors of cell proliferation and motility generated by the epithelial cells. The success of this model suggests a step-wise approach whereby additional constraints imposed by the tissue and the organism could be examined in silico and rigorously tested by in vitro and in vivo experiments, in accordance with the organicist perspective we embrace.

Keywords:  Ductal morphogenesis, mathematical models, organicism, organizational closure, acinar morphogenesis, mammary gland morphogenesis

Theory and fact are equally strong and utterly interdependent; one has no meaning without the other. We need theory to organize and interpret facts, even to know what we can or might observe. And we need facts to validate theories and give them substance. Theory and fact are equally strong and utterly interdependent; one has no meaning without the other. We need theory to organize and interpret facts, even to know what we can or might observe. And we need facts to validate theories and give them substance.
S.J. Gould (1998)

## 1. Introduction

Scientific theories provide organizing principles and construct objectivity by framing observations and experiments (Longo & Soto, this issue). On the one hand, theories construct the proper observables and on the other they provide the framework for studying them. Good theories, like Newton’s law of inertia and the conservation of momentum from his 3rd law were never abandoned but were reinstated while physics underwent further theoretical changes. Indeed, a deeper understanding of these principles was gained through E. Noether’s theorems which justify the conservation properties of energy and momenta in terms of symmetries in the state equations (van Fraassen 1989). However, a theory does not need to be “right” to guide the praxis of good experiments. Even a “wrong” theory can be useful if, when proven incorrect it is modified or even dismissed.

Here we demonstrate that the application of the principles we propose to use for the construction of a theory of organisms results in a better understanding of morphogenesis (the generation of biological form) than the common practice of using metaphors derived from the mathematical theory of information as theoretical background [for a critique see Perret and Longo in this issue, and (Longo et al. 2015)]. Our approach diverges from the biophysical methodology which is based on conservation principles and their associated symmetries, on the one hand, and with optimization principles, on the other. In contrast, biology is about an incessant breaking of symmetries (Longo et al. 2015; Longo and Montévil 2011; Longo and Soto, this issue). Taking mammary gland morphogenesis as an example, here we show that our theoretical principles are useful to provide a framework for the mathematical modelling of tissue morphogenesis. We will focus on the manner in which we propose to deal with cellular behavior.

## 2. Theoretical principles

We begin by identifying a foundational principle, the default state of cells, which is proliferation with variation and motility. This default state is a manifestation of the agency of living objects, and thus, a cause; it does not need an explanation or an external cause (Longo et al. 2015, Soto et al. 2016). The default state is what happens when nothing is done to the system. Second, we adopt the notion that organismal constraints prevent the expression of the default state. This means that constraints determine when proliferation with variation and/or motility are allowed to instantiate. Third, we consider it essential to stress that biological processes make full sense only in the context of the organism in which they take place. As stated by Claude Bernard: “The physiologist and the physician must never forget that the living being comprises an organism and an individuality. If we decompose the living organism into its various parts, it is only for the sake of experimental analysis, not for them to be understood separately. Indeed, when we wish to ascribe to a physiological quality its value and true significance, we must always refer to this whole and draw our final conclusions only in relation to its effects in the whole.” We address this aspect of biological integration using the notion of closure of constraints (see Mossio et al this issue). These constraints are considered “local invariants” because they do not change at the time scale of the process they influence. In an organism, these constraints depend collectively on each other thus attaining closure. In turn, closure provides an understanding of the relative stability of biological organizations. Fourth, organisms spontaneously undergo variation. A fundamental generator of variation is the default state. Unlike physical systems, biological ones are not framed by invariants and invariant preserving transformations. Instead, the flow of time is associated with qualitative changes of organization that cannot be stated a priori. This original feature is directly related to the historical nature of biological objects, since a specific object is the result of this unpredictable accumulation of changes (Longo et al. 2015, this issue, Montévil et al, this issue). Another notion that becomes fundamental through this principle of variation is that of contextuality. Indeed, understanding biological organization requires taking into account its interaction with the surrounding environment, both at a given time-point and through the successive environments that biological objects traverse (Soto and Sonnenschein 2005)(see Miquel and Hwang, Montévil et al and Sonnenschein and Soto in this issue). In addition to the role that constraints play in canalizing processes, they make possible the appearance of new constraints and thus changes of organization. Lastly, the framing principle states that biological phenomena should be understood as the non-identical iterations of morphogenetic processes. Biological processes iterate at all levels of organization. Organization involves iteration through the circularity of closures, but organization itself is also iterated as reproduction. Inside organisms, structures are also iterated, for example in the case of branching morphogenesis. In all cases, the principle of variation applies so that each iteration may be associated with unpredictable qualitative changes (Longo et al. 2015).

Together these principles provide a genuinely biological framework for the understanding of organismal phenomena. This framework combines the integrative viewpoint inherited from physiology, the centrality of biological variation that derives from the theory of evolution and the default state that links organismal and evolutionary biology.

## 3. The mammary gland as a model system

The mammary gland is made up of two main tissue types, namely, i) the epithelial parenchyma, its function is to produce and deliver milk, and ii) the stroma which surrounds and supports the epithelium. The stroma is composed of various cell types (fibroblasts, adipocytes, and immune cells), blood vessels, nerves, and an extracellular fibrous matrix of which the main component is collagen type-I. In the resting gland the epithelium is organized into a ductal tree. During pregnancy a second epithelial compartment, the alveoli, develop from the ducts; these are the structures that produce and secrete milk. Throughout development, reciprocal interactions between the epithelium and the stroma are responsible for the structure and function of mammary glands. Perturbations of epithelial-stromal interactions result in various pathologies including neoplasms (Soto and Sonnenschein 2011, Sonnenschein and Soto 2016).

The mammary gland undergoes morphological and functional changes throughout life. Mammary organogenesis has been studied in most detail in rodents. In mice, the mammary placodes become visible between embryonic day (E) 11 and 12 and they then develop into mammary buds by E13. At this time, several layers of mesenchyme condense surrounding the buds in a concentric fashion. In female mouse embryos the mammary bud sprouts and invades the presumptive fat pad. At E18, the mammary epithelium consists of an incipient ductal tree (Balinsky 1950; Robinson et al. 1999). Although fetal mammary gland development occurs even in the absence of receptors for mammotropic hormones suggesting that these hormones are not required at this stage, fetal mammary morphogenesis can be altered by exposure to hormonally-active chemicals (Vandenberg et al. 2007; Vandenberg et al. 2007). From the onset of puberty, the development of the mammary gland is subject to hormonal regulation. At the onset of puberty, estrogens induce the formation of club-shaped structures at the end of the ducts, called terminal end buds (TEBs). Thereafter, the epithelium begins to fill the fat pad and branches. Progesterone induces lateral branching. If pregnancy occurs, prolactin in combination with estrogen and progesterone initiates a characteristic lobuloalveolar development [reviewed in (Brisken and O’Malley 2010)]. When lactation ceases, involution of the alveolar structures occurs and the mammary gland returns to its resting state.

### 3.1. Biological models for the study of mammary gland biology

3D culture systems allow for the dynamic study of epithelial morphogenesis and the organization of the stroma. These models are intended to mimic conditions prevailing in a living organism while reducing the number of constraints present in vivo to those which theoretically and/or empirically are considered to be the most relevant ones for the subject study. When designing a 3D culture model one must first define the main characteristics of the target tissue that the model aims to reproduce, and which stage of mammary gland development the investigator is interested in reproducing in vitro. The objective of mimicking the tissue of origin is tempered by the need to make it manageable by reducing the model to a few components. This allows the researcher to infer from these results the contribution of these components to the mammary gland phenotype in vivo. The resulting model may then be compared to more complex ones resulting from the step-wise addition of relevant components, and eventually to the behavior of the gland in situ.

Epithelial–stromal interactions can be studied using 3D co-cultures of epithelial and stromal cells by analyzing matrix remodeling and epithelial morphogenesis. With regard to matrices the most biologically relevant ones are those that provide the structure and rigidity of the model tissue that allows the cellular components to attain characteristics seen in the breast.

Here we focus on collagen-based matrices since collagen is a main component of the mammary stroma that allows for breast epithelial cells to organize into structures that closely resemble those observed in vivo (Krause et al. 2012; Krause et al. 2008; Dhimolea et al. 2010; Speroni et al. 2014; Barnes et al. 2014).

#### 3.1.1. MCF10A 3D culture model

To test the mathematical model, we used data generated using the MCF10A 3D culture model described in (Barnes et al. 2014). In this 3D culture model, the proportion of acinar and ductal structures can be modified by changing the concentration of reconstituted basement membrane (Matrigel) in the extracellular matrix (ECM) (Krause et al. 2008; Barnes et al. 2014). Briefly, MCF10A human breast epithelial cells were seeded in either bovine type-I collagen matrix or in mixed matrices containing collagen and Matrigel at 5 and 50% v/v. The final collagen concentration in all gels was 1.0 mg/ml. Gels were prepared by carefully pouring 500 µl of the cell-matrix mixture into wells of a glass-bottomed 6 well plate. The gels were allowed to solidify for 30 min at 37°C before adding 1.5 ml of cell maintenance culture medium into each well. Cultures were maintained at 37°C in an atmosphere containing 6% CO2/94% and 100% humidity 5 or 7 days, and the medium was changed every 2 days. Live cell-imaging was conducted with a Leica SP5 microscope (Leica-Microsystems, Germany). Simultaneous reflectance confocal microscopy (RCM) and bright-field images were acquired using a 40 ×, 1.1 numerical aperture, water immersion objective with separate photomultiplier tubes using the Argon 488 nm laser line. The RCM signal was collected between wavelengths of 478–498 nm with a pinhole size of 57µm. For additional details see (Barnes et al. 2014).

### 3.2. Current mathematical models of mammary gland morphogenesis

Several types of mathematical models have been used for the study of mammary gland morphogenesis. These models focus on distinct aspects of this phenomenon, and address a given aspect at a particular time and space scale. These diverse models also focus on different kinds of determinants, some chemical, some mechanical, some both.

 Object studied Modeling cell proliferation Modeling motility Implicit default state Reference Collagen network remodeling Not discussed. ECM constrains movement. Motility. (Harjanto and Zaman 2013) Collagen and fibroblasts Not discussed. Cells spontaneously exert forces and move. Stress may prevent motion. Motility. (Dallon et al. 2014) Acinus in 2D Limited space prevents proliferation. Although not included in the math model, it is stated that proliferation is regulated by signals/growth factors. Not discussed. Proliferation (quiescence is invoked the discussion). (Rejniak and Anderson 2008) Acinus in 3D Only cells adjacent to the basement membrane proliferate. Cells have a proliferative potential that decreases at each division. Not discussed. Shifts from proliferation to quiescence as time elapses. Not discussed for inner cells. (Tang et al. 2011) Review on biophysical cell self-assembly No proliferation. Cells show trend to move modeled by a parameter formally similar to temperature. If this "temperature parameter" is an intrinsic property of cells, the default state is motility, otherwise quiescence. (Neagu 2006) Terminal End Buds Not causal analysis (measured or assessed indirectly). NA (Paine et al. 2016) Epithelial tree Components resulting from the action of matrix metalloproteinases are inferred to stimulate cell proliferation. Not discussed. Quiescence. (Grant et al. 2004)

Oftentimes the modelers make mathematical hypotheses on the behavior of cells without making explicit the broader biological significance of these hypotheses. Our mathematical model instead is based on general biological principles. Some models seem to implicitly rely on the default state that we propose: that is, cells spontaneously move or proliferate and the model discusses specific constraints on these behaviors. In other models, cells are quiescent without any constraints acting on them and chemicals stimulate proliferation and/or movement without removing constraints. Still, in other models even opposite behaviors co-exist and no explicit attempt to reconcile these opposites is made. In Table 1 we review the way several models deal with cell behavior.

In these models assumptions on cell behavior are largely ad hoc, varying from one model to another. Occasionally, the model and its interpretation in the Discussion are inconsistent. For example, in (Rejniak and Anderson 2008) proliferation is constrained by the available space but it is also attributed to "signals".

Various models implicitly adopt the premise that the default state of cells is proliferation or motility. Other models (Grant et al. 2004) do not adopt this premise; however they may be reinterpreted by adopting our principles. In the Grant et al model, the matrix metalloproteinases are assumed to have a positive effect on proliferation by degrading the ECM. Because the authors also infer that the default state is quiescence, they assume that this degradation produces a chemical which stimulates proliferation. A simpler hypothesis that we favor is that this degradation removes the mechanical constraint of the ECM on the default state of cell proliferation with variation and motility. To our knowledge, no model of mammary gland morphogenesis has taken both aspects of our default state into account.

The main aim of this article is to emphasize that making explicit the assumption that the default state is proliferation with variation and motility enables us to model morphogenesis on precise theoretical bases. The models reviewed above focus on different aspects of mammary gland morphogenesis. For instance, models of ECM remodeling focus on fibroblasts and do not discuss epithelial morphogenesis. Several agent-based models focus on the behavior of epithelial cells during acinus formation in conditions that preclude ductal morphogenesis. Finally, other models focus on larger scale organogenesis but do not provide a detailed account of cellular behavior. Determinants of these models are either chemical, mechanical forces or empty space. Note however that in most cases morphogenesis is assumed to be driven by chemical interactions (Iber and Menshykau 2013). This article will focus on the formation of epithelial structures, mostly ducts but also acini on the basis of mechanical interactions between cells and between cells and the ECM.

Mammary gland morphogenesis in vivo requires an interplay between the epithelial compartment and the stroma which contains fibroblasts, adipocytes and ECM. Our simplified model contains epithelial cells and a stroma devoid of cells but containing the ECM which in vivo would be a product of the stromal cells. The assumption that the ECM is the mechanically important component of the stroma opens the possibility to explore a set of minimal and manageable conditions that allow for a bi-stable determination of acini and ducts, the two structures that characterize the mammary gland parenchyma. In the rest of the text we will discuss all experimental results and mathematical assumptions under the hypothesis that the default state of cells is proliferation with variation and motility.

### 3.3. From the 3D culture model to a mathematical model

#### 3.3.1. Proliferation

As expected from the default state, breast epithelial cells proliferate maximally in serumless medium. Addition of hormone-free serum to culture medium results in a dose-dependent inhibition of cell proliferation. This is due to the effect of serum albumin, a main constraint for the proliferation of estrogen-target cells. This constraint could be lifted by lowering albumin or serum concentration, a procedure easily performed in 2D cell cultures, or by adding estrogens, which is the natural way the organism uses to neutralize the effect of serum albumin (Sonnenschein et al. 1996).

Additional constraints are those imposed by cell-cell contact, which are weaker in 2D culture than in 3D. This is because 3D structures allow for cell-cell contact in practically all directions, while in 2D these contacts are more restricted. For example, the estrogen-sensitive T47D cell line is inhibited from proliferating when placed in medium containing serum. This constraint is lifted by the addition of estrogen. When placed in 2D culture, the cell number ratio between serum plus estrogen and serum without estrogen is 3.75, while the same experiment performed in 3D culture results in a ratio of 2.5 (Speroni et al. 2014). These data indicate that the organization of cells into epithelial structures constrains cell proliferation more effectively than the one in 2D culture, where cells are not organized into closely packed epithelial structures. Additionally, in 2D cultures cells attach to the bottom of the culture dish, a surface that is exceedingly more rigid than the conditions these very cells encounter within the tissue of origin. In contrast, 3D cultures could be engineered to mimic the rigidity of the tissue of origin.

#### 3.3.2. Motility

In classical mechanics, motion is a consequence of external forces. In biology, the situation is different but compatible with mechanics. Cells need a configuration of forces to be able to move; for example, they need a support to be able to crawl on it, or fibers to which they can attach and pull in order to move. However, unlike inert objects, cells initiate movement on their own. In other words, they are autonomous agents which express their default state (Sonnenschein and Soto, 1999). Cells move unless there are constraints which prevent them from doing so. Reciprocally, a given mechanical force acting upon biological entities such as cells produces clearly different effects than forces acting on inert matter (Soto et al. 2008; Longo and Montévil 2014). For example, gravity in mechanics is just a force proportional to the mass of the object and oriented towards the center of the earth. In biology, however, gravity becomes a constant constraint, which has not been altered since the origin of life. Biological organization reacts to it in various ways. For example, swim bladders, wings, limbs and tree trunks are responsive to gravitational force but are not explained by it. The behavior of molecules in cells and the overall behavior of tissues and organs are massively impacted when in microgravity conditions (Bizzarri et al. 2014).

The constraints to motility that cells experience in the tissue environment can be modeled in a 3D culture system. In this system, the tissue environment is recreated by the matrix in which the cells are seeded. Cells exert their motility by using filopodia and pseudopodia. When they encounter a structure such as a fiber they can use it for locomotion or pull on it to attach. Breast epithelial cells seeded in a fibrilar matrix emit projections in all directions soon after seeding (Figure 1; video 1). This process allows cell elongation which precedes the formation of ducts (tubular structures) and branching (Barnes et al. 2014). Amorphous, non-fibrilar matrix proteins as well as fibers also oppose migration and cause resistance due to their relative rigidity and their lack of pores (Figure 2). Breast epithelial cells growing in a non-fibrilar matrix display limited motility and emit short projections into the matrix that retract soon afterwards. Cells rotate and divide resulting in the formation of an acinus, a sphere with a central lumen (Tanner et al. 2012) (Figure 3; video 2). Cell movement is also constrained by the pore size of the matrix, this is mostly determined by fiber alignment, fiber density and abundance of non-fibrous matrix materials. Other factors such as pore size and matrix rigidity could contribute to the morphological differences of the epithelial structures. Pore size is bigger in the fibrilar matrix than in the globular matrix (Figure 2). The globular matrix is stiffer than the fibrilar matrix, however this proved to be a minor contributor to epithelial phenotype compared to collagen fiber distribution (Barnes et al. 2014).

Cell adhesion. Cells can adhere to each other after division of their progenitor; they can also attach to any migrant cells they encounter. Once an epithelial structure, such as a duct or an acinus, is formed, adhesion will be maintained as new cells are formed or replaced. Adhesion constrains the motion of cells. During ductal morphogenesis, single cells can detach from the main structure and may also be incorporated back into the structure from which they detached (Figure 4, video 3). This phenomenon suggests an environment-sensing strategy as well as a means used by epithelial structures to modify the matrix prior to growth in a certain direction. During lumen formation, cells migrate toward the periphery of the epithelial structure leaving a space that will be filled by fluid.

#### 3.3.3. Determination of the system

The increase in cell number due to the unconstrained default state brings about re-distribution of fluids, reorganization of fibers and a certain degree of matrix compression, and/or matrix degradation. Elongation is accompanied by fiber organization into bundles projecting in the direction of the ductal tip (Barnes et al. 2014). Collagen bundles facilitate the merging of epithelial structures initially positioned at a long distance range (Guo et al. 2012)(Figure 5, Video 4).

Acellular collagen gels contain small fibers that are not seen using the classical picrosirius red/polarized light method to visualize collagen fibers. However, a similar collagen gel containing cells reveals a very different picture; twenty-four hours after seeding the cells, small yellow fibers are detected (Figure 6). This experiment shows that cells organize collagen fibers. In other words, cells exert forces upon fibers, and fibers transmit these forces for quite long distances (Guo et al. 2012).

As collagen fibers progressively organize, they constrain the proliferation and motility of cells. These constraints may be positive like those that facilitate cell migration along fibers or negative like the ones hindering migration orthogonally, and those due to “pore” formation when the collagen fibers organize into a network. Additionally, cells constrain other cells mechanically. The reciprocal interactions between the collagen and the cells are illustrated in Figure 7 which emphasizes the separation of the system into a physical component, the collagen, and a biological component, the cells. In this diagram the behavior of the cells is determined by the default state and the constraints exerted on it.

## 4. Implementing a model based on these principles

We propose a mathematical model of morphogenesis in 3D cultures that we analyze by computer simulations. In this article this model is used as a proof of concept. A mathematical description of our model is provided in Appendix A. Our model is based on the principles that we propose for the construction of a theory of organisms. Our theoretical framework restricts what is acceptable in order to model cellular behavior. For example, it is unacceptable for cells to be proliferatively quiescent without an explicit constraint keeping them in this state. Of course, in the presence of a strong constraint cells will become quiescent and will remain as such for the duration of the constraint.

### 4.1. Components of the model

The 3D culture gel is represented as a lattice, in three dimensions. The model is based on different layers that interact with each other:

1. A mechanical forces layer. Each elementary cube of our model exerts forces on the adjacent cubes. Forces propagate in space. Their orientation and propagation depends on the orientation of collagen fibers.
2. A cellular layer, which has three possible states: presence of a live cell, a dead cell (in the case of lumen formation), or no cell, in which case the cube contains ECM.
3. The collagen layer is approached at a mesoscopic scale and does not represent individual fibers. Each cube of collagen has a main orientation which will transmit forces farther, and this orientation is random in the initial conditions. This orientation is also relevant for cellular behavior (see below). Note that this layer represents the cytoskeleton in the cubes occupied by cells. Collagen tends to align with forces that are exerted on it.
4. An inhibitory layer. Cells around the lumen produce a short range inhibition of both proliferation and motility. This layer can be interpreted either as the formation of a basal membrane or as a chemical inhibitor. Both aspects are biologically relevant and could be split into two different types of inhibition in future work.
5. A “nutrient" layer governed by diffusion. Cells consume nutrients. If there are not enough nutrients, cells die leading to lumen formation. The main point here is not to produce an accurate account of lumen formation, but instead to take into account the consequences of lumen formation on morphogenesis as a lumen would disrupt the transmission of forces.

### 4.2. Cellular behavior:

Cells express four different behaviors:

1. Cells exert forces on each other and on the collagen network. The orientation of these forces is influenced by both the neighboring cells and there is also a random component. The magnitude of the force depends on the content of the neighboring cube. More precisely, the force exerted depends on the position with respect to the cell, the orientation of the cytoskeleton of the cell and the orientation of the collagen fibers in the case of an action upon a collagen fiber or fascicle. Additionally, cells tend to oppose strong mechanical stress.
2. Cells have a defined generation time and divide, except when there are constraints which prevents them from doing so. The new cell will occupy a random spot. When this spot is already occupied, the cell cannot proliferate; the cell will make another attempt at the next iteration of the simulation loop.
3. Cells move randomly unless this movement is constrained.
4. Cells die when they lack nutrients. Unlike the other behaviors, this one is an ad-hoc addition to create a lumen, since the steps involved in lumen formation are not well known. There is evidence for cell death and for cell migration, but the cue provoking lumen formation is unknown.

In this analysis, cell motility has two components: the forces exerted by cells and cell movement. Cell movement involves detachment and reattachment to other cells and to the extracellular matrix.

### 4.3. Constraints on proliferation:

1. Cells tend to proliferate along the direction of forces. The stronger the forces, the stronger the constraint becomes. In the case of three cells being aligned, the one in the middle may be unable to proliferate if there is a force in the direction of their alignment.
2. Strong mechanical stress slows down proliferation; this constraint is not required for the model to work.
3. The inhibitory layer prevents proliferation.

### 4.4. Constraints on cell movement:

1. A strong mechanical stress slows down movement; this is not required for the model to work. Also, as the number of neighboring cells increases, the stronger the effect of cell adhesion, which prevents the cell from initiating movement.
2. Movement is facilitated when it occurs along collagen fibers.
3. Movement is facilitated towards other cells.
4. The inhibitory layer prevents movement.

### 4.5. Results of implementing the mathematical model

The model that we propose exhibits a circularity that can be interpreted in terms of closure (Figure 8). This circularity concerns constraints acting on i) processes such as physical forces and ii) directly on the default state.

The implementation of this mathematical model generates biologically relevant results. When in collagen, cells form elongated structures that can be interpreted as ducts. When we remove the effects of collagen’s fibrilar structure, in order to mimic a globular matrix (Matrigel), cells form spherical structures.

#### 4.5.1. Duct formation

A single cell is surrounded by collagen. This cell starts to exert forces on the collagen and organizes it. The cell will proliferate and may move. The epithelial structure acquires additional cells through cell proliferation. Cells move but they mostly stay attached to the structure. Cells in the middle of the structure can neither move nor proliferate. The cells exert forces on each other and on the collagen. This leads to the appearance of a dominant direction in which these forces are exerted and collagen is reorganized. This direction is often that of the forces exerted initially by the first cell and also depends on initial collagen configuration. Motility and proliferation are facilitated along this main direction, while they are inhibited in the direction perpendicular to this force. As a result, the structure acquires an elongated shape. The structure grows following this dynamic until it reaches a size large enough for the lumen to form at the thickest part of the structure, which is close to the initial position of the first cell. In the context of lumen formation the inhibitor constrains proliferation and motility and inhibits the growth in the width of the structure in the vicinity of a lumen. In contrast, the tips of the elongated structure are not inhibited and elongate further without apparent restriction (Figure 9, videos 5 and 6).

In some cases, the constraints leading to duct formation can become disorganized at one of the tips. In such case, the main direction of forces exerted may change, leading to a change in the growth direction of the duct. In other cases, this disorganization leads to a bifurcation in the main direction of the constraints and to branching (videos 7 and 8). These phenomena will be the subject of further studies.

#### 4.5.2. Acinus formation:

The addition of Matrigel, to a collagen matrix changes the gel properties by coating the collagen fibers and thus hindering fiber organization (Barnes et al. 2014). In this case, collagen fibers are not accessible to the epithelial cells, which prevents the establishment of a main direction in the forces exerted by the cells. Moreover, these forces are exerted exclusively on cells rather than on fibers. As a result, the epithelial structure grows in an isotropic manner, and when the lumen is formed and the inhibitor is secreted, all cells are constrained. The acinus being formed is a smooth structure due to the constrained motility of cells, each of which is constrained by the many cells surrounding them (Figure 10, video 9).

## 5. The in vitro system and the organism

How faithfully does this in vitro system represent a phenomenon occurring inside the organism? By accepting the reciprocal relationship between the whole (organism) and its parts, it is difficult to conceive that a sub-system would operate in a biologically relevant fashion when outside the organism. However, there are examples in embryology of the relative autonomy of parts at a particular point in space-time, as for example, the autonomy of limb morphogenesis upon transplantation of the limb morphogenetic field to an ectopic location. Additionally, parceling a whole does not prevent us from ascertaining whether or not the in vitro process arrives at similar outcomes as those observed in vivo. How faithfully does this in vitro system represent a phenomenon occurring inside the organism? By accepting the reciprocal relationship between the whole (organism) and its parts, it is difficult to conceive that a sub-system would operate in a biologically relevant fashion when outside the organism. However, there are examples in embryology of the relative autonomy of parts at a particular point in space-time, as for example, the autonomy of limb morphogenesis upon transplantation of the limb morphogenetic field to an ectopic location. Additionally, parceling a whole does not prevent us from ascertaining whether or not the in vitro process arrives at similar outcomes as those observed in vivo.

Which constraints are required for a relevant model of tissue morphogenesis?

1. Biological meaning is construed by applying similar constraints to those which operate in vivo and which seem to play a role in the determination of the phenomenon. In this way, we can “reduce” the number of constraints to those necessary to answer our specific question. Deviations from expected results may potentially indicate additional constraints which could then be identified.
2. Constraints that are absolutely required to allow the cells to continue being alive (pH, nutrients, temperature) and to express their default state in conditions that replicate as much as possible the conditions present in the organism. The “optimization” of these basal conditions is done experimentally by ascertaining that the cells can proliferate as fast as possible and that the cellular phenotype that we wish to study is obtained. This is done in 2D cultures.
3. Another consideration is the historicity and specificity of biological systems. For example, fibroblasts to be used in a 3D culture of the mammary gland are isolated from human breast tissue and used within 6 passages to avoid excessive deviation from the in situ condition. Along the same lines, established epithelial cells are maintained in standardized conditions that result in the reproducibility of the phenomenon studied (i.e., duct formation) and the cell phenotype (i.e., response to a given hormone).
4. In vitro 3D models allow researchers to manipulate constraints beyond the range operating in vivo. That is, constraints are determined by the organism and its parts, while in the in vitro model the researcher also plays a direct role in modifying these constraints and parameterizing them. For example, we can manipulate the rigidity of the mammary gland model to that of bone, and learn how rigidity affects shape beyond the limits imposed by the organism (Weaver et al. 1995). This type of manipulation revealed that high rigidity inhibits lumen formation and makes epithelial structures disorganize in a way reminiscent of neoplasms (Paszek et al. 2005).
5. Specific organism level constraints induce the tissue to undergo morphological changes required for proper organ function at the right time. Hormone action on the mammary gland is an example of this. At the onset of puberty, estrogen influences the formation of TEBs, the structure at the end of the ducts that invades the stroma and guides ductal growth until the ductal tree fills the fat pad. Progesterone promotes side-branching and during pregnancy prolactin facilitates acinar development in preparation for lactation. A classical biochemical interpretation would reduce the phenomenon to a hormone-receptor interaction triggering signaling pathways inside the cell. However, this approach does not have the capacity to explain the shape changes resulting from these hormonal influences. From the tissue perspective, exposure to hormones leads to changes in collagen fiber organization which enable the cells to generate various epithelial organization patterns. In a hormone-sensitive 3D culture model, epithelial structures resulting from exposure to estrogen in combination with a progestogen or prolactin were more irregular in shape than the elongated, smooth structures resulting from exposure to estrogen alone. Consistently, combined hormone treatment resulted in higher collagen density variability within 20µm from the epithelial structure compared to E2 alone (Speroni et al. 2014) (Figure 11).

## 6. Conclusions

We posited that experimental research guided by the global theoretical approach that we are proposing would be different from that of the prevailing ones which are mostly guided by the metaphors of information, signal and program borrowed from mathematical information theories (Longo et al. 2012). Often times, modelers entering biological research treat biological objects as if they were either physical objects or computer programs. In this issue, we have presented critiques to these approaches and suggested biological singularities that have to be taken into consideration when constructing a biological theory (Longo and Soto, this issue, Longo et al. 2015). On the one hand, the metaphorical use of information pushes the experimenter to seek causality in terms of discrete structures, namely molecules, in particular DNA. This view precludes physical “constraints”, like the ones analyzed here, to causally contribute to the generation and maintenance of the organism unless they are digitally encoded as molecular signs. On the other hand, biophysical approaches involve hypotheses which are not sound in biology, such as optimization principles in the behavior of cells. The approach adopted herein is based on two of the principles proposed as foundations for a theory of organisms, namely, the default state and organizational closure. We explored whether it would be possible and informative to model epithelial glandular morphogenesis from these biological principles, rather than the usual procedure of transferring mathematical structures developed for the understanding of physical phenomena into biological ones or proposing that cells follow a program. It is worth stressing that the former represents true mathematical modeling which is based on the theoretical framework of the discipline to which the modeled phenomenon pertains, while the latter, properly described as imitation, uses principles from one discipline and applies them to another without a critical appraisal of their theoretical meaning when transported into a different theoretical context.

In this initial modeling effort, applying the two principles (default state and constraints leading to closure) were sufficient to show the formation of ducts and acini. Cells generated forces that were transmitted to neighboring cells and collagen fibers, which in turn created constraints to movement and proliferation (Figures 7 & 8). Additionally, constraints to the default state are sufficient to explain ductal and acinar formation, and point to a target of future research, namely, the inhibitors of cell proliferation and motility which in this mathematical model are generated by the epithelial cells. Finally, the success of this modeling effort performed as a “proof of principle” opens the possibility for a step-wise approach whereby additional constraints imposed by the tissue (e.g., additional cell types) and the organism (e.g., hormones) could be assessed in silico and rigorously tested by in vitro and in vivo experiments.

## Acknowledgments

This work was conducted as part of the research project “Addressing biological organization in the post-genomic era” which is supported by the International Blaise Pascal Chairs, Region Ile de France (AMS: Pascal Chair 2013). Additional support to AMS was provided by Grant R01ES08314 from the U. S. National Institute of Environmental Health Sciences and Maël Montévil was supported by Labex "Who am I?", Laboratory of Excellence No. ANR-11-LABX-0071. The funders had no role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript. The authors are grateful to Cheryl Schaeberle for her critical input and Clifford Barnes for image acquisition. The authors have no competing financial interests to declare.

## 7. References

1. Balinsky, B.I. 1950. On the prenatal growth of the mammary gland rudiment in the mouse. J. Anat. 84, 227-235.
2. Barnes, C., Speroni, L., Quinn, K., Montévil, M., Saetzler, K., Bode-Animashaun, G., McKerr, G., Georgakoudi, I., Downes, S., Sonnenschein, C., Howard, C.V., Soto, A.M. 2014. From single cells to tissues: interactions between the matrix and human breast cells in real time. PLoS ONE 9, e93325.
3. Bizzarri, M., Cucina, A., Palombo, A., Masiello, M.G. 2014. Gravity sensing by cells: mechanisms and theoretical grounds. Rend. Fis. Acc. Lincei 25 Suppl 1, S29-S38.
4. Brisken, C., O’Malley, B. 2010. Hormone action in the mammary gland. Cold Spring Harb. Perspect. Biol. 2, a003178.
5. Dallon, J.C., Evans, E.J., Ehrlich, P. 2014. A mathematical model of collagen lattice contraction. J. R. Soc. Interface 11, 20140598.
6. Dhimolea, E., Maffini, M.V., Soto, A.M., Sonnenschein, C. 2010. The role of collagen reorganization on mammary epithelial morphogenesis in a 3D culture model. Biomaterials 31, 3622-3630.
7. Gould, S.J. 1998. Leonardo’s Mountain of Clams and the Diet of Worms: Essays on Natural History, 155.
8. Grant, M.R., Hunt, C.A., Xia, L., Fata, J.E., Bissell, M.J. 2004. Modeling mammary gland morphogenesis as a reaction-diffusion process, in: Engineering in Medicine and Biology Society, JEMBS ’04 26th Annual International Conference of the IEEE, pp. 679-682.
9. Guo, C.L., Ouyang, M., Yu, J.Y., Maslov, J., Price, A., Shen, C.Y. 2012. Long-range mechanical force enables self-assembly of epithelial tubular patterns. Proc. Nat. Acad. Sci. USA 109, 5576-5582.
10. Harjanto, D., Zaman, M.H. 2013. Modeling extracellular matrix reorganization in 3D environments. PLoS One. 8, e52509.
11. Iber, D., Menshykau, D. 2013. The control of branching morphogenesis. Open Biol. 3, 130088.
12. Krause, S., Jondeau-Cabaton, A., Dhimolea, E., Soto, A.M., Sonnenschein, C., Maffini, M.V. 2012. Dual regulation of breast tubulogenesis using extracellular matrix composition and stromal cells. Tissue Eng Part A. 18, 520-532.
13. Krause, S., Maffini, M.V., Soto, A.M., Sonnenschein, C. 2008. A novel 3D in vitro culture model to study stromal-epithelial interactions in the mammary gland. Tissue Eng. Part C Methods 14, 261-271.
14. Longo, G., Miquel, P.-A., Sonnenschein, C., Soto, A.M. 2012. Is information a proper observable for biological organization? Prog. Biophys. Mol. Biol. 109, 108-14.
15. Longo, G., Montévil, M. 2011. From physics to biology by extending criticality and symmetry breakings. Prog. Biophys. Mol. Biol. 106, 340-347.
16. —— 2014. Perspectives on Organisms: Biological Time, Symmetries and Singularities. Springer, Berlin.
17. Longo, G., Montévil, M., Sonnenschein, C., Soto, A.M. 2015. In search of principles for a theory of organisms. J. Biosci. 40, 955-968.
18. Longo, G. & Soto, A.M. 2016, this issue. Why do we need theories?, Progress in Biophysics and Molecular Biology, Available online 4 July 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.06.005
19. Miquel, P.-A., Hwang, S.-Y 2016, this issue. From physical to biological individuation, Progress in Biophysics and Molecular Biology, Available online 16 July 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.07.002
20. Montévil, M., Mossio, M., Pocheville, A. & Longo, G. 2016, this issue. Theoretical principles for biology: Variation, Progress in Biophysics and Molecular Biology, Available online 13 August 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.08.005
21. Mossio, M., Montévil, M., Longo, G. 2016, this issue. Theoretical principles for biology: Organization, Progress in Biophysics and Molecular Biology, Available online 10 August 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.07.005
22. Neagu, A. 2006. Computational modeling of tissue self-assembly. Mod. Phys. Lett. 20, 1217.
23. Paine, I., Chauviere, A., Landua, J., Sreekuman, A., Cristini, V., Rosen, J., Lewis, MT. 2016. A geometrically-constrained mathematical model of mammary gland ductal elongation reveals novel cellular dynamics within the terminal end bud. PLoS Comput. Biol. 12, e1004839.
24. Paszek, M.J., Zahir, N., Johnson, K.R., Lakins, J.N., Rozenberg, G.I., Gefen, A., Reinhart-King, C.A., Margulies, S.S., Dembo, M., Boettiger, D., Hammer, D.A., Weaver, V.M. 2005. Tensional homeostasis and the malignant phenotype. Cancer Cell 8, 241-254.
25. Perret, N. & Longo, G. 2016, this issue. Reductionist perspectives and the notion of information, Progress in Biophysics and Molecular Biology, Available online 22 July 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.07.003
26. Rejniak, K.A., Anderson, A.R. 2008. A computational study of the development of epithelial acini: II. Necessary conditions for structure and lumen stability. Bull. Math. Biol. 70, 1450-1479.
27. Robinson, G.W., Karpf, A.B.C., Kratochwil, K. 1999. Regulation of mammary gland development by tissue interaction. J. Mammary Gland Biol. Neoplasia 4, 9-19.
28. Sonnenschein C. and Soto A.M. The Society of Cells: Cancer and Control of Cell Proliferation. New York: Springer Verlag; 1999.
29. Sonnenschein, C., Soto, A.M., Michaelson, C.L. 1996. Human serum albumin shares the properties of estrocolyone-I, the inhibitor of the proliferation of estrogen-target cells. J. Steroid Biochem. Molec. Biol. 59, 147-154.
30. Sonnenschein, C., Soto, A.M., Carcinogenesis explained within the context of a theory of organisms, Progress in Biophysics and Molecular Biology, Available online 3 August 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.07.004
31. Soto, A.M., Longo, G., Montévil, M., Sonnenschein, C. 2016, this issue. The biological default state of cell proliferation with variation and motility, a fundamental principle for a theory of organisms, Progress in Biophysics and Molecular Biology, Available online 2 July 2016, ISSN 0079-6107, doi: 10.1016/j.pbiomolbio.2016.06.006
32. Soto, A.M., Sonnenschein, C. 2005. Emergentism as a default: cancer as a problem of tissue organization. J. Biosci. 30, 103-118.
33. Soto, A.M., Sonnenschein, C. 2011. The tissue organization field theory of cancer: A testable replacement for the somatic mutation theory. BioEssays 33, 332-340.
34. Soto, A.M., Sonnenschein, C., Miquel, P.-A. 2008. On physicalism and Downward Causation in Developmental and Cancer Biology. Acta Biotheoretica 56, 257-274.
35. Speroni, L., Whitt, G.S., Xylas, J., Quinn, K.P., Jondeau-Cabaton, A., Georgakoudi, I., Sonnenschein, C., Soto, A.M. 2014. Hormonal regulation of epithelial organization in a 3D breast tissue culture model. Tissue Eng. Part C Methods 20, 42-51.
36. Tang, J., Enderling, H., Becker-Weimann, S., Pham, C., Polyzos, A., Chen, C.-Y., Costes, S.V. 2011. Phenotypic transition maps of 3D breast acini obtained by imaging-guided agent-based modeling. Integrative Biology 3, 408-421.
37. Tanner, K., Mori, H., Mroue, R., Bruni-Cardoso, A., Bissell, M.J. 2012. Coherent angular motion in the establishment of multicellular architecture of glandular tissues. Proc. Nat. Acad. Sci. USA 109, 1973-1978.
38. van Fraassen, B.C. 1989. Laws and Symmetry. Oxford University Press, Oxford.
39. Vandenberg, L.N., Maffini, M.V., Wadia, P.R., Sonnenschein, C., Rubin, B.S., Soto, A.M. 2007. Exposure to environmentally relevant doses of the xenoestrogen bisphenol-A alters development of the fetal mouse mammary gland. Endocrinology 148, 116-127.
40. Weaver, V.M., Howlett, A.R., Langton-Webster, B., Petersen, O.W., Bissell, M.J. 1995. The development of a functionally relevant cell culture model of progressive human breast cancer. Seminars in Cancer Biology 6, 175-184.

## Appendix

### A. Mathematical description of the model

In this appendix, we describe the mathematical model of epithelial morphogenesis in collagen gels that we propose as a proof of concept for the use of the notion of proliferation with variation and motility as the default state of cells.

A gel is described as a three dimensional material whose properties change over time.

In our simulations, collagen and more generally the system is approximated as a three dimensional lattice with different dimensions corresponding to the different scalar or vector fields relevant to the system (collagen orientation, chemical concentrations, presence of cells, etc). The program has a main loop which corresponds to the update of cellular behaviors. The updates of collagen or chemicals are performed by sub-loops. This implementation is justified by the assumption that the characteristic speed of processes described by the sub-loops is faster than the cellular changes. All coefficients are given with respect to the time scale of the main loop.

Initially, every element of the array has a random collagen orientation with a uniform distribution and every element is independent of the other. The nutrient layer is set uniformly to $1$ and the inhibitory layer to $0$. A single cell is typically put in position $\left(25,25,25\right)$ in a $50×50×50$ array.

#### A.1. Collagen and forces

We only take into account the mean orientations of the fibers. This mean orientation provides a vector for every element of the lattice, thus defining a vector field $\stackrel{⃗}{C}\left(x,y,z,t\right)$. Collagen orientation vectors have norm $1$. In order to lighten the notation we will keep the normalization of the vectors implicit in this text. Note that for all intend and purpose, vectors $\stackrel{⃗}{C}\left(x,y,z,t\right)$ and $-\stackrel{⃗}{C}\left(x,y,z,t\right)$ are equivalent because we assume that collagen fibers are not oriented; thus, this symmetry is respected in all equations.

Forces are represented by several vector fields ${\stackrel{⃗}{F}}_{i}\left(x,y,z,t\right)$, with $i=x$, $y$, or $z$. ${\stackrel{⃗}{F}}_{x}\left(x,y,z,t\right)$ corresponds to the force exerted by the element in position $\left(x-1,y,z\right)$ on the element in position $\left(x,y,z\right)$ at time $t$. Note that according to the principle of reaction of classical mechanics, $-{\stackrel{⃗}{F}}_{x}\left(x,y,z,t\right)$ is then the force exerted by the element in position $\left(x,y,z\right)$ on the element in position $\left(x-1,y,z\right)$. Collagen is considered as an anisotropic elastic material at short time scales. Collagen orientation is altered at a larger time scales. We simulate the collagen and force propagation by the finite difference method.

The anisotropy of collagen is modeled by a Young’s modulus which depends on collagen orientation:

 $\left(\left(1-\alpha \right)\stackrel{⃗}{C}\left(x,y,z,t\right)+\alpha \right)\left(1-d\right)$

(we use $d=0.1$ and $\alpha =0.5$). The representation of collagen that we use is a crude macroscopic representation, but we feel that it is a computationally “lightweight” approach that is sufficient to understand the morphogenesis of the epithelial structures of interest. Note also that both the mathematical and the biological models that we discuss aim ultimately to understand morphogenesis in vivo. 3D cell cultures are already an experimental model of this phenomenon. As a result, there is no reason here to aim for a highly detailed representation of collagen.

Cells transmit forces in the same manner than collagen (orientation is then the orientation of the cytoskeleton). Cells may also exert forces; this is discussed below. However, elements which describe the lumen do not transmit forces.

Collagen changes its orientation when forces are exerted on it, more precisely

 $\stackrel{⃗}{C}\left(x,y,z,t+1\right)=𝜖\stackrel{⃗}{C}\left(x,y,z,t\right)+\frac{\mathit{\Delta t}}{{F}_{1}{\tau }_{c}}\stackrel{⃗}{F}\left(x,y,z,t\right)$

This equation is valid modulo a normalization and a sign factor. Here, $𝜖=±1$ so that the two vectors have the same orientation (a positive dot product) which account for the fact that collagen orientation is defined up to a factor $-1$ ,as mentioned above. We use a relatively small $\mathit{\Delta t}∕{F}_{1}{\tau }_{c}$ ($0.01$) because remodeling is a relatively slow process in comparison with mechanical equilibrium. This equation corresponds to an averaging of the orientation of collagen and of the orientation of the forces at each time step with a weight of the force that depends on its magnitude.

In the case of acini formation in Matrigel, we assume that cells cannot exert forces on the extracellular matrix. In our model, this assumption is equivalent to the assumption that the extracellular matrix cannot transmit forces.

Cells will respond to the stress exerted on it:

$\begin{array}{c}\stackrel{⃗}{F}\left(x,y,z,t\right)={\stackrel{⃗}{F}}_{x}\left(x,y,z,t\right)+{\stackrel{⃗}{F}}_{y}\left(x,y,z,t\right)\\ +{\stackrel{⃗}{F}}_{z}\left(x,y,z,t\right)+{\stackrel{⃗}{F}}_{x}\left(x+1,y,z,t\right)\\ +{\stackrel{⃗}{F}}_{y}\left(x,y+1,z,t\right)+{\stackrel{⃗}{F}}_{z}\left(x,y,z+1,t\right)\end{array}$

Note that the different axes of this stress may be compressive or tensile independently. Also the expression above means that we do not distinguish shear stress from axial stress for the cellular response which is a simplifying approximation.

We also assume that the magnitude of this vector is a relevant quantity for cells:

 ${F}_{m}\left(x,y,z,t\right)=‖\stackrel{⃗}{F}\left(x,y,z,t\right)‖$

#### A.2. Chemical layers

We consider two chemical concentration scalar fields. The first is the “nutrient” layer which is involved in lumen formation. We model this layer by a diffusion equation with fixed boundary conditions (set to $1$), and we simulate it by finite differences. Cells consume these nutrients. When the local concentration in nutrient is below a threshold, the cell dies which leads to lumen formation at the center of structures.

The second chemical component is the inhibitory layer discussed in the text. It is simulated as a chemical that decays rapidly and diffuses. Its effects are thus confined to the neighborhood of its sources. The sources are cells adjacent to a lumen. Because of these specific features, extracellular matrix deposition in the same region would have approximately the same distribution pattern. When this layer is above a threshold in a given element occupied by a cell, both proliferation and motility are impossible for this cell. Further work will aim to distinguish the properties of a chemical inhibitor and the properties of the basement membrane components secreted by epithelial cells.

#### A.3. Cellular proliferation and motility

##### A.3.1. Cell Movement

We model cell-cell contact interactions in a macroscopic manner. We consider ${n}_{\mathit{nei}}\left(x,y,z,t\right)$ the number of neighbors of a cell at position $\left(x,y,z\right)$. Note that we normalize the contribution of every neighbor by distances, so that a cell that is a neighbor with a relative position $\left(1,1,0\right)$ counts for $1∕\sqrt{2}$. We also use the distribution of probability ${p}_{x,y,z,t}$ with

$\begin{array}{c}{p}_{\left(x,y,z,t\right)}\left(i,j,k\right)\propto \delta \left(x+i,y+j,z+k,t\right){\left({n}_{\mathit{nei}}\left(x+i,y+j,z+k,t\right)+\frac{0.5}{\sqrt{{i}^{2}+{j}^{2}+{k}^{2}}}\right)}^{2}\end{array}$

where $\left(i,j,k\right)\in {\left\{-1,0,1\right\}}^{3}-\left(0,0,0\right)$ and $\delta \left(x+i,y+j,z+k,t\right)$ is $1$ if there are no cell in this position or $0$ otherwise. Note that the counts here ignore the cell at position $\left(x,y,z\right)$ since this cell is the object under influence of the other cells in our uses of this distribution. This distribution means that the choice of a direction is isotropic for an isolated cell. However, when there are other cells next to it, a cell will tend to select a direction towards the other cells even though the opposite is still possible.

In the absence of the inhibitor discussed above, cells initiate motion with a probability

 $\mathrm{exp}\left(-\frac{{F}_{m}\left(x,y,z,t\right)}{{f}_{m}}-\frac{{n}_{\mathit{nei}}\left(x,y,z,t\right)}{{n}_{\mathit{inhib}}}\right)$

Here, ${f}_{m}$ is the characteristic magnitude of the inhibition of cell movement by physical forces (we use ${f}_{m}=40$). The inhibition of motility by cell to cell contact is described by ${n}_{\mathit{inhib}}$ (we use ${n}_{\mathit{inhib}}=3$). Then, the possible motion follows a random direction $\stackrel{⃗}{r}$ with the probability distribution ${p}_{\left(x,y,z,t\right)}$. This motion leads to the new position $\left({x}^{\prime },{y}^{\prime },{z}^{\prime }\right)$. The motion is performed if the following condition is met:

 $\frac{{n}_{\mathit{nei}}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)+0.5}{{n}_{\mathit{nei}}\left(x,y,z,t\right)+0.5}\prod _{\stackrel{⃗}{v}}{\left(\frac{|\stackrel{⃗}{r}.\stackrel{⃗}{v}|+0.25}{||\stackrel{⃗}{v}||+0.25}\right)}^{1∕3}-X>0$

where $X$ is a random variable with uniform distribution over $\left[0,1\right]$ and $\stackrel{⃗}{v}$ is $\stackrel{⃗}{C}\left(x,y,z,t\right)$,$\stackrel{⃗}{C}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)$ and $\stackrel{⃗}{F}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)$. We thus assume that movement is facilitated towards positions where the number of neighboring cells is larger than in the initial position. Positions which are easier to access correspond to a displacement in the direction of collagen orientation and local forces. We combine these last effects by considering their geometric average.

##### A.3.2. Proliferation

Every cell has an internal variable $𝜃$ which corresponds to the progress of its cell cycle. Cells may only proliferate when $𝜃$ reaches $1$. In the lack of the inhibitor mentioned above, $𝜃$ is incremented by $\mathit{\Delta 𝜃}=\mathit{\Delta t}∕{\tau }_{p}$ (we use $\mathit{\Delta 𝜃}=1∕6$) where $\mathit{\Delta t}$ is the duration of an iteration of the program and ${\tau }_{p}$is the time required for proliferation when there are no constraints (maximal proliferation).

When $𝜃$ reaches $1$, the cell attempts to proliferate with a probability $\mathrm{exp}\left(-{F}_{m}\left(x,y,z,t\right)∕{f}_{p}\right)$ where ${f}_{p}$ is a characteristic magnitude of the inhibition of proliferation by physical forces (we use ${f}_{p}=80$). The direction in which the cell tries to proliferate is given by a combination of $F\left(x,y,z,t\right)$ and a random vector $\stackrel{⃗}{r}$ generated on the basis of the probability distribution ${p}_{x,y,z,t}$ described for motility. This random vector leads to the position $\left({x}^{\prime },{y}^{\prime },{z}^{\prime }\right)$ that is a neighbor of the initial position. More precisely the direction of the proliferation attempt is given by:

 ${𝜖}_{1}\stackrel{⃗}{F}\left(x,y,z,t\right)+{𝜖}_{2}{w}_{c}\stackrel{⃗}{C}\left(x,y,z,t\right)+{w}_{r}{n}_{\mathit{nei}}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)\stackrel{⃗}{r}$

Let us recall that $C$ in an element occupied by a cell represents the orientation of the cytoskeleton. ${w}_{c}$ is the impact of the cytoskeleton orientation on the direction of proliferation; we use ${w}_{c}=1$. ${w}_{r}$ is the weight of the random factor. Of course, all these weights are relative to the magnitude of the force $F$. ${𝜖}_{i}$ are random coefficients with values which are either $1$ or $-1$. If the resulting direction in which the new cell would appear is occupied, the process is aborted and a proliferation attempt is performed at the next time step. Cells never stop trying to proliferate, but they will stop proliferating when strongly constrained.

#### A.4. Cellular forces

Every cell has an internal variable $\stackrel{⃗}{\varphi }\left(t\right)$ (a normalized vector) which corresponds to the retention of the direction of the former forces it exerted. This vector corresponds to the orientation in which the cell “chooses” to preferentially pull or push. The new orientation is:

$\begin{array}{c}\stackrel{⃗}{\varphi }\left(t+1\right)=10\stackrel{⃗}{\varphi }\left(t\right)+2X+\stackrel{⃗}{F}\left(x,y,z,t\right)\\ -0.5{F}_{m}{\left(x,y,z,t\right)}^{2}\stackrel{⃗}{F}\left(x,y,z,t\right)+\sum _{\stackrel{⃗}{v}\in A}\stackrel{⃗}{v}.\stackrel{⃗}{C}\left(x,y,z,t\right)\stackrel{⃗}{v}\end{array}$

where the sum is over the set $A$ of directions to neighbor positions which are occupied by cells. $X$ is a random vector of norm $1$. Note that cells react to forces by amplifying small forces and opposing strong forces which we model by the polynomial response above.

The force that a cell at position $\left(x,y,z\right)$ can exert on an adjacent element at $\left(x,y,z\right)+\stackrel{⃗}{v}=\left({x}^{\prime },{y}^{\prime },{z}^{\prime }\right)$ is:

 $\left(|\stackrel{⃗}{v}.\stackrel{⃗}{C}\left(x,y,z,t\right)|+|\stackrel{⃗}{C}\left(x,y,z,t\right).\stackrel{⃗}{C}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)|\right)F\stackrel{⃗}{v}$

where $F$ is either ${F}_{\mathit{cell}}$ if the adjacent position $\left({x}^{\prime },{y}^{\prime },{z}^{\prime }\right)$ is a cell or ${F}_{\mathit{col}}$ if this position is occupied by collagen. We use, for example, ${F}_{\mathit{cell}}=4$ and ${F}_{\mathit{col}}=3$.

Then the force that the cell attempts to exert is finally:

 $\left(|\stackrel{⃗}{v}.\stackrel{⃗}{C}\left(x,y,z,t\right)|+|\stackrel{⃗}{C}\left(x,y,z,t\right).\stackrel{⃗}{C}\left({x}^{\prime },{y}^{\prime },{z}^{\prime },t\right)|\right)F\stackrel{⃗}{v}.\stackrel{⃗}{\varphi }\left(t\right)\stackrel{⃗}{v}$

The forces exerted by cells are added to the forces present in the collagen gel, and the new mechanical system is then simulated. Note that the force effectively exerted may be different because the mechanical simulation of the system leads to the elimination of systemic inconsistencies. For example in our model it is not possible for the cell to exert forces on the lumen.

### Video

#### Video 1. Projections of breast epithelial cells seeded in a fibrilar matrix.

Cells emit projections in all directions soon after seeding. These cell projections are involved in collagen organization.

#### Video 2. Breast epithelial cells forming an acinus in a non-fibrilar matrix at day 4.

Cells display limited motility and emit only short projections into the matrix. Cells rotate and divide resulting in the formation of an acinus, a sphere with a central lumen.

#### Video 3. Branching duct at day 7 of culture.

A cell detaches from the main structure and is incorporated back into the structure.

#### Video 4. Collagen fibers and breast epithelial structures after 6 days in culture.

Cells organize collagen in a collagen only matrix and the collagen bundles (green) facilitate the merging of epithelial structures.

#### Video 5 and 6. Formation of a duct by the mathematical model.

Cells are represented in green and the ductal lumen in orange in 3 dimensions. Cells organize collagen over time by exerting forces; in turn, collagen constrains cell proliferation and motility. These interactions lead to the emergence of a main direction of growth.

#### Video 7 and 8. Formation of branching ducts by the mathematical model.

Cells are represented in green and the ductal lumen in orange in 3 dimensions. Ducts branch spontaneously in our model (see also the end of video 5 and 6). For these simulations, we reduced the range of the inhibitor which increases the odds of branching.

#### Video 9. Formation of an acinus by the mathematical model.

Cells are represented in green and the lumen in orange in 3 dimensions. Here, the interactions between the ECM and the cells were removed in order to mimic the effect of Matrigel. In this condition, cells proliferate and move in an isotropic manner leading to the formation of a rounded structure.