Proteolytic and non-proteolytic regulation of collective cell invasion: tuning by ECM density and organization

Cancer cells manoeuvre through extracellular matrices (ECMs) using different invasion modes, including single cell and collective cell invasion. These modes rely on MMP-driven ECM proteolysis to make space for cells to move. How cancer-associated alterations in ECM influence the mode of invasion remains unclear. Further, the sensitivity of the two invasion modes to MMP dynamics remains unexplored. In this paper, we address these open questions using a multiscale hybrid computational model combining ECM density-dependent MMP secretion, MMP diffusion, ECM degradation by MMP and active cell motility. Our results demonstrate that in randomly aligned matrices, collective cell invasion is more efficient than single cell invasion. Although increase in MMP secretion rate enhances invasiveness independent of cell–cell adhesion, sustenance of collective invasion in dense matrices requires high MMP secretion rates. However, matrix alignment can sustain both single cell and collective cell invasion even without ECM proteolysis. Similar to our in-silico observations, increase in ECM density and MMP inhibition reduced migration of MCF-7 cells embedded in sandwich gels. Together, our results indicate that apart from cell intrinsic factors (i.e., high cell–cell adhesion and MMP secretion rates), ECM density and organization represent two important extrinsic parameters that govern collective cell invasion and invasion plasticity.


Pilot simulations for calibrating values of perimeter constraint (λ p ) and intrinsic motility (µ 0 )
The value of perimeter constraint (λ p ) was assumed in the model, and a value of 0.5 was chosen after several simulations so as to minimize excessive cell deformations (data not shown). This value of λ p was then used for calibrating the value of cell motility strength µ 0 . The strength of intrinsic motility µ 0 was optimized such that (i) the cells showed a rotational motion, similar to, 1 when confined within a non-degradable circular geometry, and, (ii) the cells showed a translational motion when placed in an extracellular matrix with already available free path. Similar to the observations of, 1 cells, when placed inside a non-degradable confinement, exhibit rotational motion (Supp. Video V1). After running several simulations, a value of µ 0 = 50 was determined to satisfy the above constraints. To determine the approximate value of 1 MCS in our simulations, we run the simulations in hindrance free environment (i.e., without ECM fibres) and calculated the cell speed in units of µm.MCS −1 . This value of speed ( 425 1800 ≈ 0.23 µm.MCS −1 ) was then equated to the experimentally reported values of cell speed (≈ 23 µm.h −1 ) 2 to estimate the value of 1 MCS. From this, 1 MCS was estimated to be approximately 36 seconds.

Model implementation 2.1 CC3D based cellular Potts model
The open source framework CompuCell3D (CC3D) was used to develop our multiscale computational model. 3 The package is based on highly modular design and provides the facility to incorporate the cell based cellular Potts model (or Graner-Glazier-Hogeweg(GGH)) Model) 4, 5 and PDE based chemical reaction diffusion dynamics in a single framework. CC3D plugins, written in C++, were used to incorporate different energies including cell-motility, adhesion energies, area constraint, and perimeter constraint. Additionally, while all the plugins used in our model are part of the standard CC3D package, we wrote two steppables (python routines to add MCS level events) to incorporate aspects of ECM degradation and cell-polarity evolution and 2 steppables exclusively for visualization and data logging. Brief of these steppables is given below.
PolarityEvolutionSteppable: This steppable was developed to update the cell polarity after every MCS. As discussed in the material & method section and used in previous studies , 1, 6 the average displacement of cell during last 10 MCS (D avg,10 ) was quantified and then the new polarity of the cells was set equal to the direction of D avg,10 . Initial polarities of all cells were set to be random. To integrate this polarity with cell-motility, the direction was multiplied by strength of motility (i.e., µ 0 ) and this vector was then supplied to the 'ExternalPotential' plugin (from python code only) which then implements the cell-motility ( 3 and CC3D manual). This steppable also logs the centroid of all cells to CSV file for further analysis.
ECMDegradationSteppable: This steppable was developed to implement MMP mediated ECM degradation and was called after every MCS. This plugin checks the MMP concentration at all pixel positions and removes the fibre pixel (if any) from a pixel if the MMP concentration at that pixel location is ≥ 1.
IdFieldVisualizationSteppable: This steppable defines a cell level scaler field and assigns a random number to each cell. This number can then be used to visualize an individual cell.
MMPSecretionSteppable: This steppable logs the temporal profile of total MMP in the lattice.

