Segregation of co-cultured multicellular systems: review and modeling consideration

Abstract

Cell segregation caused by collective cell migration (CCM) is crucial for morphogenesis, functional development of tissue parts, and is an important aspect in other diseases such as cancer and its metastasis process. Efficiency of the cell segregation depends on the interplay between: (1) biochemical processes such as cell signaling and gene expression and (2) physical interactions between cells. Despite extensive research devoted to study the segregation of various co-cultured systems, we still do not understand the role of physical interactions in cell segregation. Cumulative effects of these physical interactions appear in the form of physical parameters such as: (1) tissue surface tension, (2) viscoelasticity caused by CCM, and (3) solid stress accumulated in multicellular systems. These parameters primarily depend on the interplay between the state of cell–cell adhesion contacts and cell contractility. The role of these physical parameters on the segregation efficiency is discussed on model systems such as co-cultured breast cell spheroids consisting of two subpopulations that are in contact. This review study aims to: (1) summarize biological aspects related to cell segregation, mechanical properties of cell collectives, effects along the biointerface between cell subpopulations and (2) describe from a biophysical/mathematical perspective the same biological aspects summarized before. So that overall it can illustrate the complexity of the biological systems that translate into very complex biophysical/mathematical equations. Moreover, by presenting in parallel these two seemingly different parts (biology vs. equations), this review aims to emphasize the need for experiments to estimate the variety of parameters entering the resulting complex biophysical/mathematical models.

Introduction

The process of cell segregation via collective cell migration (CCM) is an integral part of morphogenesis and cancer metastasis (Batlle and Wilkinson, 2012; Carey et al., 2013; Barriga and Mayor, 2019; Devanny et al., 2021). Since these tumor spheroids are very heterogeneous cellular systems, to shed light on them we focus on different breast cell subpopulations with different degrees of mesenchymal character, which can be cultured together. Note that the mesenchymal character of cells is characterized by elongated cell shape, increased migratory cell ability, establishment of front-rear cell polarity, and weakening of cell–cell adhesion. In contrast, the epithelial character of cells is characterized by cuboidal shape, reduced cell mobility, apical-basal polarity, and establishment of strong E-cadherin-mediated cell–cell adhesions (Gandalovičová et al., 2016). The epithelial-to-mesenchymal transition (EMT) has been recognized as one of the hallmarks of breast cancer metastasis (Wang and Zhou, 2013). However, this EMT is not always complete, with cells inside the tumor spheroids having various degrees of mesenchymal character. For example, in a recent study, Devanny et al. (2021) studied six breast cell lines (five of which being cancerous and one being benign) and ordered them from the more epithelial cells (MCF-10A, ZR-75-1) to the more mesenchymal cells (MDA-MB-436, MDA-MB-231). These experimental results suggest the need to better understand the interactions between different cells with different mesenchymal character if we want to understand the segregation process.

While the epithelial-like cells establish strong E-cadherin-mediated adherens junctions (AJs), mesenchymal-like cells establish weak AJs (Devanny et al., 2021). Epithelial-like cells do not establish focal adhesions (FAs) with the extracellular matrix (ECM) already present within the spheroid, while the mesenchymal-like cells establish cell-ECM FAs (Devanny et al., 2021). Changes in the strength of cell–cell adhesion contacts provoke various biochemical mechanisms which influence gene expression and have a feedback impact on cell movement (Barriga and Mayor, 2019). The strength of cell–cell and cell–matrix adhesion contacts is also influenced by transport of the products of cellular metabolism through glycocalyx between neighbor cells (Oberleithner and Hausberg, 2007; Oberleithner et al., 2011). The glycocalyx is a dense gel-like structure made by glycoproteins, glycosaminoglycans, proteoglycans, and associated plasma proteins. This gel-like structure can promote integrin clustering and influence membrane bending (Chighizola et al., 2022). The interactions of cadherin and integrin with the glycocalyx have a feedback on the cellular metabolism. The strength of cell–cell adhesion contacts also influences the mode of cell movement. While epithelial-like cells migrate in the form of strongly connected cell clusters, mesenchymal-like cells migrate in the form of weakly connected cell streams (Foty et al., 1996; Clark and Vignjevic, 2015; Pajic-Lijakovic and Milivojevic, 2021a). Consequently, the movement of mesenchymal-like cells is primarily dissipative and corresponds to viscoelastic liquids (Foty et al., 1996). In contrast to the mesenchymal-like cells, the movement of epithelial-like cells corresponds to viscoelastic solids and induces additional cell residual stress accumulation (Pajic-Lijakovic and Milivojevic, 2021a, 2021b). Accumulated cell residual stress within the epithelium induces an increase in cell packing density which reduces further cell movement (Trepat et al., 2009; Pajic-Lijakovic and Milivojevic, 2021b). The cell movement reduction leads to a change in the state of viscoelasticity and can result in the cell jamming state transition (Nnetu et al., 2013; Notbohm et al., 2016). While epithelial-like cells undergo the jamming state transition, mesenchymal-like cells avoid cell jamming (Grosser et al., 2021). Consequently, during the process of cell segregation, only a part of epithelial cells migrates, while the other part stays in the resting state caused by the cell jamming (Nnetu et al., 2013; Pajic-Lijakovic and Milivojevic, 2021a). Migrating and resting parts of epithelial-like subpopulation show significant difference in the physical parameters such as cell stiffness and the surface tension. Contractile (migrating) epithelial-like cells are much stiffer than noncontractile (resting) ones due to an accumulation of contractile energy (Kollmannsberger et al., 2011; Schulze et al., 2017; Pajic-LIjakovic and MIlivojevic, 2022b). Besides cell stiffness, migrating epithelial-like cells have higher surface tension than the resting epithelial-like cells (Devanny et al., 2021). It is in accordance with the fact that cell contractility enhances the strength of E-cadherin-mediated AJs (Devanny et al., 2021). In contrast to the epithelial-like cells, the contractility induces repulsions among mesenchymal-like cells, which leads to a decrease in the surface tension (Devanny et al., 2021). Based on significant difference in the physical parameters such as the cell stiffness and the surface tension between contractile (migrating) and noncontractile (resting) epithelial cells, the epithelial-like subpopulation can be treated as two-phase systems (Pajic-Lijakovic and Milivojevic, 2019, 2020a). However, the mesenchymal-like cells can be treated as a monophase system.

Batlle and Wilkinson (2012) pointed to three types of mechanisms which underlie the process of cell segregation and formation of the biointerface between the different cell subpopulations. These mechanisms are related to cell signaling, strength of cell–cell adhesion contacts, and cell activation by inducting the actomyosin contractility. When multiple cell populations are co-cultured, signaling from one cell subpopulation influences assembly and contraction of the biointerface-associated actomyosin in the adjacent cells of the other subpopulation (Batlle and Wilkinson, 2012; Lee et al., 2012; Heine et al., 2021; Senigagliesi et al., 2022). A good example is the widely examined co-cultured breast MCF-10A/MDA-MB-231 system, where proteins such as vinculin, laminin-5, and fibronectin secreted by cells are involved in heterotypic interactions between MCF-10A and MDA-MB-231 cells that influence cell adhesion, movement, and can induce EMT of epithelial MCF-10A cells (Bateman et al., 2010; Nikkhah et al., 2011). The strength of cell–cell adhesion contacts is related to the difference in cadherin-mediated adhesion. Cadherins are also involved in a variety of cellular processes, including polarity and gene expression (Barriga and Mayor, 2019). Cumulative effects of homotypic cell–cell interactions along the monocultured cellular surfaces in contact with surrounding liquid medium influence the macroscopic tissue surface tension, while the cumulative effects of heterotypic cell–cell interactions along the biointerface established between two cell subpopulations within co-cultured systems influence the macroscopic interfacial tension (Pajic-Lijakovic and Milivojevic, 2023; Pajic-Lijakovic et al., 2023a). These two physical parameters, tissue surface tension and interfacial tension, depend on the type and strength of cell–cell adhesion contacts, cell contractility (Devanny et al., 2021), and extension or compression of multicellular system caused by CCM (Guevorkian et al., 2021). While macroscopic tissue surface tension has been measured by cell aggregate uniaxial compression between parallel plates (Mombach et al., 2005) and by cell aggregate micropipette aspiration (Guevorkian et al., 2021), the macroscopic interfacial tension between two cell subpopulations has not been measured yet. Deeper insight into interplay between tissue surface tensions of cell subpopulations accompanied by the interfacial tension between them, as well as the viscoelasticity caused by CCM is necessary for understanding the cell segregation within co-cultured systems.

