Role of senescent cells in the motile behavior of active, non-senescent cells in confluent populations

Characteristics of cell migration in a confluent population depend on the nature of cell-to-cell interactions as well as cell-intrinsic properties such as the directional persistence in crawling. In addition, biological tissues (or cell cultures) almost always carry anisotropies and they too can significantly affect cell motility. In the light of this viewpoint, the emergence of cellular senescences in a confluent population of active cells raises an interesting question. Cellular senescence is a process through which a cell enters a permanent growth-arrest state and generally exhibits a dramatic body expansion. Therefore, randomly emerging senescent cells transform an initially homogeneous cell population to a “binary mixture” of two distinct cell types. Here, using in vitro cultures of MDA-MB-231 cells we investigate how spatially localized cellular senescence affect the motility of active cells within a confluent population. Importantly, we estimate the intercellular surface energy of the interface between non-senescent and senescent MDA-MB-231 cells by combining the analysis on the motile behaviors of non-senescent cells encircling senescent cells and the result of extensive numerical simulations of a cellular Potts model. We find that the adhesion of normal cells to senescent cells is much weaker than that among normal cells and that the ‘arclength’ traveled by a normal cell along the boundary of a senescent cell, on average, is several times greater than the persistence length of normal cell in a densely packed homogeneous population. The directional persistent time of normal cell during its contact with a senescent cell also increases significantly. We speculate that the phenomenon could be a general feature associated with senescent cells as the enormous expansion of senescent cell’s membrane would inevitably decrease the density of cell adhesion molecules.

www.nature.com/scientificreports/ Earlier, it was reported that embryonic cells of two different types that are initially dissociated and randomly mixed can even spontaneously sort to reestablish homogeneous binary phase separations, even when the involved cells do not actively crawl around 20 . Subsequently, the difference in the intercellular adhesive property of the two cell subpopulations, in connection with a Cellular Potts Model (CPM), was attributed to this sorting phenomenon 21 . In the simulation of a CPM, the final configuration of cell populations approaches to a state having the global minimum of the total energy associated with the system, of which intercellular surface energy (which determines the "stickiness" of interface) is one of the main factors affecting the tissue (or confluent cell population) dynamics. In fact, it was suggested that different intercellular surface energies for different interfaces in a binary mixture of cells that led to a complete cell sorting 20,22 . In other cases, different levels of surface energies could lead to different flocking behavior of homogeneous cell populations 19,23,24 .
In this paper, we carefully estimate the intercellular surface energy E sn between normal and senescent MDA-MB-231 cells by combining the result of numerical simulations of a (previously validated) CPM and experimental analysis on the angular speeds ωs of normal cells that move along the boundary of a senescent cell. We find that E sn (~ − 20) is much higher than the intercellular surface energy between normal cells (E nn ~ − 65); that is, normal cells form a much weaker adhesion to senescent cells than normal cells themselves. Consequently, each encounter of a normal cell with a senescent cell body provides a small, curvy trajectory, during which the directional persistent time (and length) has increased significantly.