Custom tool for generating fibrillar matrix
To generate ECM lattice with varying fibre density (ξ ) and orientation (ψ), a customized C++ program was written. Line drawing algorithm was used to generate ξ fibres, each of them 30 − 40µm long and aligned at an angle ψ. The utility outputs the generated lattice in *.PIFF formate which then serves as input to the CC3D model using the 'PIFInitializer' steppable (part of standard CC3D package).

MMP activity of MCF-7 increases with increase in collagen density
MCF-7 cells (50 × 10 3 cells per cm 2 ) were grown in 500 µl of media on 3D collagen gels of varying densities (0.5, 1 and 2 mg/ml, respectively). After 36 hours of incubation, media were collected and centrifuged at 5000 rpm for 10 minutes to remove cells and gel particles (if any) and subjected to gelatin zymography. 7 In brief, equal volumes of media were mixed with 6X Laemmli non-reducing buffer, incubated for 15 minutes at room temperature. Samples were then electrophoresed on 8% SDS-PAGE gel containing 1 mg/ml gelatin. After electrophoresis, the gels were washed twice with renaturation buffer containing 2.5% triton X-100 for 30 minutes each, followed by rinsing 3 times for 30 minutes with reaction buffer containing 50 mM Tris-HCl buffer (pH 7.6) with 5 mM CaCl 2 , 200 mM NaCl and incubated overnight at 37 o C. The gels were stained with 0.5% Coomassie brilliant blue G-250 solution containing 10% acetic acid and 20% methanol for 30 min and destained with 7.5% acetic acid solution containing 10% methanol. Areas of gelatinase activity were detected as clear bands against the blue-stained gelatin background, which was quantified by densitometric analysis using ImageJ software. On analysis, we observed an increase in MMP (MMP-1, MMP-2 and MMP-9) activity with increase in collagen density ( Figure S1).  Figure S8. Effect of ECM hindrance on cell translocation. Simulations were performed with J ce = J cf , i.e., cells do not have any preference to attach to ECM fibres. Proteolysis was turned off in these simulations (i.e., λ = 0). Translocation of invading cells for varying cell-cell adhesion (J cc ) and ECM density (ξ ).

Cell invasion in the presence of cell division
To study any potential influence of cell division on cell invasion, cells were allowed to divide. In all these simulations, cell division time was considered to be 24 hours (= 2400 MCS). To incorporate cell division into our model, 'MitosisSteppable' available with CompuCell3D was used. 3 Algorithmically, every cell was assigned an internal clock (age) to keep track of time (in hours) passed after the last division ( Figure S9A). At t = 0, age of all cells were initialized to a random integer value age ∈ [0, 24]. After each 100 MCS (∼1 hour), the value of age for each cell was incremented by one. Once the age of a cell reaches 24 hours, the cell was divided into two daughter cells and the age of daughter cells were reset to '0'. Random initialization of age was used to make sure that (i) some cells divided within 1800 MCS (otherwise they need atleast 2400 iterations), and (2) cell division happens asynchronously (i.e., all cells need not to divide at the same point of time). Due to this asynchronous cell proliferation, population size increased from 69 cells at t = 0 to ∼120 cells at the end of simulations (t = 18 hours).
Simulations were performed for varying cell-cell adhesions (i.e., J cc = 1, 16 and 40), moderate ECM density (i.e., ξ = 600), moderate ECM proteolysis (i.e., λ = 0.05 s −1 ) and randomly aligned matrices (i.e., ψ = 'Random') (similar to Figure 2). For each value of J cc , 10 simulations were performed. Radius of gyration (RoG), percentage ECM degradation and total MMP secretion were quantified after 1800 MCS. Since each cell division event introduces two new cells and removes the parent cell during the course of simulation, not all cells were present during the entire duration of simulation (1800 MCS). Hence, neither cell translocation (D) after 18 hours (total simulation duration) nor the total distance traveled by a cell (d) in 18 hours could be quantified. We, therefore, relied on radius of gyration (RoG) for quantifying cell scattering. Upon quantification of RoG, no drastic differences were observed when cell division was taken into account for J cc = 1, 16; slight increase in RoG was observed for J cc = 40 ( Figure S9B). Mild increase in percentage ECM degradation and MMP secretion was observed in the presence of cell proliferation ( Figure S9C and S9D). Since cell proliferation increases total cell population size, this mild increase in total MMP secretion and percentage ECM degradation was expected. Thus, these results suggest that cell proliferation does not appreciably change the attributes of cell invasion observed without taking cell division into account.