Despite extensive research devoted to study the segregation of various co-cultured systems, we still do not understand the role of tissue surface tension, interfacial tension, and viscoelasticity in cell segregation. One reason is that the co-cultured systems are simplifications of the complex in vivo systems. Nevertheless, the co-cultured systems, under in vitro conditions, can be used to extract some important variables for the biophysical models and eventually parameterize these models. To this end, we focus on tumor breast cell spheroids since: (1) they are the most experimentally investigated tumor multicellular units, usually co-cultured as a free-floating aggregate (although it can be co-cultured also with ECM); (2) they can be considered as toy systems for investigating mixtures of cells with different biophysical properties, which can also be used to parameterize the corresponding biophysical models. It is necessary to identify the forces which influence cell segregation within co-cultured breast cell spheroids, and how they depend on the physical parameters such as surface tension of each subpopulation, interfacial tension between them, interfacial tension gradients, viscoelasticity caused by CCM, and solid stress accumulated within the spheroid core region and (2) show how we can include all these forces and physical parameters in a biophysical model. A more complex spherical multicellular unit is represented by the tumor organoids, which require ECM and a cocktail of growth factors for culture, and which resembles the original tumor tissue (Gunti et al., 2021). However, before we understand the segregation of different cells with epithelial and mesenchymal properties inside the more complex tumor organoids, we need to better understand cell segregation inside the simpler tumor spheroids. For this reason, throughout this study we ignore any discussion about the ECM and growth factors, and focus only on cell–cell interactions.

The focus of this review is to discuss the interrelations between viscoelasticity caused by CCM and surface properties of epithelial-like and mesenchymal-like subpopulations in contact, which influence cell segregation. To this end, we start by discussing the morphology of monocultured breast cell spheroids, and then we focus on co-cultured cell spheroids and the segregation process that takes place between different cell types during CCM. Particular attention is given to the role of different types of mechanical stresses to cell responses (at cell-level and spheroid-level – during collective cell movement). In addition, we point out the importance of accounting for the viscoelasticity in advancing biological physics research and discuss potential opportunities that can be addressed with these tools. We conclude this discussion by presenting a new mesoscopic multiphase model developed by combining some old and new physical models which describe the interplay between viscoelasticity and tissue surface tension in the segregation of co-cultured systems.

Morphology of monocultured breast cell spheroids

Breast cell lines can be characterized based on the degree of their mesenchymal character (Devanny et al., 2021). High degree of mesenchymal character represents the characteristic of the triple negative cancer cell lines (Dai et al., 2017). Devanny et al. (2021) ranged breast cell lines from more epithelial (MCF-10A, ZR-75-1) to more mesenchymal (MDA-MB-436, MDA-MB-231) in character. These different breast cell lines also have different morphological group characteristics based on their cell packing density, the type of cell adhesion contacts, and cell shapes. Strength of cell–cell adhesion contacts and cell packing density influence the surface characteristics of multicellular systems and viscoelasticity caused by CCM which are relevant for further modeling consideration.

These group characteristics can be summarized as: (1) mass, (2) grape-like, and (3) stellate (Kenny et al., 2007; Devanny et al., 2021), which are shown in Figure 1. For example, mass morphology (MCF-10A, MCF-7, ZR-75-1) is characterized by cell close packing, while stellate morphology (MDA-MB-157, MDA-MB-231, MDA-MB-436) corresponds to loose packing structures (Devanny et al., 2021). Grape-like morphology (MDA-MB-468, MDA-MB-361, MDA-MB-453) is between them. The cell types which form mass cell spheroids are capable of establishing strong cadherin-mediated cell–cell AJs. In contrast, the main characteristics of cells which establish stellate morphology is $ {\beta}_1 $-integrin-mediated FAs between cells and ECM. Cadherin-mediated AJs and integrin-mediated FAs have cooperative interrelation (Zuidema et al., 2020). The cell types that form grape-like cell spheroids have low levels of E-cadherin and can establish AJs accompanied by the $ {\beta}_1 $-integrin-mediated FAs (Kenny et al., 2007). The cell types that form stellate spheroids rearrange primarily using $ {\beta}_1 $-integrin-mediated FAs (Kenny et al., 2007). Cell types that form stellate morphology are elongated, while the cell types which correspond to mass and grape-like morphologies are more round (Devanny et al., 2021).

Figure 1. Morphologies of monocultured breast cell spheroids.

Besides cell shape, cell packing density, and the state of cell–cell and cell–matrix adhesion contacts, these morphological groups can be characterized by various mechanisms of cell migration. The mechanism of cell migration influences cell rearrangement and cell response under stress conditions. Two migration mechanisms were observed within breast cell collectives: pressure-driven bleb-like protrusions and filopodia-like protrusions (actin-polymerization-based protrusions). Srivastava et al. (2020) pointed out that bleb-like mechanism of cell movement often exists in tissues in which cells are under mechanical stress. The stress induces disconnection between the cell membrane and the cytoskeleton. Compressive stress of 100 Pa is sufficient to favor bleb-like movement of Dictyostelium cells (Srivastava et al., 2020). Compressive stress generated during cell movement through the collagen I porous structure is enough to intensify bleb-like movement of cancer cells which form grape-like morphology (Guzman et al., 2020). The MDA-MB-231 cells used blebs to pass through nonadhesive confluent environment (Riehl et al., 2021). However, these cells rather perform filapodia-like movement through adhesive confluent environment such as porous collagen I structure (Guzman et al., 2020). In contrast, MDA-MB-468 cells perform movement with blebs through porous collagen I gel (Guzman et al., 2020).

Note that in this review (and in particular for the biophysical modeling of these co-cultured biological systems) we ignore the biochemical aspects of the in vivo tumors (e.g., hormones, growth factors) or the heterogeneous structure of the breast tissue. These aspects will have to be integrated with the biophysical aspects in a second stage, once we clarify these biophysical aspects.

CCM in co-cultured multicellular spheroids and cell residual stress accumulation

Next, we discuss the segregation of co-cultured multicellular spheroids caused by CCM (CCM). The spheroids are widely used to model 3D cell systems that undergo segregation. CCM occurs under stress conditions. Solid stress is accumulated within a spheroid core region caused by cell division and interaction between the spheroid and the surrounding medium (Kalli and Stylianopoulos, 2018). The accumulated solid stress within a core region of multicellular spheroids is primarily compressive and corresponds to a few kPa (Kalli and Stylianopoulos, 2018). The growth-induced compressive stress within breast, colon, pancreatic, and brain tumors under in vivo conditions corresponds to 0.21–20 kPa (Kalli and Stylianopoulos, 2018). In humans, normal physiological blood flow induces shear stress of 1–5 Pa (Baeyens et al., 2016). Higher shear stress, caused by blood flow, can induce various vascular diseases such as atherosclerosis and aneurysm formation (Cunningham and Gotlieb, 2005). Frictional shear stress of several tens of Pa can be generated at the biointerface between epithelial cells and soft medical implants (Pitenis et al., 2018). The CCM itself also induces the generation of stress, normal, and shear. Compressive stress is generated within a migrating cell cluster and during a collision of cell clusters caused by the uncorrelated motility, while the shear stress is generated at the biointerface between migrating cell clusters and surrounding cells (Pajic-Lijakovic and Milivojevic, 2020a, 2020b, 2021a). The generated stress caused by CCM is an order of magnitude lower than the solid stress. Tambe et al. (2013) measured the distribution of cell normal and shear residual stresses accumulated during free expansion of Madin–Darby canine kidney type II (MDCK) cell monolayers. Both stresses were in the range of 100–150 Pa. Notbohm et al. (2016) pointed out that the maximum accumulated normal stress during the rearrangement of confluent MDCK cell monolayers corresponds to 300 Pa. The accumulated stress caused by 3D CCM could be much higher (Pajic-Lijakovic and Milivojevic, 2019, 2020b). Cells tolerate well normal stress up to a few kPa, while the shear stress of only a few Pa can induce cell shape changes, gene expression, cytoskeleton softening, remodeling of cell–cell and cell–matrix adhesion contacts, and EMT (Flitney et al., 2009).

Response of cells under various stress conditions