Motility of normal cells in contact with a senescent cell in experiments.
For a given high-density MDA-MB-231 cell culture, senescent cells were observed to emerge spontaneously from the initial seeding in a random fashion both in space and time (Fig. 1a). Then, the cell bodies (including nuclei) of senescent cells expanded enormously on 2D culture substrate over a few days, and accordingly, their morphologies were often referred to as "fried eggs": See the huge empty spaces, which are guided by yellow dashed lines in Fig. 1a. Few years ago we analyzed this phenomenon in detail and found a surprising role of senescent cells actively recruiting newly replicated non-senescent cells, which show mitotic cell-rounding, towards their body centers to form a small cell cluster (as marked by a white arrow in Fig. 1a) 14 . With the emergence of these senescent cells, we may view the monoclonal cell culture as a binary mixture of two different cell populations: that is, senescent vs. non-senescent, normal cell populations. Apart from their obvious phenotypic differences, we hypothesized that physical as well as biochemical interactions between senescent cells and normal cells would be different from those among the normal cells. Alternatively, the scatters of senescent cells could be viewed as large inhomogeneities embedded in the confluent population of normal cells. Then, a question naturally arises: How the presence of these large inhomogeneities would affect the motility of normal cells, which were known to be super-diffusive 19 . Four exemplary moving traces of normal cells are superimposed in the snapshot image Fig. 1a, and two (orange and yellow) of them clearly show a quite long period of contact time (marked by red color) with the senescent cells in their vicinity.
For the crawling normal cells, the territories occupied by senescent cell bodies were like an island to which normal cells could not invade (except for adjacent replicating cells showing mitotic cell-rounding 14 ): See, for example, the orange-colored cell trace which "wraps" around the senescent cell territory demarcated by a red solid line in Fig. 1a. Interestingly, normal MDA-MB-231 cells appeared to have some preference to stay in contact with senescent cells, and this propensity was rather obvious in a low-density cell culture where cells prefer to form small clusters. Two different examples are highlighted by two (red and blue) bounding boxes in Fig. 1b: Two non-senescent normal cells in contact stay together forming a bounded pair within the small red boxed area, while several normal cells are migrating along the huge senescent cell's boundary inside the blue boxed area. The bounded pair of normal MDA-MB-231 cells generically exhibit a quite robust steady rotation (with a period ~ 100 min) as in the representative example illustrated in Fig. 1c, and the underlying mechanism for the pair rotation was discussed in our (KL and HL) recent work 19 . Briefly, the phenomenon requires (1) an active crawling of the cells with a proper degree of directional persistence (strength and memory) and (2) a good degree of cell-cell adhesion, among other properties. Similarly, a normal cell migrating along the boundary (yellow dashed line) of a large senescent cell is given in Fig. 1d (see Supplementary Video 1); the cell's moving velocity fluctuates during the time course of migration but is quite steady as we will discuss more shortly. The representative example given Fig. 1d clearly suggests that the normal MDA-MB-231 cells have a preference to stay in contact with senescent cells. Then, the main question is how much the preference of cell-cell adhesion has changed due to cellular senescence. With a careful analysis of many crawling trajectories of normal cells moving around senescent cells, which have a wide range of sizes, and those of computer simulations of a Cellular Potts Model, we provide an answer to this question.
A total of 60 cases with a similar configuration (one or a few normal cells moving along the boundary of a senescent cell) was quantified, and the trajectories of normal cells were traced and plotted in Fig. 2a with the centroids of senescent cells positioned at the origin). On average, the motility of non-senescent cells was very steady with only a very few reverse turns observed. The contact between the normal cell and senescent cell forming a pair rotation lasted a minimum of 400 min, with some cases reaching 1800 min. We have decomposed the trajectories into (R, θ) as shown in Fig. 2b. Since the senescent cell body is not an ideal disk defined with a fixed radius, R(t), the distance between the centroid of the normal cell and that of the senescent cell, fluctuates significantly in time.
Nevertheless, a small systematic increase in R over a day or so can be seen in Fig. 2b (top). The angular position θ of a normal cell encircling is also quite noisy but exhibits a clear systematic trend over time (see Fig. 2b, bottom frame). The local slope is nothing but the instantaneous angular speed ω: It also fluctuates a lot in time but has a characteristic value for different senescent cell bodies having a different R as illustrated in Fig. 2c  www.nature.com/scientificreports/ fitting constant of A = 1.79 ± 0.10 μm/min, which is nothing but tangential speed. The wide range of R covered in Fig. 2c simply reflects the diverse sizes of the senescent cells pairing with a non-senescent cell; as such we do not have any control over the range of R. Shown at the bottom of Fig. 2c are the histograms of ω and R of n = 60 different cases. Later, we run our model simulation based on these average values of ω and R . The fairly nice fit of ω ∼ 1 R (black solid line) in Fig. 2c suggests that the instantaneous tangential speed v θ (= Rω = A = 1.79 μm/ min) of moving normal cell along the boundary of a senescent cell is rather independent of the properties (such as size or age) of the senescent cell. Indeed, the ensemble average of v θ is measured to be 1.89 ± 0.53 μm/min which is off by 5.59% (see Fig. 2d). In other words, we would expect no significant curvature effect at play influencing the value of v θ , meaning that v θ is very much an inherent motile property of non-senescent cells. Subsequently, we hypothesize that the large senescent cell bodies are like a wall that does not form a strong physical bonding but attracts non-senescent cells only very weakly. A cellular Potts model (see "Methods" for details) was used to recapitulate these experimental results and validate our hypothesis. Importantly, the ensemble average of ω (and R) was used in the CPM simulations to estimate the interfacial energy associated with the interface between senescent cells and normal cells. successively recreate the behavior of an active normal cell moving around a senescent cell as observed in experiments (see Fig. 3a and Supplementary Video 2). Recapitulation of the experimental results was carried out by assuming the same set of parameter values, which was validated in a previous study investigating MDA-MB-231 cell motility 19 (see "Methods" for further detail). The given challenge in this work was to estimate a proper value of the interaction energy E sn for the interface between two different types of cells: a large, non-motile (S = 0) senescent cell and a small, motile (S = 2.8), normal cell. The overall sizes of the two cells forming a pair were controlled by the parameters V normal target and V sen target , where V normal target was fixed to be 3600 μm 2 and V sen target was scanned as one of two control parameters. Figure 3b plots the distance Rs (top frame) between the centroids of two cells forming a pair for n = 60 independent runs made with different initial conditions. Also shown in Fig. 3b (bottom frame) are the matching (unwrapped) rotation angles of the normal cell about the centroid of a senescent cell. Even for a fixed set of CPM parameter values, different runs would produce a small range of different rates of θ increase since there were stochastic events of reverse turns which temporally slowed down the overall rotation (see Fig. 3a and Supplementary Video 2). For the variable θ plotted in Fig. 3a, all such reverse turns were "unwrapped" so that the overall rate of θ (i.e., ω ) increase could be visually more evident.
A phase diagram of ω is drawn with two control parameters: V sen target , the targeted volume of the host senescent cell, and the interfacial energy E sn , controlling the level of adhesiveness at the interface of two cells forming a pair (see Fig. 3c). The solid (dashed) white line(s) in Fig. 3c represent the level curve(s) of the experimentally obtained ω (ω ± �ω) . Likewise, we have constructed a phase diagram of R using the same set of control parameters as shown in Fig. 3d. For an obvious reason, the average distance ( R ) between the two (senescent vs. non-senescent) cells forming an attached pair has a strong (weak) dependence on V sen target ( E sn ): In fact, R is mainly determined (c) Phase diagram of (time as well as ensemble) average angular velocity given in a 2D plane of E sn and the target volume of the hosting senescent cell V sen target . The white solid line is the level-curve of the experimentally measured value of ω = 3.9 × 10 -3 rad/min, and the dashed white lines mark the ranges of �ω = 2.6 × 10 -3 rad/ min. Shown in the right-hand side of the phase diagram is a scan of ω along the vertical direction (red Another important feature of the cell pair shown in Fig. 3a is the shape of the normal cell crawling along the boundary of a senescent cell. Shown in Fig. 3e is a phase diagram of p-value (= perimeter √ V normal ), a dimensionless measure of cell shape which has been considered as an important characteristic of cell crawling behavior. On top of it, we have placed a level curve of (ensemble and time averaged) experimental p-value (4.47 ± 0.65, solid blue line), which is quite close to that (4.28) of the CPM parameter set of the "red dot" ( V sen target = 3.4 × 10 4 , E sn = −20 ). But we should point out that the perimeter length, measured with experimental image data, can vary significantly depending on the way images are processed; obviously, a low-pass filtering will reduce the perimeter length (while not affecting the area much), subsequently leading to a smaller p-value. For the specified p-value (of the level curve), we have used only unprocessed raw images without any low-pass filtering. After all, ω and R are considered as a more reliable benchmark than p-value.
Shown in Fig. 3f are two-level curves (one for R and the other for ω ) intersecting at the position marked by a red dot. So, we have identified a particular value of interfacial energy E sn (= −20) which is experimentally most relevant for the inhomogeneous pair of MDA-MB-231 cells. Also illustrated alongside the level curves are six different configurations of the cell pairs, each of which corresponds to a square-marked position on the phase diagram. According to Fig. 3d, it is clear that R is almost linearly proportional to V sen target but their proportionality changes as a function of E sn ; also, we see that the curvature (or V sen target ) of the boundary of the senescent cell barely affects the shape of the adherent normal cell, whereas E sn significantly does so. The graph on the righthand side of Fig. 3c shows that the value of E sn also affects ω greatly (when the value of V sen target is fixed). Finally, Fig. 3g shows a phase diagram of (time-averaged) instantaneous tangential speed v θ , and one thing stands out clearly: v θ barely changes as a function of V sen target , which is consistent with the experimental result. In fact, it is E sn that affects not only the shape of the moving non-senescent cell but also its tangential speed: The more the cell is elongated with a stronger adhesion the faster it migrates.
Actively moving normal cells in the vicinity of a senescent cell in confluent population. The first application of CPM in biology was to explain the phenomenon of "phase separation" between two different kinds of cell populations that had different intercellular surface energies for different types of interfaces. Here, we have a very similar issue: How might the motility of non-senescent cells be modified when there are large senescent cell bodies scattered in their confluent population. Having figured out the experimentally relevant set of parameter values ( E sn = −20 , V target = 3.4 × 10 4 μm 2 ), we are now ready to run a realistic simulation as illustrated in Fig. 4a (also see Supplementary Video 3). A large senescent cell body, which can move only passively, is placed in the middle and barely moves. Four exemplary moving traces of the color-marked normal cells in Fig. 4a are given in Fig. 4b. All four marked cells in confluency have come across with the senescent cell in the middle and stayed in contact with it for some duration until they departed for elsewhere. The residence time upon each random visit to the senescent cell varied significantly from one to the other cases and its probability distribution fits well to a gamma distribution (α = 3.17, β = 62.00) as illustrated in Fig. 4c. The mean residence time was measured to be 196 min. For the same set of parameter values, but for the case of a single normal cell encircling a senescent cell, its residence time is almost infinite (not shown). And the significant reduction of the residence time in a confluent population is, of course, a consequence of normal-to-normal cell interaction which shows a much stronger adhesion with E nn ∼ −65 19 . As in the CPM simulation results of Fig. 4c, we have measured the residence time τ exp contact of normal cells (n = 30) that are in contact with a senescent cell in experiments, and its distribution is given in Fig. 4d. Like its CPM counterpart of Fig. 4c, the pdf fits a gamma function reasonably well (α = 4.22, β = 69.75) with a mean of 294 min, which is about 50% larger than that of its CPM counterpart; and we speculate that senescence-associated secretory phenotypes (SASPs) could be attributed to this increment as discussed in 14 .
We have also made CPM simulations mimicking the confluent binary mixture of cells in a larger domain to see how the motility of actively crawling normal cells would be altered in the presence of many large senescent cells (see Fig. 5a, left frame). As in the cases of heterogeneous cell pairs of Fig. 2a, the normal cells have a clear tendency to encircle the senescent cell on each encounter (see Fig. 5a, right frame): The blown-up insets shown in Fig. 5b are two exemplary cases. Clearly, the directional persistence of the cells moving along the perimeters of senescent cells are enhanced significantly as shown in Fig. 5c [persistence times: τ no_contact = 43.42 (min) τ contact = 123.50 (min)]. In the meanwhile, instantaneous speed has no difference for the two different cases, and this seems to be a reminiscence of the fact that v θ does not depend on R for the cases of heterogeneous pair of cells. Consequently, each encounter of a normal cell with a senescent cell body adds a small (curved) trajectory; yet it does not cause a significant 'packing density' fluctuation in time 25 (along the moving trajectory) as τ contact is only about three times that of τ no_contact . As expected, both τ contact and l contact (contact length along the boundary of senescent cell) increase monotonically as a function of E sn (and V sen target ). That is, the less sticky the interface is (or the larger the senescent cell body gets), the longer τ contact and l contact become (see Supplementary Fig. S1).

Summary and discussion
To understand the effects of cellular senescence on the motile behavior of MDA-MB-231 cells in confluency, we followed a bottom-up approach where heterogeneous systems with a much simpler configuration, namely, a cell pair formed by a normal cell and a senescent cell, were compared to the previously studied homogenous counterpart, that is, a cell pair formed by two normal cells 19 . After carefully analyzing the trajectories of two different types of cells during many free encounters in the experiment, a CPM was set up to recapitulate the results of experimental observations. The same CPM simulation was repeated for many different sets of parameter values, and a suitable value for the interfacial surface energy relevant for the interface between the two different www.nature.com/scientificreports/ types of cells was estimated. A similar strategy was applied to confluent populations having inhomogeneities (i.e., senescent cells), with special attention given to the residence time (i.e., contact time duration) of the normal cell, during which it stayed in contact with the hosting senescent cell. Unlike the pairwise rotation of normal cell pairs, the movement of the normal cell around a hosting senescent cell was non-mutual. Moreover, the characteristics (size, shape) of the involved senescent cells varied significantly from one to another due to the inherent heterogeneity of the cultured sample. Importantly, our analysis showed a strong 1/R functional dependence for the normal cell's angular speed ω , while its tangential speed along the senescent cell boundary v θ remained constant (i.e., invariant over a significant curvature variation).
Based on the average values of ω and R obtained in experiments, the adhesion energy E sn between senescentnon-senescent cell pair interaction was estimated to be much higher (E sn = − 20) than that (E nn = − 65) of equivalent non-senescent homogeneous pair interaction reported in our recent study 26 . Then, the significantly large E sn (compared to E nn ) should be responsible for the nearly ballistic movements of non-senescent cells that are in contact with senescent cells. After all, for the normal MDA-MB-231 cells the boundaries of hugely expanded senescent cells merely acted like an impenetrable wall which was far less sticky than peer normal cells in the population. Although there are previous studies collectively suggesting a role for aberrant adhesion protein (e.g., integrin) signaling in senescence, but how cellular senescence modifies cell-matrix adhesion is poorly understood 27 . In the future, it will be important to identify the biochemical mechanism underlying the reduced interfacial surface energy caused by cellular senescence.
Within a confluent population on a 2D substrate, non-senescent MDA-MB-231 cells exhibit a "superdiffusive" motility, where the directional superdiffusivity is even stronger than that of freely crawling cells 26 . Yet, the www.nature.com/scientificreports/ neighbor-enhanced persistence length was measured to be only 72 μm, which was about 1.2 times a cell diameter (cell diameter is assumed to be 60 μm), and the moving trajectory of a normal MDA-MB-231 cell could be viewed as a sequence of uncaging (or position swapping) events or position swapping of two adjacent cells. For the actively crawling normal cells, moving along the less sticky boundary of senescent cells would be energetically more favorable than joining into a densely-packed, stickier population: Therefore, the arclength along the senescent cell boundary covered over residence times (see Fig. 2a) was naturally far greater than the diameter of a normal cell. Then, the detachment of a normal cell from the boundary of the senescent cell could be incurred by a stochastic position swapping event with neighboring normal cell(s) as discussed in ref 19 .

Conclusion
With the help of CPM simulations we found that senescent MDA-MB-231 cells, which emerged spontaneously and randomly in space and time in cell culture, act like large physical obstacles providing only a low level of stickiness for actively crawling normal MDA-MB-231 cells in confluent populations. Consequently, each encounter of a normal cell with a senescent cell body adds a small but significant piece of circular trajectory on top of already super-diffusive moving trajectories of normal cells in confluency. Although there are several different versions of CPMs incorporating different mechanisms, for example, for the directionally persistent motility, the CPM model that we have employed in this investigation was validated earlier in a previous study modelling the super-diffusive nature of MDA-MB-231 cells in confluency 19 . This CPM includes neither detailed biochemistry orchestrating the cytosolic actin dynamics nor the biochemical interactions mediating the interfacial adhesiveness; therefore, it is quite surprising and fascinating that this simple model could recapitulate the observed phenomena in a very quantitative manner. Lastly, we should point out that senescent cells, in general, are known to excrete various SASPs 16,17,28 which can act as a chemo-attractant (or -repellant); a property that was attributed to the formation of cell clusters at the core of each senescent cell body in our previous study. In this work, we did not incorporate the potential chemotactic force that could have facilitated the adhesive property of the senescent cell boundary; and we speculate that it could have been an additional factor that made the average τ exp contact about 50% larger than the average τ contact of CPM simulation. Finally, we should comment that although the 2D phenomenon that we www.nature.com/scientificreports/ discuss in this paper is conceptually important, most biological tissues and organs are three-dimensional. As such, we have a plan to investigate the same issue in more physiological 3D environments.

Methods
Cell culture. High-density culture of MDA-MB-231 cells (~ 70% confluency) growing in DMEM culture medium (10% FBS) on a 3.5 mm petri dish was trypsinized with Trypsin-EDTA solution 10 × (Cat. No. 59418C, Sigma-Aldrich) diluted to 5 × with DPBS (Cat. No. D8537, Sigma-Aldrich) for 3 min after removing the DMEM medium. The solution was then gently stirred, centrifugated, and the supernatant discarded. The remaining cells were washed with DPBS solution and centrifugated again; this washing procedure was carried out twice before the addition of 1 ml of DMEM solution to the final remaining cells. Following cell counting, two different samples that either have a medium initial confluence (60%) or low initial confluence (20%) were prepared. All samples were stored in a laboratory incubator at 37 °C and 5% CO 2 perfusion (for about one day for low-density samples and 4-5 days for high-density samples) until they were transferred to the microscope incubator chamber for observation. MDA-MB-231 cells were purchased from Korean Cell Line Bank (KCLB No 30026).
Live-cell imaging. A small, thin cylindrical, temperature-controlled chamber incubator was designed and lab-built for live-cell imaging over long periods. The inside of the chamber was kept at a constant temperature of 36 ± 0.1 °C with a heating coil controlled by a temperature controller (SDM9000, Sanup, Korea) equipped with a PT100 temperature probe. Two ITO-coated optical windows of the incubator along the imaging axis were also electrically heated to avoid condensation. A mixed gas of CO 2 (5%) and air (95%) was continuously perfused to the incubator. The incubator was positioned and fixed on the x-y stage of an inverted microscope (Olympus IX710). A sample petri-dish containing the designated samples of MDA-MB-231 cells in a proper medium (89% RPMI serum, 10% Fetal Bovine serum, and 1% Penicillin-Streptomycin) was then loaded into the chamber. Time-lapse phase-contrast images were acquired with time intervals of 1 min using ProgRes MFcool CCD Camera (Jenoptik, Germany) with 4 × (NA 0.13) phase-contrast objective lenses and Micromanager software. The camera has a sensor resolution of up to 1360 × 1020 pixels with a measurement precision of 2 µm. A mechanical shutter was used to block the light in between image acquisitions. Imaging processing. The crawling MDA-MB-231 cells in our experiments show a significant activity of pseudopods, whose existence can make them look more elongated and ramified than those of CPM simulations. Since these pseudopod activities seem irrelevant to the normal cell-senescent cell interface formation, (which is relevant to E sn for the case of CPM), they are excluded in the cell-shape analysis. On the other hand, with the gray-scale phase-contrast images and the separation process relying on a (gray scale) threshold alone, it is not possible to clearly distinguish pseudopods from cell bodies as the boundaries between them are sometimes unclear. On top of the innate cellular inhomogeneity of MDA-MB-231 cells, this could be a significant factor contributing to the uncertainty in the measured p-value. Perhaps more importantly, the perimeter length can vary significantly depend on the way image was taken and processed. In this sense, ω and R are considered as a better benchmark than p-value, when they were compared with those of CPM simulations. The cellular potts model simulation. Cellular Potts Model (CPM) is a lattice-based model that represents a given biological cell as a group of connected lattice sites. Each cell is allowed to interact with the cell surroundings along its boundary (e.g., other cells, physical substrates and barriers, and medium). The lattice domain evolutions (e.g., cell-cell interaction) are governed by the system's free energy. CPM evolves in time and space based on a Monte Carlo method, in which (1) one of the lattices is picked randomly; (2) subsequently, we select a 'neighbor site' randomly among all 'surrounding sites' of the chosen lattice site; (3) the probability P accept of the initially chosen site accepting the selected neighbor site's cell-type is based on Metropolis probability [ P accept = 1 (if �H < 0); e −�H/T (if �H > 0)], where H and T represent the change in the total Hamiltonian and the level of cell membrane fluctuations 29 , respectively. The system would evolve towards the state of minimum energy. The H includes the following three components: www.nature.com/scientificreports/ where A i and P i are the current area and perimeter of the cell indexed by i , and E τ i,τ j is the interfacial (adhesion) energy. area , perimeter , A target , P target , and S are fixed parameters. H adhesion models the change in the energy conferring the persistence of moving cells. For the chosen set of parameter values that we used for Figs. 4 and 5, the size of �H P constraints = perimeter � i (P target − P i ) 2 [ �H A constraints = area � i (A target − A i ) 2 ] is most (least) significant and H persistence comes second as for updating the physical state of each cell (see Supplementary  Fig. S2 for further details). p i is the polarity vector of i'th cell, and x i the hypothetical displacement of the cell's centroid if the neighbor's lattice site is copied to the target lattice site. At each Monte Carlo step, p i evolves with the addition of r i , which is the spatial displacement of the centroid at the current time, and the factor reflecting directional memory loss −p/τ , where τ represents the time constant for the memory loss. If τ → ∞ , p i (t) = p i (t − 1) + �r i , which means the current polarity vector depends on the entire history of the passage taken by the cell. If τ → 1 (i.e., the minimum time step), p i+1 (t) = �r i .
All (non-senescent) model cells have directional persistence in their motility and some degree of "stickiness" with neighboring cells; so, in a confluent population, they can be viewed as super-diffusive particles interacting with one another via interfacial interaction 26 . On the other hand, senescent cells are deprived of self-persistent motion (i.e., the parameter associated with self-propulsion is zero), since the cell bodies of senescent MDA-MB-231 cells are in general very huge and they barely move in comparison with non-senescent tumor cells 14 . In this investigation, we have varied only two parameters in efforts to investigate the effect of the senescent cells' overall size (or the curvature of cell boundary in 2D space) and adhesiveness between the senescent and non-senescent cells: namely, senescent cell's target volume V sen target [and also P sen target in proportion to (V sen target ) 1/2 ] and adhesion energy E sen−cell between senescent and non-senescent cell boundaries. All other parameter values are fixed as following: τ = 4 , area = 1 , surface = 1 , A target = 100 , P cell = P senescent = 0.9 , E cell−cell = −65, E medium−medium = 0 , E cell−medium = 0 , and the parameter T was set to 10. Here, the subscript (or superscript) ' sen ' and ' cell ' refer to the senescent and normal cells, respectively. In the simulation, 1 MCS corresponds to updating sites the number of times equal to the number of total sites in the system and it was estimated to be 2 min in earlier work of ours 19 . All simulations were done inside a 600 × 600 square lattice space with periodic boundary conditions for over 10 4 Monte Carlo steps (MCS) using a CPM simulation platform Morpheus 30 .