Cell invasion in the presence of cell supply
To study the effect of continuous supply of cells on cellular invasion, simulations were performed where cells were continuously added to the centre of the grid (where initial population was present). In this setup, the initial population represents a 2D slice of a 3D tumour where cells are continuously supplied from upper and lower layers of the 3D tumour. The rate of cell addition was kept proportional to the free space in the tumour centre (where initial population was present) ( Figure S10A). Precisely, after every 10 MCS, N FreePixels a 0 cells were added to the lattice centre, where N FreePixels represents the number of free pixels and a 0 represents the number of target pixels (or target area in pixels) per cell (100 in our case). This continuous supply of cells caused a drastic increase in total number of cells in the population (> 200). Cell scattering (RoG) and percentage ECM degradation were quantified after 1000 MCS (instead of 1800 MCS), as beyond this time, due to their large number, cells started reaching the boundaries and moving along the boundaries. Temporal MMP profile was also traced for 1000 MCS. As expected, compared to the case of constant number of cells, higher levels of MMP secretion and percentage ECM degradation were observed for the case of continuous cell supply ( Figure S10B and S10C). However, the behaviour of their their dependence on J cc remained the same. While the dependency of RoG on J cc differed between the two cases, it is difficult to make a direct comparison between the two cases. This is because for the case of continuous cell supply, cells are added at the centre thereby creating a bias in their spatial distribution, and hence influence the RoG measure differently ( Figure S10D).

Supplementary Videos
Supplementary Video V1. Rotational motion of cells confined in non-degradable circular geometry. Video shows that when cells trapped in a circular confinement, they exhibit rotational motion. Simulation time is ≈ 18 hours. The circular geometry was designed using CellDraw utility available in CompuCell3D(CC3D) package.
Supplementary Video V2 -V4. Cell invasion for varying cell-cell adhesion. Video show the cell invasion with high cell-cell adhesion (i.e., J cc = 1, Supp. Video V2), medium cell-cell adhesion (i.e., J cc = 16, Supp. Video V3) and low cell-cell Supplementary Video V5 -V8. Proteolysis independent cell invasion in random matrices. Videos show cell invasion when cells with different cell-cell adhesion abilities were allowed to invade through matrices composed of randomly aligned fibres at two different ECM densities. ECM degradation was turned off. Video V5: High cell-cell adhesion and low ECM density (i.e., J cc = 1 and ξ = 400). Video V6: Low cell-cell adhesion and low ECM density (i.e., J cc = 40 and ξ = 400). Video V7: High cell-cell adhesion and high ECM density (i.e., J cc = 1 and ξ = 1000). Video V8: Low cell-cell adhesion and high ECM density (i.e., J cc = 40 and ξ = 1000). Switching off of ECM proteolysis caused a significant drop in translocation for all combinations of cell-cell adhesion and ECM density.
Supplementary Videos V9 -V12. Proteolysis independent cell invasion in aligned matrices. Videos show cell invasion when cells with different cell-cell adhesion abilities were allowed to invade through aligned matrices at two different ECM densities. ECM degradation was turned off. Video V9: High cell-cell adhesion and low ECM density (i.e., J cc = 1 and ξ = 400). Video V10: Low cell-cell adhesion and low ECM density (i.e., J cc = 40 and ξ = 400). Video V11: High cell-cell adhesion and high ECM density (i.e., J cc = 1 and ξ = 1000). Video V12: Low cell-cell adhesion and high ECM density (i.e., J cc = 40 and ξ = 1000). A significant increase in cell translocation was observed in aligned matrices compared to random matrices (Supp. Videos V5 -V8). The exact increase in translocation was observed to be cell-cell adhesion dependent with maximum at highest cell-cell adhesion (i.e. J cc = 1) ( Figure 5C, Figure S7).
Supplementary Videos V13 and V14. Migration of MCF-7 cells in sandwich collagen gels. Videos show migration 11/12