The response of cells under various mechanical stress conditions accounts for the interplay between various subcellular processes which influence the rate of cell spreading. Moreover, this response is a multi-scale temporal process. The cytoskeleton remodeling and the change in the state of cell adhesion contacts that lead to a cell shape change occurs on a time scale of minutes (Wottawah et al., 2005). For example, it is known that the epithelial MCF-10A cells in the monocultured spheroid core region are rounder, while the cells in the surface region are elongated (Devanny et al., 2021), and this time scale for changes in cell shape corresponds to cell polarization (Alert et al., 2019) and cadherin turnover time (Lee and Wolgemuth, 2011). A time scale of several tens of minutes corresponds to gene expression (Petrungaro et al., 2019) and cell persistence time (Mc Cann et al., 2010), while the CCM takes place on a time scale of hours. Cell division as a possible cause of cell rearrangement is neglected at the time scale of hours, because it occurs on much longer time scales, that is, days. Consequently, the movement of cells, the resulted mechanical strains (volumetric and shear), and the accumulation of normal and shear residual cell stresses occur at a time-scale of hours, while the relaxation of these stresses occurs at a time-scale of minutes (Pajic-Lijakovic and MIlivojevic, 2020b, 2020c). Cell movement induces successive short-time mechanical stress relaxation cycles under constant strain per cycle, while the strain change occurs at a time-scale of hours. Schematic presentation of short-time stress relaxation cycles caused by CCM, which results in the residual stress accumulation during movement of epithelial-like collectives, is presented in Figure 2.

Figure 2. Cell residual stress accumulation within migrating epithelial-like collectives: schematic presentation.

While mechanical stress enhances the movement of some cell types, it has no effect or reduces the movement of others. For example, the compressive stress of 773 Pa suppresses the movement of MCF-10A and MCF-7 cells (Tse et al., 2012). This stress corresponds to the compressive stress generated within breast tissue by cell growth (Kalli and Stylianopoulos, 2018). The normal mechanical stress of a several hundreds of Pa can be induced by collective movement of epithelial cells (Tambe et al. (2013)). In contrast, this stress is capable of enhancing the movement of highly aggressive 4T1 and MDA-MB-231 cells, as well as, 67NR cells (Tse et al., 2012). Riehl et al. (2020) considered and compared responses of MDA-MB-231 and MDA-MB-468 cancer cells, as well as MCF-10A epithelial cells under shear stress of 1.5 Pa. The shear stress can be induced by interstitial flow (Riehl et al., 2020) or by collective movement of epithelial cells (Tambe et al. (2013)). While the shear stress stimulates movement of the MDA-MB-231 cells along the flow direction, this stress has no impact on the movement of MDA-MB-468 cells, and even reduces the movement of MCF-10A epithelial cells. Higher compressive stress can suppress movement of epithelial MCF-10A cells and induce the cell jamming state transition (Grosser et al., 2021; Pajic-Lijakovic and Milivojevic, 2021b). Consequently, the behavior of co-cultured multicellular systems should be considered in the context of cell mechano-sensitivity.

Several parameters have been discussed in the context of cell mechano-sensitivity such as: (1) single-cell stiffness (Tse et al., 2012), (2) level of E-cadherin (Mohammed et al., 2021), and (3) the mechanism of cell movement (Guzman et al., 2020; Devanny et al., 2021). Tse et al. (2012) reported that stiffer cells are less mechano-sensitive. Rudzka et al. (2019) and Riehl et al. (2021) pointed to a correlation between cancer cell stiffness and their invasiveness. The average Young modulus decreases from 1.05 kPa for (less invasive) MCF-7 cells to 0.94 kPa for (more invasive) T47D cells and 0.62 kPa for (the most invasive) MDA-MB-231 cells (Omidvar et al., 2014). The stiffness of single cells is a product of cell mechanical and biochemical interactions with their surroundings which provoke various internal molecular mechanisms within single cells themselves responsible for their adaptation such as the rearrangement of their cytoskeletons and change of the strength of cell–cell adhesion contacts. Yousafzai et al. (2016) found that neighboring cells significantly alter cell stiffness. The MDA-MB-231 cells become stiffer when they are in monocultured systems, while HBL-100 and MCF-7 exhibit softer character.

Mohammed et al. (2021) pointed to E-cadherin-mediated cell–cell adhesion contact as the main cause of the MCF-10A epithelial cell movement reduction and the jamming state transition under the compressive stress. While MDA-MB-231 cells keep filopodia-like movement through dense adhesive surrounding such as collagen I gel, MDA-MB-468 cells prefer movement with blebs (Guzman et al., 2020). Enhanced motility of MDA-MB-231 cells under stress could be connected with their mechanism of movement and the adaptability of FAs, while the movement with blebs makes MDA-MB-468 cells more resistant (Riehl et al., 2021).

To understand the rearrangement of co-cultured cellular systems, in addition to cell response to various stresses it is also necessary to emphasize the role of surface tension of cell pseudo-phases (i.e., various cell types in contact) accompanied by the interfacial tension between them in cell rearrangement.

Surface tension: the driving force for segregation of co-cultured cellular systems

Tissue surface tension depends on the interplay between single-cell contractility and the state of AJs (Stirbat et al., 2013; Pajic-Lijakovic et al., 2023c). Stronger cadherin-mediated cell–cell adhesion contacts, characteristic for the epithelial-like cells, lead to the establishment of higher tissue surface tension in comparison with the mesenchymal-like cells (Devanny et al., 2021).

Cell contractility influences the state of AJs and FAs as well as their crosstalk (Zuidema et al., 2020). Devanny et al. (2021) revealed that contractility plays a fundamentally different role in the cell lines in which the rearrangement is driven primarily by integrins (MDA-MB-468, along with MDA-MB-231 and MDA-MB-436 cells) versus by cadherins (MCF-10A). They emphasized that the cell contractility suppression reduces the tissue surface tension of the epithelial MCF-10A cells. The monocultured MCF-10A cell spheroids treated by 20μM blebbistatin (a molecule capable of suppressing cellular contractility) lose the smooth surface and became more irregular in shape (Devanny et al., 2021). Intensive contractility of surface cells enhances the strength of E-cadherin-mediated adhesion contacts and on that base increases the surface tension. In contrast, enhanced contractility of mesenchymal-like cells induces an increase in cell–cell repulsions, which reduce the surface tension (Devanny et al., 2021). The surface tension of epithelial MCF-10A spheroids is $ 45\pm 18\;\frac{mN}{m} $, while the surface tension of carcinoma MCF10DCIS.com spheroids is lower and equal to $ 21\pm 9\;\frac{mN}{m} $ (Nagle et al., 2022). The surface tension of carcinoma F9 WT cell spheroids is significantly lower and equal to $ 4.74\pm 0.28\;\frac{mN}{m} $.

Besides the tissue surface tension, the segregation of co-cultured multicellular systems depends also on the interfacial tension between cell pseudo-phases.

The interfacial tension as a product of the interactions between various cell types

Interfacial tension between cell subpopulations in co-cultured systems depends on the surface tensions of the pseudo-phases in contact and the heterotypic cell–cell interactions. These heterotypic interactions account for biochemical and mechanical interactions between cells. Mechanical interactions have an impact on the stress generation at the biointerface between the cell subpopulations depending on the relative velocity between them (Pajic-Lijakovic et al., 2023b). Biochemical interactions include cell signaling and gene expression which have a feedback impact on the state of cell–cell and cell-ECM adhesion contacts, migration speed and persistence. Based on our knowledge, an interfacial tension between the subpopulations has not been measured yet. This parameter could be measured by applying some noninvasive technique such as the resonant acoustic rheometry. This technique, which uses the surface capillary waves, has been used for measurement of the interfacial tension within various soft matter systems (Hobson et al., 2021). We will provide qualitative analysis of the interfacial tension in comparison with the surface tensions of the subpopulations within the biophysical model.

Interactions between cells in co-cultured multicellular systems depend on the way of the system co-culture. Various methods have been used such as hanging drop, nonadherent surface, spinning flask, and rotating vessel methods, as well as, microfluidic, acoustic, water-in water emulsion, and 3D printing methods (Froehlich et al., 2016; Raghavan et al., 2016; Bowers et al., 2020; Chae et al., 2021). The hanging drop method has been widely used and ensure cell aggregation within droplets under the influence of surface tension and gravitational forces, while the nonadherent surface method includes cell seeding on the nonadherent scaffolds such as polyacrylamide hydrogel or agarose hydrogel (Chae et al., 2021). Two cell types can be seeded on the same scaffold (direct co-culture) or on separate scaffolds and then cultivated together (in-direct co-culture) (Jo et al., 2021). The MCF-10A cells undergo the EMT when they are surrounded by the MDA-MB-231 cells within a directly co-cultured cellular system (Jo et al., 2021). The EMT in this case is stimulated by the inability of MCF-10A cells to establish AJs with the same type of cells (Jo et al., 2021). The EMT induces weakening of AJs which leads to a decrease in the surface tension of MCF-10A cell subpopulation that has a feedback impact on cell segregation. The consequence of AJs weakening within mesenchymal cell phenotype is a more intensive cell movement (Barriga and Mayor, 2019). However, in in-directly co-cultured MCF-10A/MDA-MB-231 cellular systems, the MCF-10A cells keep their epithelial phenotype (Jo et al., 2021).

Signaling of one cell type in co-cultured cellular systems influences the movement of the other. The MCF10A cells secrete laminin-5 and fibronectin, which stimulate the FAs of cancer cells such as the MDA-MB-231 cells, which could also stimulate their movement (Nikkhah et al., 2011). The movement of the MDA-MB-231 cells is enhanced when they are surrounded by the MCF10A cell (Lee et al., 2012). Heine et al. (2021) revealed that the migration of MDA-MB-231 cells is increased by the presence of the MDA-MB-436 cells. Small extracellular vesicles derived from MDA-MB-231 promote proliferation and movement in MCF7 cells (Senigagliesi et al., 2022).

In the following, we briefly discuss some experimental studies that connect the surface tension and the various types of segregation phenomena that can be observed in multicellular systems.

Segregation of co-cultured multicellular systems: various scenarios

Several scenarios can be considered in the context of the rearrangement of co-cultured cellular systems: complete segregation, partial segregation, and mixed segregation, as shown in Figure 3:

Figure 3. Segregation of co-cultured breast cell spheroids: various scenarios.

The scenarios of cell segregation can result by two mechanisms: interfacial tension between the subpopulations accompanied by the interfacial tension gradients and the viscoelasticity caused by CCM. The MCF-10A cells form compact spheroids in monocultured systems, while in co-cultured cellular systems these cells perform partial or complete segregation depending on the surface tension of cancer pseudo-phase as well as the interfacial tension between them. The MDA-MB-436/MCF-10A co-cultured cellular systems, as well as the MDA-MB-231/MCF-10A co-cultured cellular systems, perform complete segregation such that MCF-10A cells, which have higher surface tension, reach out the spheroid core region, while the MDA-MB-231 (or MDA-MB-436) cells with much lower tissue surface tension cover the spheroid surface region (Carey et al., 2013; Devanny et al., 2021). However, the MCF-10A cells perform partial segregation in a co-culture with the MDA-MB-468 cells (Devanny et al., 2021). These results can be understood when we keep in mind that the MDA-MB-468 cells rather than the MDA-MB-436 cells and the MDA-MB-231 cells are able to establish to some extent cadherin-mediated AJs, which could be additionally stimulated by the presence of the MCF-10A cells in their surroundings (Kenny et al., 2007). Co-culture of the MDA-MB-468 cells with other cancer cell types, such as the ZR-75-1 cells, also leads to partial segregation (Devanny et al., 2021). The ZR-75-1 cells are capable of establishing both $ {\beta}_1 $ integ0rin and E-cadherin and on that base establish higher surface tension in comparison with MDA-MB-468 cells (Devanny et al., 2021). The MDA-MB-468/MDA-MB-157 co-cultured cellular systems perform partial segregation. It is in accordance with the fact that the MDA-MB-157 lacks E-cadherin but can establish cell–cell adhesion contacts by using other cadherin types, and on that base establish higher surface tension (Devanny et al., 2021). The MDA-MB-468 cells form mixed spheroids with the MDA-MB-436 cells (Devanny et al., 2021).

The segregation of co-cultured cellular systems: modeling consideration

In the literature, there are various biophysical models for cancer spheroids (Alexandri et al., 2013; Chen and Zou, 2018; Urcun et al., 2021), all of which incorporate in various forms some interplay between mechanics and cell rearrangement. However, none of these models explicitly consider the viscoelasticity caused by CCM; in particular, the epithelial cells behave as viscoelastic solids, whereas the mesenchymal cells can be treated as viscoelastic liquids; thus, these different types of cells need to be modeled differently. These biomechanical details are important not only for modeling and engineering tissues, but in the long term they are important for cancer treatment.

In the following, we propose a new biophysical model that takes into consideration all previous biological and biomechanical aspects that seem to be important for the segregation of co-cultured breast cancer spheroids. The aim of this modeling consideration is to make suggestions about the possible interplay of various physical parameters on the cell segregation based on experimental data from the literature rather than providing quantitative or qualitative results. It is in accordance with the fact that exact values of some parameters for considered co-cultured systems do not exist. To this end, we focus on two cell subpopulations: epithelial-like cells and mesenchymal-like cells. The epithelial-like phenotype accounts for various breast cell types capable of establishing stronger E-cadherin-mediated AJs, such as the MCF-10A, MCF-7, or ZR-75-1 cells (Devanny et al., 2021). The mesenchymal-like phenotype accounts for various cell types which establish weak AJs and $ {\beta}_1 $ integrin-mediated FAs, such as MDA-MB-157, MDA-MB-231, and MDA-MB-436 cells (Devanny et al., 2021). Since the epithelial-like cell subpopulation has higher surface tension compared to the mesenchymal one, we can treat the epithelial-like subpopulation as a dispersed subpopulation, while the mesenchymal subpopulation can be treated as a continuous cell subpopulation. Due to differences in the surface tension and stiffness between different cells, the epithelial subpopulation can be considered as a two-phase system formed of migrating cells (one phase) and resting cells (another phase). The presence of resting cells in the epithelial-like subpopulation represents a consequence of the jamming state transition discussed before.

The viscoelasticity caused by CCM is the key factor responsible for this difference in rearrangement of epithelial-like and mesenchymal-like subpopulations, which will be discussed below more in detail. Migrating and resting parts of epithelial-like subpopulation can be characterized by different values of physical parameters such as surface tension and cell stiffness. The surface tension of migrating epithelial-like cells (denoted by $ {\gamma}_e^m $) is higher than the one for the resting epithelial-like cells ( $ {\gamma}_e^r $), that is, $ {\gamma}_e^m>{\gamma}_e^r $ (Devanny et al., 2021). It is in accordance with the fact that contractility of epithelial like cells enhance the strength of E-cadherin-mediated AJs (Devanny et al., 2021). The surface tension of mesenchymal-like cells (denoted by $ {\gamma}_c $) is the lowest, that is, $ {\gamma}_c\ll {\gamma}_e^r<{\gamma}_e^m $.

Besides the difference in the surface tension between the migrating and resting parts of the epithelial-like subpopulation, there are also differences in cell stiffness. Migrating epithelial-like cells are much stiffer than resting ones due to an accumulation of the contractile energy (Kollmannsberger et al., 2011; Lange and Fabry, 2013; Pajic-Lijakovic and Milivojevic, 2019). Schulze et al. (2017) revealed that the Young’s modulus of contractile MDSC cell monolayer is $ \sim 33.0 \pm 3.0\; kPa $, while the modulus of noncontractile cells is significantly lower and equal to $ \sim 15.6\pm 5.5\; kPa $. Lange and Fabry (2013) reported that muscle cells can change their elastic modulus by over one order of magnitude from less than 10 kPa in a relaxed (resting) state to around 200 kPa in a fully activated (contractile) state. Mesenchymal-like breast cancer cells are softer than resting epithelial cells (Lekka, 2016).

To model the segregation of co-cultured cancer systems, and to understand the impact of biophysical forces on such cell segregation, it is necessary to account for the rearrangement of epithelial and mesenchymal cells due to cell–cell mechanical and biochemical interactions under simplified in vitro conditions. Understanding all mechanisms behind the segregation of different cell types in vitro is a first step toward the understanding of the complex mechanisms behind the segregation of cell types in vivo. Consequently, the mesoscopic model formulated below describes the rearrangement of three pseudo-phases: (1) migrating, epithelial-like pseudo-phase, (2) resting, epithelial-like pseudo-phase, and (3) mesenchymal-like pseudo-phase. These phase assumptions are in accordance with the fact that only a migrating part of epithelial cells actively contribute to the process of segregation, while the resting part of epithelial cells is in the jamming state.

Cell segregation: the mesoscopic modeling consideration

The phase model is formulated to describe the role of biophysical factors in cell segregation occurred in heterogeneous intra-tumor and inter-tumor cellular systems by considering simplified in vitro model systems such as co-cultured spheroids. The inter- and intra-tumor heterogeneity accounts for the coexistence of cell subpopulations that differ in their genetic, phenotypic, or behavioral characteristics within a given primary tumor, and between a given primary tumor and its metastasis (Martelotto et al., 2014). The phenomenon can be caused by (1) biochemical factors such as genetic and epigenetic factors, as well as, fluctuation in signaling pathways (Marusyk et al., 2012) and (2) physical factors such as viscoelasticity and biomechanics of cell subpopulations, as well as, the dynamics at the biointerface between them (Runel et al., 2021; Pajic-Lijakovic et al., 2023b). The phase model, which describes the rearrangement of these three-phase cellular systems, is formulated based on the biological assumptions discussed previously in Sections “Introduction, Morphology of monocultured breast cell spheroids, CCM in co-cultured multicellular spheroids and cell residual stress accumulation” and it is briefly presented in Figure 4. Graphical presentation of the biophysical model is necessary, in order to emphasize clearly how the different components of the model discussed above interact with each other. The blue arrows indicate model parameters which are included into balances of the volume fractions of the pseudo-phases and force balances. These two parts of the model represents a system of equations. The dark red arrow points out that the cell residual stress is included into the force balances. The cell residual stress is complex parameter which consists of several contributions. The one of the contributions represents a cell residual stress generated by CCM, which is discussed from the rheological standpoint and indicated by the purple arrow. The meaning of the variables and parameters is as introduced throughout this Section “The segregation of co-cultured cellular systems: modeling consideration” and Appendices A and B. For the purpose of formulation, the phase model, it is necessary to distinguish the pseudo-phases and characterize their radial distribution within the spheroid. One phase represents the mesenchymal-like cell subpopulation, while the epithelial-like subpopulation consists of migrating and resting cell pseudo-phases. The pseudo-phases are inhomogeneously distributed within the spheroid volume. The volume fractions of these different phases (migrating epithelial $ {\phi}_{em} $, resting epithelial $ {\phi}_{rm} $, and mesenchymal $ {\phi}_c $) satisfy the following:

(1) $$ {\phi}_{em}\left(r,\tau \right)+{\phi}_{rm}\left(r,\tau \right)+{\phi}_c\left(r,\tau \right)=1 $$

where $ r $ is the radial distance in the spheroid and $ \tau $ is the long-time scale (a time scale of hours). The segregation of epithelial and mesenchymal subpopulations caused by CCM occurs on a time scale of hours $ \tau $. On this time scale, the spheroid can be treated as a canonical ensemble such that the total number of cells inside the spheroid is constant. These cell volume fractions relate to the spheroid volume through the following equations:

• • Migrating epithelial cells: $ {\phi}_{em}\left(r,\tau \right)=\frac{1}{dV}{\sum \limits}_{i=1}^{{N_r}^{\ast }}\Delta {V}_{em\;i} $, where $ dV=4{r}^2\pi dr $ is the volume increment of the spheroid, $ \Delta {V}_{em\;i} $ is the volume of the ith migrating epithelial cluster, $ {N_r}^{\ast } $ is the number of migrating epithelial clusters within the spheroid’s volume increment $ dV $.

• • Resting epithelial cells: $ {\phi}_{rm}\left(r,\tau \right)=\frac{1}{dV}{\sum \limits}_{l=1}^{{N_r}^{\ast \ast }}\Delta {V}_{er\;l} $, where $ \Delta {V}_{er\;l} $ is the volume of the lth resting epithelial cluster, $ {N_r}^{\ast \ast } $ is the number of resting epithelial clusters within the spheroid’s volume increment $ dV $.

• • Mesenchymal cells: $ {\phi}_c\left(r,\tau \right)=\frac{1}{dV}{\sum \limits}_{k=1}^{{N_r}^{\ast \ast \ast }}\Delta {V}_{ck} $, where $ \Delta {V}_{ck} $ is the volume of the kth cluster of mesenchymal cells, $ {N_r}^{\ast \ast \ast } $ is the number of mesenchymal clusters within the spheroid’s volume increment $ dV $.

Figure 4. The schematic presentation of the biophysical model.

Consequently, the spheroid volume increment can be expressed as: $ dV={\sum \limits}_{i=1}^{{N_r}^{\ast }}\Delta {V}_{em\;i}+{\sum \limits}_{l=1}^{{N_r}^{\ast \ast }}\Delta {V}_{er\;l}+{\sum \limits}_{k=1}^{{N_r}^{\ast \ast \ast }}\Delta {V}_{ck} $. Three biointerfaces are taken into consideration: (1) the “c-me” biointerface (i.e., the biointerface between migrating epithelial pseudo-phase and mesenchymal pseudo-phase), (2) the “c-re” biointerface (i.e., the biointerface between resting epithelial pseudo-phase and mesenchymal pseudo-phase), and (3) the “re-me” biointerface (i.e., the biointerface between migrating and resting epithelial pseudo-phases).

Accordingly, with the fact that only migrating epithelial cells contribute to the cell segregation, it is necessary to formulate the phase model which accounts for the jamming state transition of epithelial cells and vice versa influenced by the epithelial cell residual stress accumulation and presence of mesenchymal-like cells in their surroundings. The schematic presentation of the biophysical model is presented in Figure 4.

The model we developed consists of three interconnected parts. The main part of the model describes the dynamics of cell segregation at the spheroid level, while the dynamics at the level of cell clusters is formulated in the Appendices A and B. The main part of the model describes the interplay between the changes of the pseudo-phase local volume fractions and force balance equations for the migrating epithelial-like pseudo-phase and mesenchymal pseudo-phase depending on cell signaling, surface characteristics of the pseudo-phases, and the cell residual stress accumulation. The second part of the model formulates the normal and shear cell residual stresses accumulated within the pseudo-phases. Cell normal stress per pseudo-phase consists of isotropic and deviatoric parts. The isotropic part represents the consequence of the surface properties of the pseudo-phases, while the deviatoric part of stress is induced by CCM. Cell shear residual stress represents the product of natural and forced convections. While natural convection is caused by the interfacial tension gradient established at the biointerfaces between the pseudo-phases, the forced convection is induced by CCM (Pajic-Lijakovic et al., 2023a, 2023b). The normal and shear residual stress per pseudo-phases are formulated in the Appendix A. The part of cell residual stress (normal and shear) caused by CCM depends on the state of single cells and cell rearrangement. The state of single cells accounts for cell contractility and the state of AJs, while the cell rearrangement depends on cell packing density and cell velocity. These parameters influence the constitutive behaviors of migrating cell collectives. The cell residual stress caused by CCM is formulated in the Appendix B. The surface tensions of the pseudo-phases and interfacial tensions between them also depends on the strength of AJs and cell contractility. The strength of AJs and cell contractility, which represent a product of cell adaptation to microenvironmental conditions, depend on homotypic and heterotypic cell–cell and cell–matrix interactions, cell mechanotransduction, gene expression, and transport of cell metabolites through glycocalyx (Oberleithner and Hausberg, 2007; Oberleithner et al., 2011; Barriga and Mayor, 2019; Devanny et al., 2021). Consequently, biochemical and physical mechanisms cooperate together and influence the segregation process. Schematic presentation of the role of physical parameters in the segregation process is shown in Figure 5 and discussed in the context of formulated modeling equations.

Figure 5. Schematic presentation of the role of physical parameters in the segregation process.

The main part of the model is formulated based on modified phase model C proposed for thermodynamic systems far from equilibrium (Ala-Nissila et al., 2004) combined with phase models for viscoelastic phase transition proposed by Tanaka (1997). Accordingly, the cellular interactions are expressed by coupling two scalar fields, that is, the volume fraction of migrating epithelial-like pseudo-phase and volume fraction of mesenchymal-like pseudo-phase. The change of the volume fraction of migrating epithelial-like pseudo-phase in the presence of mesenchymal-like pseudo-phase is expressed as:

(2)

where $ {D}_e $ is the effective dispersion coefficient of migrating epithelial-like pseudo-phase, $ F\left({\phi}_{em},{\phi}_c\right) $ is the Ginsburg-Landau free energy functional equal to $ F\left({\phi}_{em},{\phi}_c\right)={F}_e\left({\phi}_{em}\right)+{F}_c\left({\phi}_{em},{\phi}_c\right) $, $ {F}_e\left({\phi}_{em}\right) $ is the part of free energy functional which describes biochemical and mechanical interactions within the migrating epithelial-like pseudo-phase during contact with resting epithelial-like pseudo-phase, $ {F}_c\left({\phi}_{em},{\phi}_c\right) $ is the part of free energy functional which describes the contribution of mesenchymal-like pseudo-phase to the movement of epithelial cells, $ \frac{\delta F\left({\phi}_{em},{\phi}_c\right)}{\delta {\phi}_{em}} $ is the functional derivative equal to $ \frac{\delta F\left({\phi}_{em},{\phi}_c\right)}{\delta {\phi}_{em}}={\mathit{\lim}}_{\varepsilon \to 0}\frac{F\left[\left({\phi}_{em}\left({r}^{\prime },\tau \right),{\phi}_c\right)+\varepsilon \delta \left({r}^{\prime }-r\right)\right]-F\left[\Big({\phi}_{em}\left({r}^{\prime },\tau \right),{\phi}_c\right]}{\varepsilon } $, $ \varepsilon $ is an increment of $ {\phi}_{em} $. (Ala-Nissila et al., 2004), $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{emR}}\left(r,\tau \right) $ represents the residual stress accumulation within migrating epithelial pseudo-phase equal to $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{emR}}\left(r,\tau \right)=\frac{dV}{N_r^{\ast }}{\sum \limits}_{i=1}^{N_r^{\ast }}{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{emR}\;\boldsymbol{i}}\delta \left(r-{r}_i\right) $, $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{emR}\;\boldsymbol{i}} $ is the residual stress accumulated within the ith migrating epithelial cluster within the spheroid volume increment $ dV $, $ \delta \left(\bullet \right) $ is the Dirac delta distribution function, and $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}} $ is the solid stress. The residual stress accumulation within the ith migrating epithelial cluster is formulated in the Appendix A, while the part of the stress caused by CCM is formulated in the Appendix B. Note that the constitutive models in Appendices A and B, which discuss several viscoelastic regimes, depend on three other parameters: $ {n}_j $, the cell packing density of resting epithelial clusters which corresponds to the cell packing density at the jamming state; $ {n}_{mei} $, the packing density of ith migrating epithelial cell clusters; and $ {n}_c $, the cell packing density of mesenchymal cells. Since the jamming state transition induced by the cell residual stress accumulation leads to an increase in cell packing density, we have $ {n}_{mei}<{n}_j $. Moreover, since the mesenchymal-like cells keep moving and avoid cell jamming (Grosser et al., 2021), their local packing density satisfies $ {n}_c<{n}_j $.

Returning to Eq. (2), the total stress within the migrating epithelial pseudo-phase, that is, can suppress movement by decreasing the volume fraction of migrating epithelial cells and can induce the jamming state transition. The free energy $ {F}_e\left({\phi}_{em}\right) $ has been expressed based on modified model by Cohen and Murray (1981) in the form:

(3) $$ {F}_e\left({\phi}_{em}\right)=\int \left[{f}_e\left({\phi}_{em}\right)+\frac{1}{2}{k}_e^{m-r}{\left(\overrightarrow{\nabla}{\unicode{x03D5}}_{em}\right)}^2\right]{d}^3r $$

where $ {f}_e\left({\phi}_{em}\right) $ represents the energetic effect of the volumetric rearrangement of epithelial cells driven by the surface tension of epithelial cells which is simplified as: $ {f}_e\left({\phi}_{em}\right)\approx \frac{1}{2}{k}_e{\phi_{em}}^2 $, $ {k}_e $ is the parameter which accounts for interactions between migrating epithelial cells within a spheroid, and the other terms represent the gradient energy contributions, while $ {k}_e^{m-r} $ is the interaction parameter which accounts for physical interactions between migrating and resting epithelial-like pseudo-phases at the biointerface. The interfacial tension between migrating and resting epithelial-like pseudo-phases can be expressed as: $ {\gamma}_e^{m-r}\left(\tau \right){A}_e^{m-r}\left(\tau \right)=\int \frac{1}{2}{k}_e^{m-r}{\left(\overrightarrow{\nabla}{\phi}_{em}\right)}^2{d}^3r $ (where $ {A}_e^{m-r}\left(\tau \right) $ is the interfacial area between the pseudo-phases). The biointerface is finite with the thickness which is an order of magnitude larger than the size of single cells (Pajic-Lijakovic and Milivojevic, 2020a). Consequently, the gradient of the interfacial tension at the r-m biointerface influence cell movement along the interface from the resting epithelial pseudo-phase to the migrating one, which is expressed by the Marangoni flux (Pajic-Lijakovic and MIlivojevic, 2022c).

The free energy $ {F}_c\left({\phi}_{em},{\phi}_c\right) $ can be expressed as:

(4) $$ {F}_c\left({\phi}_{em},{\phi}_c\right)={\displaystyle \begin{array}{l}\int \Big[{f}_c\left({\phi}_{em},{\phi}_c\right)+\frac{1}{2}{k}_{ec}^m{\left(\overrightarrow{\nabla}{\phi}_{em}\right)}^2\\ {}+\hskip2px \frac{1}{2}{k}_{ec}^r{\left(\overrightarrow{\nabla}\Big(1-{\phi}_{em}-{\phi}_c\right)}^2\Big]{d}^3r\end{array}} $$

where $ {f}_c\left({\phi}_{em},{\phi}_c\right) $ accounts for the cumulative effects of mesenchymal-like cell signaling which influences the volumetric rearrangement of epithelial cells and the second term represents the gradient energy contribution, while $ {k}_{ec}^m $ is the interaction parameter which accounts for physical interactions between migrating epithelial-like pseudo-phase and mesenchymal-like pseudo-phase at the biointerface. The interfacial tension between migrating epithelial-like pseudo-phase and mesenchymal-like pseudo-phase can be expressed as: $ {\gamma}_{ec}^m\left(\tau \right){A}_{ec}^m\left(\tau \right)=\int \frac{1}{2}{k}_{ec}^m{\left(\overrightarrow{\nabla}{\phi}_{em}\right)}^2{d}^3r $ (where $ {A}_{ec}^m\left(\tau \right) $ is the interfacial area between these pseudo-phases), while the interfacial tension between resting epithelial-like pseudo-phase and mesenchymal-like pseudo-phase is $ {\gamma}_{ec}^r\left(\tau \right){A}_{ec}^r\left(\tau \right)=\int \frac{1}{2}{k}_{ec}^r{\left(\overrightarrow{\nabla}\left(1-{\phi}_{em}-{\phi}_c\right)\right)}^2{d}^3r $ (where $ {A}_{ec}^r\left(\tau \right) $ is the interfacial area between these pseudo-phases). The biointerface is finite and the gradient of interfacial tension established at the c-me biointerface stimulates movement of mesenchymal-like pseudo-phase along the biointerface toward to the migrating, epithelial-like pseudo-phase (Pajic-Lijakovic and Milivojevic, 2020a).

The energy $ {f}_c\left({\phi}_{em},{\phi}_c\right) $ can be expressed as:

(5) $$ {f}_c\left({\phi}_{em},{\phi}_c\right)\approx \frac{1}{2}{k}_c{\phi_{em}}^2{\phi_c}^2-\frac{1}{2}{\beta}_c{\phi_c}^2 $$

where $ {k}_c $ is the parameter which accounts for interactions between migrating epithelial cells and mesenchymal cells which arise as a product of cell signaling and $ {\beta}_c $ is the measure of cumulative effects of tractions of mesenchymal-like cells. While cell signaling is capable of enhancing the movement of both pseudo-phases (the first term of the right-hand side of Eq. (5)), the traction of mesenchymal cells can reduce movement of mesenchymal-like cells (the second term of the right-hand side of Eq. (5)). This reduction is pronounced when mesenchymal-like cells establish stronger FAs with the ECM already present within the multicellular spheroid.

The change of the volume fraction of mesenchymal-like pseudo-phase in the presence of migrating, epithelial-like pseudo-phase is expressed as:

(6) $$ \frac{\partial {\unicode{x03D5}}_c\left(r,\tau \right)}{\partial \tau }+\overrightarrow{\nabla}\left({\phi}_c{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{R}}\right)=\overrightarrow{\nabla}\left[{D}_c\left({\unicode{x03D5}}_c\overrightarrow{\mathbf{\nabla}}\frac{\delta F\left({\unicode{x03D5}}_{em},{\phi}_c\right)}{\delta {\unicode{x03D5}}_c}\right)+\overrightarrow{\nabla}\bullet \left({\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}}-{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}}\right)\right] $$

where $ {D}_c $ is the effective dispersion coefficient which quantifies spreading of mesenchymal pseudo-phase, $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{R}} $ is the relative velocity between migrating epithelial and mesenchymal cell pseudo-phases equal to $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{R}}\left(r,\tau \right)={\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right)+{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right) $, $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $ is the local velocity of mesenchymal-like cells equal to $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right)=\frac{dV}{N_r^{\ast \ast \ast }}{\sum}_{k=1}^{N_r^{\ast \ast \ast }}{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}\boldsymbol{k}}\delta \left(r-{r}_k\right) $, $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{ck}} $ is the average velocity of the center of mass of the kth cluster of mesenchymal cells, and $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right) $ is the local velocity of epithelial-like cells equal to equal to $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right)=\frac{dV}{N_r^{\ast }}{\sum}_{i=1}^{N_r^{\ast }}{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}\boldsymbol{i}}\delta \left(r-{r}_i\right) $, $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{ei}} $ is the average velocity of the center of mass of the ith cluster of epithelial cells. This formulation of the relative velocity is in accordance with the fact that the velocities $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}} $ and $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $ have an opposite directions. While mesenchymal cells migrate from the spheroid core region toward its surface driven by the solid stress, epithelial cells migrate from the spheroid surface region toward its core region driven by the surface tension. The velocity of cancer cells $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $ satisfies the condition $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right)\ge {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right) $ (Clark and Vignjevic, 2015), while the velocity of epithelial cells satisfies two conditions: (1) $ 0<{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right)\le {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $ for migrating cells and (2) $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right)=0 $ for resting cells under jamming. Consequently, the relative velocity is $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{R}}\left(r,\tau \right)>{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $. The first term on the right-hand side of Eq. (6) describes interactions between mesenchymal-like and migrating, epithelial-like pseudo-phases which enhance movement of mesenchymal-like cells while the second term describes an influence of the accumulated stress within the spheroid core region equal to $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}}-{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}} $ on movement of the mesenchymal-like cells (where $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}} $ is the residual stress caused by collective movement of mesenchymal-like pseudo-phase). The stress $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}} $ is expressed as $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}}\left(r,\tau \right)=\frac{dV}{N_r^{\ast \ast \ast }}{\sum}_{k=1}^{N_r^{\ast \ast \ast }}{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}\;\boldsymbol{k}}\delta \left(r-{r}_k\right) $, $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{cR}\;\boldsymbol{k}} $ is the residual stress accumulated within the kth mesenchymal cluster within the spheroid volume increment $ dV $, and $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}} $ is the solid stress. The residual stress accumulation within the kth mesenchymal cluster is formulated in the Appendix A, while the part of the stress caused by CCM is formulated in the Appendix B. In order to understand the process of cell segregation, it is necessary to estimate change in velocity for mesenchymal cell pseudo-phase $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right) $ and migrating epithelial pseudo-phase $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right) $.

The force balance for the rearrangement of migrating epithelial pseudo-phase and mesenchymal pseudo-phase

CCM causes mechanical waves generation in the form of oscillatory change in cell velocity, resulted strain and corresponding cell residual stress (Serra-Picamal et al., 2012; Notbohm et al., 2016; Pajic-Lijakovic and Milivojevic, 2020c, 2022a). Oscillatory change of epithelial cell velocity (i.e., effective inertia), in the case of cell segregation, arises as the result of competition between interfacial tension force and the mixing force against viscoelastic force. The viscoelastic force is capable of suppressing movement of epithelial cells (Grosser et al., 2021; Pajic-Lijakovic and Milivojevic, 2023). The interfacial tension force is expressed by modified model proposed by Pajic-Lijakovic and Milivojevic (2020c) as:

(7) $$ {\phi}_{em}\left\langle {n}_{em}\right\rangle {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{st}}}^{\boldsymbol{e}}={\phi}_{em}\left\langle {n}_{em}\right\rangle \left[-{S}_e^{r-m}\hskip0.24em {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{e}}}^{r-m}-{S}_e^{c-m}{{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{e}}}^{c-m}\right] $$

where $ \left\langle {n}_{em}\right\rangle =\frac{{\sum \limits}_{i=1}^{{N_r}^{\ast }}{N}_{em i}}{{\sum \limits}_{i=1}^{{N_r}^{\ast }}\Delta {V}_{em\;i}} $ is the average packing density of single migrating epithelial cluster within the spheroid volume increment $ dV $, $ {N}_{emi} $ is the number of cells within the ith migrating epithelial cluster, $ {S}_e^{r-m}={\gamma_e}^r-\left({\gamma_e}^m+{\gamma_e}^{r-m}\right) $ is the spreading coefficient of migrating epithelial cells in contact with resting epithelial cells (such that $ {S}_e^{r-m}<0 $), $ {S}_e^{c-m}={\gamma}_c-\left({\gamma_e}^m+{\gamma_{ec}}^m\right) $ is the spreading coefficient of migrating epithelial cells in contact with mesenchymal-like cells (such that $ {S}_e^{c-m}<0 $), $ {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{e}}}^{r-m} $ is the local displacement field caused by movement of epithelial clusters near the re-me biointerface, $ {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{e}}}^{c-m} $ is the local displacement field caused by movement of epithelial clusters near the c-me biointerface. The ratios $ {X}_1=\frac{\gamma_c+{\gamma}_{ec}^m}{\gamma_e^m}<1 $ and $ {X}_2=\frac{\gamma_c+{\gamma}_{ec}^r}{\gamma_e^r}<1 $. Maximizing the ratios $ {X}_1 $ and $ {X}_2 $, such that $ {X}_1,{X}_2\to 1 $ by enhancing cell–cell adhesion contacts within the mesenchymal-like pseudo-phase can reduce the spreading of mesenchymal cells. Interfacial tension force acts to reduce the interfacial area between epithelial and mesenchymal cell subpopulations. The viscoelastic force is resistive force and acts to reduce movement of epithelial cells. The mixing force $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{m}}}^{\boldsymbol{e}} $ can be expressed as: $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{m}}}^{\boldsymbol{e}}=-{\phi}_{em}\overrightarrow{\nabla}\frac{\delta F\left({\phi}_{em},{\phi}_c\right)}{\delta {\phi}_{em}} $. The viscoelastic force represents a consequence of inhomogeneously distributed cell residual stress, and is expressed by Murray et al. (1988) and Pajic-Lijakovic et al. (2023a; 2023c) as: $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{Tve}}}^{\boldsymbol{e}}=\overrightarrow{\mathbf{\nabla}}\bullet \left({\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{e}\boldsymbol{mR}}+{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}}\right) $. The corresponding force balance is expressed by Pajic-Lijakovic et al. (2023a; 2023c) as:

(8)

where $ {\left\langle m\right\rangle}_e $ is the average mass of a single epithelial cell and $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}} $ is the epithelial cell velocity equal to $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\left(r,\tau \right)=\frac{d{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{e}}}{d\tau} $, $ \frac{D{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}}{D\tau}=\frac{\partial {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}}{\partial \tau }+\left({\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}}\bullet \overrightarrow{\mathbf{\nabla}}\right){\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{e}} $ is the material derivative (Bird et al., 1960).

While the viscoelastic force reduces movement of epithelial cells, this force represents a driving force for migration of mesenchymal like cells. It is in accordance with the fact that accumulated cell stress enhances movement of cancer cells and suppresses movement of epithelial cells (Tse et al., 2012; Riehl et al., 2020, 2021). Viscoelastic force in this case is expressed by modified model proposed by Pajic-Lijakovic and Milivojevic (2020c) as: $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{Tve}}}^{\boldsymbol{c}}=\overrightarrow{\mathbf{\nabla}}\bullet \left({\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{SD}}-{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{c}\boldsymbol{R}}\Big)-{\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{R}\;\boldsymbol{ECM}}\right) $ (where $ {\overset{\sim }{\boldsymbol{\sigma}}}_{\boldsymbol{R}\;\boldsymbol{ECM}} $ is the residual stress accumulated within ECM). Besides the viscoelastic force, the mixing force and the interfacial tension force also drive extension of mesenchymal-like cells. The mixing force can be expressed as: $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{m}}}^{\boldsymbol{c}}={\unicode{x03D5}}_c\overrightarrow{\mathbf{\nabla}}\frac{\delta {F}_c\left({\unicode{x03D5}}_{em},{\phi}_c\right)}{\delta {\unicode{x03D5}}_c} $. The interfacial tension force is equal to $ {\phi}_c\left\langle {n}_c\right\rangle $ $ {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{st}}}^{\boldsymbol{c}}={\phi}_c\left\langle {n}_c\right\rangle \left[{S}_c^{c-m}\hskip0.24em {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{c}}}^{c-m}+{S}_c^{c-r}{{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{c}}}^{c-r}\right] $ (where $ \left\langle {n}_c\right\rangle =\frac{{\sum \limits}_{k=1}^{{N_r}^{\ast \ast \ast }}{N}_{ck}}{{\sum \limits}_{k=1}^{{N_r}^{\ast \ast \ast }}\Delta {V}_{c\;k}} $ is the average packing density of single mesenchymal cluster within the spheroid volume increment $ dV $, $ {N}_{ck} $ is the number of cells within the kth mesenchymal cluster, $ {S}_c^{c-m}={\gamma_e}^m-\left({\gamma}_c+{\gamma_{ec}}^m\right) $ is the spreading coefficient of the mesenchymal like cells toward the migrating epithelium (such that $ {S}_c^{c-m}>0 $), $ {S}_c^{c-r}={\gamma_e}^r-\left({\gamma}_c+{\gamma_{ec}}^r\right) $ is the spreading coefficient of the mesenchymal like cells toward the resting epithelium (such that $ {S}_c^{c-r}>0 $), $ {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{c}}}^{c-m} $ is the local displacement field caused by movement of mesenchymal-like clusters near the c-me biointerface, and $ {{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{c}}}^{c-r} $ is the local displacement field caused by movement of mesenchymal-like clusters near the c-re biointerface. The ratios $ {X}_3=\frac{\gamma_e^m+{\gamma}_{ec}^m}{\gamma_c}>1 $ and $ {X}_4=\frac{\gamma_e^r+{\gamma}_{ec}^r}{\gamma_c}>1 $ such that $ {X}_3>{X}_4 $. Minimizing the ratios $ {X}_3 $ and $ {X}_4 $ by enhancing cell–cell adhesion contacts within the mesenchymal-like pseudo-phase can reduce the spreading of mesenchymal cells. The oscillatory change of cell velocity in the case of movement of mesenchymal-like cells represents the consequence of the competition between: the viscoelastic force, mixing force, interfacial tension force, against the traction force. The traction force as a consequence of established FAs influences movement of mesenchymal cells, while the epithelial cells do not establish FAs within a spheroid (Devanny et al., 2021). This force is capable of reducing movement of mesenchymal cells depending on the strength of FAs (Fuhrmann et al., 2017). The traction force is expressed as: $ \rho {{\overrightarrow{\boldsymbol{F}}}_{\boldsymbol{tr}}}^{\boldsymbol{c}}=\rho k{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{ECM}} $ (where $ k $ is an elastic constant of single FA, $ \rho $ is the number density of FAs, and $ {\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{ECM}} $ is the displacement field of ECM caused by movement of cancer cells) (Murray et al., 1988). The corresponding force balance can be expressed by modifying the model proposed by Pajic-Lijakovic and Milivojevic (2020c) for 2D CCM as:

(9)

where $ {\left\langle m\right\rangle}_c $ is the average mass of a single cancer cell and $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}} $ is the epithelial cell velocity equal to $ {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\left(r,\tau \right)=\frac{d{\overrightarrow{\boldsymbol{u}}}_{\boldsymbol{c}}}{d\tau} $, $ \frac{D{\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}}{D\tau}=\frac{\partial {\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}}{\partial \tau }+\left({\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}}\bullet \overrightarrow{\mathbf{\nabla}}\right){\overrightarrow{\boldsymbol{v}}}_{\boldsymbol{c}} $ is the material derivative (Bird et al., 1960). While movement of mesenchymal cells corresponds to the convective regime, movement of epithelial-like cells changes from convective regime through conductive regime to the damped-conductive regime.

The efficient process of cell segregation can be postulated based on the proposed modeling consideration. In this purpose, following parameters are introduced, such as: (1) volume fraction of epithelial cells in the resting state $ {\phi}_{er}\to 0 $ and (2) the surface tension ratios $ {X}_1,{X}_2,{X}_3,{X}_4\to 1 $.

An increase in the degree of mesenchymal character of epithelial-like subpopulation leads to a decrease in the volume fraction of cells in the resting (jamming) state, corresponding surface tension, as well as, the interfacial tensions between the pseudo-phases. The surface tension of epithelial-like pseudo-phase corresponds to an order of magnitude from several $ \frac{mN}{m} $ to tens of $ \frac{mN}{m} $ (Mombach et al., 2005), while the surface tension of mesenchymal-like cells is significantly lower (Devanny et al., 2021). Interfacial tensions between the subpopulations are lower than the surface tension of migrating epithelial cells, but correspond to a same order of magnitude with it. The volume fraction of epithelial-like cells in the resting state during the segregation process could be even larger than 15% of whole epithelial-like subpopulation and placed primarily within the spheroid core region. The cell normal and shear residual stress accumulation caused by CCM corresponds to several tens of Pa (Tambe et al., 2013) while the solid stress within the spheroid core region is significantly larger and corresponds to several kPa (Kalli and Stylianopoulos, 2018).

The formulated biophysical model, pointed to the interrelation between the physical parameters responsible for the cell segregation process by accounting for the viscoelasticity of the pseudo-phases and effects along the biointerfaces between them as was presented graphically in Figure 5. Some parameters such as the residual stress accumulation caused by movement of epithelial collectives, solid stress accumulated in the spheroid core region, and tissue surface tension have been already measured, while the others such as the interfacial tension between two cell subpopulations in direct contact and interfacial tension gradient have not been measured yet. Gsell et al. (2023) recently confirmed cell movement along the multicellular surface in contact with liquid medium driven by the tissue surface tension gradient, that is, the Marangoni effect. However, the surface tension gradient itself has not been measured. The maximum cell residual stress caused by collective movement of epithelial monolayers corresponds to a few hundreds of Pa (Serra-Picamal et al., 2012; Notbohm et al., 2016). The solid stress corresponds to a few kPa (Kalli and Stylianopoulos, 2018). Various experimental techniques have been used for the measurement the cell stress caused by CCM such as: (1) the monolayer stress microscopy (Tambe et al., 2013), (2) microbead/droplet-based stress sensors (Campàs et al., 2014; Dolega et al., 2017). Some of these techniques, such as monolayer stress microscopy, are suitable for the measurement of the stress in 2D, while the others such as using droplet-based stress sensors can be used for the measurement of anisotropic normal stresses only. Development of the suitable measuring technique is a prerequisite for the improvement of our knowledge about the cell rearrangement. The tissue surface tension has been measured under equilibrium conditions only (i.e., the static tissue surface tension). Various experimental techniques have been used for the measurement of the static tissue surface tension, such as cell aggregate compression between parallel plates (Mombach et al., 2005; Marmottant et al., 2009), cell aggregate micropipette aspiration (Guevorkian et al., 2021), and magnetic force tensiometer (Nagle et al., 2022). The experimental values of the static tissue surface tension vary from a few $ \frac{mN}{m} $ to several tens of $ \frac{mN}{m} $ depending on cell type and measuring technique. However, the tissue surface tension and the interfacial tension between the two subpopulations are time dependent parameters. This modeling consideration is an effort to stimulate further biological research in this field.

Conclusion

The segregation of co-cultured cellular spheroids made by breast cells which have various degrees of mesenchymal character is considered and discussed in the context of physical parameters such as tissue surface tension, viscoelasticity caused by CCM, and accumulated solid stress within the spheroid core region. Three scenarios are possible: (1) complete segregation, (2) partial segregation, and (3) mixed segregation. To account for these scenarios, we proposed a mesoscopic phase model which accounts for biochemical and physical interactions between two cell subpopulations: epithelial-like subpopulation (i.e., the breast cells with low degree of mesenchymal character) in contact with the mesenchymal-like subpopulation (i.e., the breast cells with high degree of mesenchymal character). The model describes the change of volume fractions of migrating and resting epithelial-like pseudo-phases, as well as, the change in volume fraction of mesenchymal-like pseudo-phase, during the segregation process, as a product of their interactions. These interactions are characterized by interplay between: (1) the surface tensions of pseudo-phases, (2) interfacial tensions between them, (3) interfacial tension gradients, (4) cell residual stress accumulation, and (5) the solid stress accumulated within the spheroid core region. The choice of model parameters is in accordance with the physical and biological properties of the subpopulations in contact.

Contractile (migrating) epithelial-like pseudo-phase has the highest surface tension, while the mesenchymal-like pseudo-phase has the lowest surface tension, which enables their intensive spreading toward the epithelium. Accordingly, with the fact that only migrating, epithelial pseudo-phase actively contributes to the segregation process, the efficiency of the segregation process can be improved by minimizing the volume fraction of epithelial cells in the resting state. The volume fraction of epithelial cells in the resting state can be minimized by reducing the compressive and shear stresses accumulated within the epithelial subpopulation. The main parameters which influence compressive and shear residual stress accumulation within the epithelium are: (1) the interfacial tension between mesenchymal and migrating epithelial pseudo-phases, (2) the corresponding interfacial tension gradient, and (3) characteristics of epithelial cell rearrangement related to the local cell packing density and cell velocity which have a feedback impact on the constitutive behavior of migrating epithelium.

The biophysics multiphase model proposed in this review is only the first step in the modeling of complex multicellular spheroids. Additional experiments are needed in order to measure interfacial tensions among the pseudo-phases as well as the interfacial tension gradients and to correlate them with cell residual stress accumulation. Acquisition of experimental data to parameterize these types of multi-scale and multiphase models is another step in the modeling process. This is particularly difficult since it requires multiple experiments to measure various parameters across different scales: from the local number density of FAs or the elastic constants of single FAs (at molecular level), to the surface tensions for different cell types (at the level of cell clusters), and local displacement fields caused by the movements of different cell types (at tissue level).