In brief

Leineweber and Fraley develop an approach to measure multiple types of interactions between individual cells and their surrounding 3D matrix while also tracking migration. By modeling these multiscale data, they accurately predict migration heterogeneity and identify a link between cell-matrix adhesion and the uncoupling of cell speed and persistence.

INTRODUCTION

Cell migration is a complex behavior that emerges from biophysical and biochemical interactions between thousands of molecular parts within and between cells and their environment. Comprehensively measuring migration machinery across space and time is not currently possible; therefore, pairing experiments with modeling efforts is crucial to advancing our understanding of cell migration. Most studies have been limited to 2D cell migration on flat surfaces, where cells flatten and become easy to image due to the absence of a 3D extracellular matrix (ECM). In these studies, the complexities of the migration machinery have been abstracted into three predominant subprocesses that run concurrently and are spatially coordinated: protrusion, adhesion, and contraction. 1 Many features of cell morphodynamics and migration on 2D substrates can be explained by treating these subprocesses as functional modules connected in a simplistic circuit [2][3][4] originally predicted by Abercrombie, which goes as follows: (1) protrusion of the leading edge, (2) formation of strong adhesions at the leading edge, (3) aging adhesions at the trailing edge that are (4) released by contraction of the cell leading to forward movement. 5 More recent modeling efforts have focused on gaining a molecular-level understanding of each subprocess. 6 Still, this conceptual framework does not account for the major role that the ECM plays in confining and resisting cell movement.

Navigation through 3D tissue-like environments is physiologically relevant for many migratory cells; however, the technical challenges inherent to studying 3D migration are formidable.

Particular challenges include microscopy and image analysis limitations in 3D and the increased complexity brought about by additional modes of migration in 3D vs. 2D. [7][8][9][10] Predictions of cell behavior in 3D based on data acquired in 2D are often unreliable because signaling and mechanical parameters do not always directly translate from 2D to 3D. [11][12][13][14][15][16][17] Moreover, an additional subprocess is necessary to consider for migration through 3D environments: ECM remodeling. [18][19][20][21] In vivo interstitial pore sizes range from $0.025 to 0.1 mm in diameter, 22 and those of the basement membrane range from $0.6 to 3.85 mm. 23 These pore sizes are smaller than most cell bodies (10-100 mm in diameter) and cell nuclei (3-7 mm in diameter when deformed during matrix metalloproteinases (MMP)-independent migration). 21 Therefore, cells migrating through tissues must remodel the matrix by either physical or biochemical mechanisms. Many 3D cell migration studies have used low density matrices that have much larger pore sizes than the ECM in vivo, 21,[24][25][26][27][28] which enable cells to migrate without remodeling. Indeed, it has been demonstrated that MMP-independent migration is a function of pore size. 21 In confining ECM, in which pore sizes are representative of in vivo tissues, it is not well understood how the processes of matrix remodeling, protrusion, adhesion, and contraction are integrated to produce different modes of 3D migration. Additionally, the relationships between these processes are unclear in the context of 3D cell migration.

Without an integrated framework for whole-cell 3D migration in confining ECM built on the four key subcellular processes, it remains difficult to distinguish the origins of migration heterogeneity. For example, cells of the same type can display significantly

ll

OPEN ACCESS different migration behaviors, even in ostensibly the same 3D environment. 18,[28][29][30][31] Does cell migration heterogeneity arise from fundamentally different interrelations between the subprocesses, or can it be explained by a threshold or shift in activity of one or several of the subprocesses? Experts argue that answering this question may be the most important contribution toward unraveling the mechanisms of migration. 32 An integrated framework would also help contextualize seemingly contradictory results from experimental perturbations. 1 The pathways controlling migration are often nonlinear and redundant. They may also feedback on each other, making it exceedingly difficult to address questions through single-plex molecular experiments or mechanochemical models. For example, several different combinations of physical mechanisms can be fitted to explain the same experimental behaviors, but directly testing them can be impossible. Obtaining sufficient physical and biochemical constants required for accurate models in 3D systems can be difficult to impossible. An alternative approach for modeling how cells process subcellular information into wholecell behaviors is to use a data-driven methodology based on the quantitative effects of perturbations to key high-level process modules. [33][34][35][36][37][38] Here, we present a data-driven model of 3D cell migration in confining ECM built on integrated measurements of protrusion, adhesion and contraction (traction), and matrix remodeling in single migrating cells. Single-cell tracking data revealed three modes of coupling between cell speed and persistence, which were linked to distinct combinations of cell-ECM interactions. This work represents an advancement in our understanding of how heterogeneous migration can arise and provides actionable insights into engineering cell migration behavior.

RESULTS

Cell migration is heterogeneous and less common behaviors can be enriched by perturbing biophysical processes that dictate cell-ECM interactions

To determine whether different levels of migration subprocess activity or different interrelations between the subprocesses account for heterogeneous migration behaviors within a given cell population (Figure 1A), it is necessary to study ''average'' cells, as well as ''rare'' cells that display fewer common behaviors. The distributions of migration behavior for cells migrating in confining 3D collagen type I matrices can be characterized by the persistent random walk (PRW) model, 13,39 which uses the mean squared displacement (MSD) of cells to attribute values of cell speed (S) and persistence time (P) to migrating cells (Equation 1). The other parameters in this equation include the time lag (t), dimensionality of the tracking (n), and the positioning error (SE).

For example, MDA-MB-231 (MDA) display a wide range of total displacements (10.32-183.72 mm), speeds (0.022-0.395 mm/ min), and persistence times (0.431-994.846 min). Across multiple cell types, including MDAs, HT1080s, and HFF-1s, both total cell displacement and cell speed are logarithmically distributed (Figures 1B and 1C), and persistence shows a bimodal logarith-mic distribution (Figure 1D). This suggests a general distribution of cell migration behavior in confining 3D collagen matrices and shows that capturing less common migration behaviors requires extensive sampling.

We hypothesized that inhibitors targeting the individual processes of contractility (rho-associated protein kinase inhibitor, [ROCKi]), matrix remodeling (MMPi), cytoskeletal protrusion (F-ACTINi), and adhesion (ITGB1i) could enrich different regimes of migration behaviors within the general distribution. As anticipated, these treatments shifted the peak of the distributions of displacement, speed, and persistence compared with vehicle control while remaining within the general ranges (Figures 1E-1H). For example, ROCKi, MMPi, and F-ACTINi increased the population of low-persistence cells and shifted the cell speed distribution toward lower values. ITGB1i inhibition, on the other hand, increased the number of highly persistent cells and shifted the cell speed distribution toward faster cells. To confirm that the effects of the small molecule inhibitors have similar effects to those of more targeted knockdowns, we also assayed cells with ROCK1 knocked down (Figures S1A-S1E) and ITGB1 knocked down (Figures S1F-S1I). The effects of the more specific knockdowns were more subtle but trended in similar directions as the inhibitor treatments. Cumulatively, these results suggest that changes in the level of activity of one or more processes could potentially account for naturally occurring migration heterogeneity. However, it is not clear whether shifts in activity are also accompanied by different interrelations between the core subprocesses.

Cell-ECM interaction measurements capture whole-cell biophysical behavior

To determine the interrelations between the core subprocesses and migration outcomes, each process must be simultaneously measured in individual cells as they migrate through the ECM.

Because the activity of a given process is not well described by any single measurement, we integrated nine imaging-based measurements of cell-ECM interactions to read out various aspects of subprocess activity. These measurements were captured by time-lapse z stack imaging in three fluorescent channels for MDAs in 3D collagen as follows: blue fluorescent matrix-embedded beads enabled measurement of percent (%) bead movers, maximum bead displacement, and instantaneous bead displacement as readouts of cellular contractility against the matrix; green fluorescence from dye-quenched collagen (DQ) measured the remodeling of the matrix by cells; red fluorescent protein expression in cells allowed us to measure instantaneous cell displacement, protrusion rate, lifetime, and max length as readouts of cytoskeletal activity; and the ratio of the average instantaneous bead displacement over the average instantaneous cell displacement gave a measure of the coupling between a cell and the surrounding ECM (Table 1). We will subsequently refer to this imaging platform and resulting measurements as ''biophysical imaging.'' 3D volume view time-series demonstrate the technique's ability to capture heterogeneous subprocess activity between individual cells, where the differences in cell protrusive activity (red), bead movement (blue), and matrix remodeling (green) are easily visualized (Figures 2A and 2B). The cell in Figure 2A remains rounded, with minimal protrusive, contractile, or matrix remodeling activity. In contrast, the cell in Figure 2B retracts a protrusion, resulting in displacement of the microbeads, and has a strong DQ signal indicating matrix remodeling. In some instances, we found cells remodeling the matrix by protruding into it and retracting (Video S1). In other cases, matrix degrada-tion was localized at the cell body near the neck of longer protrusions (Video S2). These results demonstrate the utility of this imaging approach to measure heterogeneous cell-ECM interactions at the single-cell level. Next, we measured how cell-ECM interactions changed in response to inhibition of each of the core subprocesses using the biophysical imaging platform. Because each inhibitor treatment shifted the distribution of cell migration behavior (Figures 1F-1H), we hypothesized that they would also shift cell-ECM interaction distributions such that distinct cell states could be identified. Inhibitors were used as opposed to molecular knockdowns to obtain more robust effects on the individual biophysical processes, while still resulting in migration behavior in the range of vehicle cells (Figures 1F-1H). 3D reconstructions representing an average cell from each inhibitor treatment are shown (Figure 2C) and distributions of cell-ECM interaction responses are plotted (Figures 2D-2L). Compared with vehicletreated cells, F-ACTINi cells (green) displaced the matrix less (Figures 2D and 2F), moved slower (Figure 2G), remodeled the matrix less (Figure 2H), and extended few protrusions (Figures 2I-2K). ROCKi (blue) decreased the extent of matrix displacement (Figures 2D and 2F), cell movement (Figure 2G), and matrix remodeling (Figure 2H), but increased cytoskeletal protrusion activity (Figures 2I-2K). MMPi cells (red) were slower on average than vehicle cells (Figure 2G), but surprisingly did not remodel the matrix to a significantly lesser extent (Figure 2H). MMPi treatment targeted the main family of ECM collagenases, so the insignificant decrease in the DQ measurement compared with control cells was unexpected. As an additional check, we confirmed that the DQ signal is an accurate readout of matrix degradation by using a degraded-collagen hybridizing peptide (CHP) as a secondary measurement (Figure S2). Finally, ITGB1i cells increased their instantaneous cell speed compared with the vehicle control. In total, each inhibitor treatment differentially regulated cell-ECM interactions and shifted cell migration distributions in distinct ways.

Principal component analysis (PCA) on this cell-ECM interaction dataset clustered cells in the same treatment conditions relatively well, despite the heterogeneity observed (Figure 2M; Table 2). The separation of the data along PC1 was driven by fairly equal contributions from instantaneous cell displacement, protrusion rate, instantaneous bead displacement, and the displacement ratio (Figure 2N). Interestingly, PC1 also seemed to order treatment groups from least to most migratory (compare Figure 2E with Figures 1F and 1G). PC2 separation was driven by protrusion length and lifetime in the positive direction, and instantaneous bead displacement and displacement ratio in the negative direction. This axis helped to separate the ROCKi group from the F-ACTINi and MMPi populations. To test the generalizability of this approach, we also performed biophysical imaging on HT1080 cells and mapped these measurements onto the principal components (Figure S3A). HT1080 cells clustered together and were most similar to the ITGB1i treatment group (Figures S3B-S3K), and we confirmed via western blot that HT1080s have less ITGB1 protein than MDAs. These results highlight this imaging platform's capacity to detect molecular differences between cells. Together, these results show that the nine cell-ECM measurements obtained from biophysical imaging can discriminate between cell states associated with distinct migration distributions arising from intra-and intercellular heterogeneity, as well as molecular inhibition of migration processes.

Cell trajectories are well modeled by the PRW model when speed and persistence are coupled Coupling between cell speed and persistence is not universal in 3D migrating cells With our biophysical imaging method established, we next sought to understand how cell-ECM interactions map to long-term cell migration trajectories. We devised a two-stage experimental approach in which cells were tracked for 8 h, followed by biophysical imaging for 1 h. Direct comparisons could then be made between long-term migration behavior and snapshots of the cell-ECM interaction state within the same cell. In this integrated experimental protocol, cell trajectories (Figure 3A) followed similar trends as those from our earlier cell tracking experiments (Figure 1E), where fluorescent beads and DQ were not included in the matrix. ITGB1-inhibited cells were again the most migratory, whereas F-actin inhibited cells were mostly non-migratory. MMP and ROCK inhibition decreased migration compared with the vehicle condition.

Because two migration parameters, speed and persistence, are sufficient to describe and predict cell migration trajectories, we fit the MSD of individual cell trajectories to the PRW model (Figure 3B) to extract values for cell speed (S) and persistence time (P) (Equation 1). We found that many of the analyzed cells, though not all, followed a migration regime wherein speed and persistence were coupled (Figure 3C). This coupling fits a general form of the ''universal coupling between speed and persistence'' (UCSP) equation that has been previously reported (Equation 2). 40 (Equation 2)

A subpopulation of cells, mostly from the ITGB1-inhibited population, did not follow this coupling law, and in fact appear to display a negative relationship between speed and persistence (Figure 3C, data outside of circle). Interestingly, we do not see uncoupling of S and P at low P. These results show that cell speed and persistence time are often, although not universally, coupled in migrating cells in confining 3D matrices. Because inhibition of ITGB1 produced the most cells in the decoupled migratory regime, this suggests that adhesion plays an important role in determining persistence time and in coupling persistence to speed. Cell-ECM measurements predict cell speed and persistence in coupling regime Because PCA, a linear data transformation method, clustered cell states (Figure 2E) and ordered them along PC1 similar to the order of their migration distributions (Figures 1F-1H), we reasoned that a simple linear regression could model the relationship between cell-ECM interactions and migration. We also considered that cells following the UCSP law may rely on a fundamentally different configuration of cell-ECM interactions than the cells that do not obey this relationship. Therefore, we used partial least squares regression (PLSR) to ask which cell-ECM interactions are the best predictors of S and P for cells that follow the UCSP behavior. All 511 possible combinations of the nine cell-ECM interaction measurements were tested as the independent variable matrix to fit the dependent variable S or P and leave-one-out cross validation was performed to account for over-fitting. The resulting goodness of fit (R 2 ) and predictive ability (Q 2 ) of every possible regression model was calculated, and the model with the highest Q 2 was chosen as the optimal model. Applying this approach with P as the dependent variable yielded an optimal model with an R 2 = 0.524 and a Q 2 = 0.449 (Figure 3D). The optimal model for P accounts for about 52.5% of the variance using only two PLS components (Figure 3E), and consists of bead-cell displacement ratio, protrusion length, instantaneous bead displacement, maximum bead displacement, DQ, and percent bead mover measurements (Figure 3F). The observed vs. fitted values of z-normalized P values shows that the model performs well across the different inhibitor treatments (Figure 3G).

PLSR analysis was also performed for S (Figures 3H-3K), which revealed that S is best predicted by a model consisting of the same cell-ECM measurements used for P, but without needing bead displacement or percent bead movers (R 2 = 0.529, Q 2 = 0.469). The loading for bead-cell displacement ratio correlates with low speed and persistence, whereas high maximum bead displacement, percent bead movers, and DQ are most correlated with enhanced speed and persistence (Figures 3F and 3J). The finding that these interaction measure-ments are the most predictive of speed and persistence update the UCSP model that previously reported only actin flow was required to predict migration behavior. 40 Using only the three protrusion measurements representing actin cytoskeletal activity in the regression model (Figure S4) achieved less predictive power for P (R 2 = 0.305, Q 2 = 0.201) and S (R 2 = 0.292, Q 2 = 0.180), and was suboptimal compared with the more comprehensive regressions (Figures 3G and 3K). In the updated model, the inclusion of bead-cell displacement ratio, DQ, and multiple bead measurements indicate that the processes of adhesion, matrix remodeling, and contractility are key determining factors along with actin-based protrusion dynamics for 3D cell migration.

To test the performance of these models, we simulated PRW migration trajectories using S and P values predicted from cell-ECM interactions (Figure 3L). We then compared these with simulated PRW trajectories using the S and P values extracted from cell MSDs (Figure 3M) and to the original real trajectories of the cells (Figure 3A). Trajectories predicted from cell-ECM interactions recapitulate the global effects of the inhibitor treatments and also the heterogeneity within each population. This generalizability of the model was next tested by incorporating HT1080s, which predominantly followed the S and P coupling regime (Figures S5A and S5B). The optimal model to predict HT1080 speed and persistence used similar biophysical measurements as that of the MDAs (Figures S5C and S5D). Importantly the model generated by MDAs and agnostic to HT1080s predicted HT1080 speed and persistence nearly identically as the model generated using both MDAs and HT1080s (Figures S5E and S5F). The cell-ECM measurements required for the optimal model demonstrate the necessity of capturing the coordination between cytoskeletal protrusions, contractility, matrix remodeling, and adhesion for this mode of cell migration.

Trajectories of cells whose speed and persistence are not well coupled are modeled by a distinct set of cell-ECM interactions The PLSR model that predicts cell migration in the UCSP regime (Figures 3D-3L) does not achieve a good fit for cells outside of this regime, suggesting that these globally uncoupled cells occupy a distinct cell-ECM interaction state. These cells had longer persistence times than those in the coupled range, but their cell speeds were within the range of the other measured cells (Figure 4A). Because these cells are highly persistent, we reasoned that they may be better modeled using the anisotropic PRW (APRW) model. APRW takes into account cells with preferred migration directions by deconvolving the migration into a primary direction and an orthogonal non-primary direction, each with a persistence and speed (P1,S1 and P2,S2, respectively). 39 When the cells with globally uncoupled speed and persistence were fit using the APRW model, a subset of these cells displayed coupling in the primary direction of migration (Figure 4B). Coupling was also observed for some cells along the non-primary direction axis (Figure 4C). Therefore, cells that do not display global coupling between speed and persistence, as modeled by PRW, can still display coupling behavior along primary or non-primary directions of migration, which is captured by the APRW model. Cell-ECM interaction state in anisotropic S vs. P coupling We next asked whether a distinct set of cell-ECM interactions were associated with globally uncoupled (APRW) cell migration. We further delineated cells that couple speed and persistence along the primary direction of APRW migration (red squares) from those that do not (blue triangles). PLSR of cells with speed and persistence coupling along the primary migration axis (P1 and S1) were well modeled by only a few combinations of cell-ECM interactions (Figure S6A). The optimal model for persistence in the primary direction consisted of instantaneous cell displacement, DQ, protrusion length, and percent bead movers (Figures S6A-S6D; Table S1). DQ, protrusion length, and percent bead movers, were also contributors to persistence in the PRW model, though in the APRW model the coefficient for protrusion length becomes negative, indicating that anisotropic persistent cells tend to have shorter protrusions (Table S1). Additionally, P1 relied on only one bead measurement, indicating a decreased reliance on matrix displacement in determining persistence. P2 was predicted by maximum bead displacement, protrusion rate, DQ, and protrusion lifetime.

Speed in the primary direction (S1) of APRW coupled cells is best predicted by cell-ECM interactions (Figure S6G) that are quite distinct from those that predict speed in PRW coupled cells (Figure 3J), sharing only DQ as a component. S1 instead needed maximum bead displacement, protrusion rate, protrusion lifetime, and percent movers for optimal predictive performance. S2 model components overlap strongly with S1, only substituting instantaneous bead displacement and protrusion rate for maximum bead displacement and protrusion lifetime (Figure S6O). APRW trajectories predicted from cell-ECM interactions were very similar to those simulated from MSD fits of the real trajectory data (Figure 4D).

Cell-ECM interaction state of cells with predominantly uncoupled S vs. P Cells that did not display global coupling between S and P in the PRW model or in the primary direction in the APRW model could be regressed using a distinct set of cell-ECM interactions. Only four combinations of measurements had R 2 vs. Q 2 values to predict primary persistence. This consisted of DQ, protrusion lifetime, protrusion length, instantaneous cell displacement, and bead-cell displacement ratio (Figures S7A-S7D). Cells in this regime tend to have highly linear trajectories, which were predicted well by our model (Figure 4E). High protrusion rate paired with low bead speed, instantaneous cell displacement, DQ, and protrusion lifetime predicted high P1 (Figure S7C). The PLSR model for S1 (Figures S7E-S7H) was similar to that for the APRW coupled model (Figures S7E-S7H). This indicates that in both subpopulations of the globally uncoupled cells, cell speed determination is independent from the factors that dictate cell persistence. Thus, the status of S and P coupling defines three different cell migration modes that rely on distinct combinations of cell-ECM interactions to achieve their persistence and speed (Figures 4F and 4G). Only matrix remodeling (DQ) was required to predict speed and persistence for all three populations, demonstrating the importance of this measurement when analyzing 3D migration modes.

Matrix remodeling coordination with adhesion and cell protrusion define distinct modes of cell migration

To visualize key changes between different migration modes, cell-ECM interactions were plotted for cells undergoing PRW coupled, APRW coupled, or APRW uncoupled migration. We found no significant differences in any individual measurement between these groups, although the two APRW groups displayed slight differences compared with the PRW group due to their enrichment with predominantly ITGB-1-inhibited cells (Figures 5A-5I). Analysis of the correlations among cell-ECM measurements in each group of cells revealed key changes in how subprocesses are coordinated (Figures 5J-5L). Strikingly, PRW coupled cells have strong positive correlations between multiple bead and protrusion measurements, as well as DQ and instantaneous cell displacement (Figure 5J). These contrast with the more heterogeneous correlations in the two APRW groups (Figures 5K and 5L). Notably, the strong positive correlations between protrusion rate and DQ ll OPEN ACCESS Article Developmental Cell 58, 1414-1428, August 7, 2023 1421

becomes weaker in the APRW coupled cells, and eventually negatively correlated in the APRW uncoupled cells (Figure 5M). Thus, the randomness of the walk is associated with the balance between cell protrusive activity and matrix degradation. Through this correlation analysis, distinct differences in subprocess coordination also emerged between the two APRW groups. There was a significant switch in the relationship between DQ and protrusion length (Figure 5N). APRW coupled cells showed a positive correlation between these measurements, but APRW uncoupled cells showed a negative correlation. Representative micrographs of cells from these two groups exemplify these relationships (Figure 5O). Cells that couple speed and persistence in the primary direction of migration (APRW coupled) show more matrix remodeling with more protrusive activity. Cells with uncoupled speed and persistence, on the other hand, tend to be more rounded if they are remodeling the matrix and more elongated if they are not. These patterns are enriched for ITGB1-inhibited cells. These findings prove that cells rely on different combinations of cell-ECM interactions to migrate, defined by three random walk modes of migration. Therefore, it is possible to use cell tracking data to infer cell-ECM interaction and, by extension, the overall biophysical state of the cell.

DISCUSSION

Cell migration trajectories predict associated cell-ECM interaction state using S and P coupling as a classifier In this study, we used a data-driven approach to discover how cell-ECM interactions produce heterogeneous migration behaviors. This was achieved through an integration of cell-ECM interaction reporters and cell tracking, which revealed relationships between motility subprocesses and whole-cell behavior. Although previous studies applied different walk models to describe heterogeneous migration behavior in 3D environments, 41 ours connects the walk behavior to underlying cellular subprocesses coordination. We show that confined 3D migration can be classified into three modes of speed and persistence coupling, and we developed a predictive model for each mode that links cell-ECM interaction state to cell trajectories. For all migration modes, matrix remodeling was essential, and yet no one subprocess measurement could achieve high predictive capability in isolation. This modeling approach could serve as an important tool toward integrating disparate conceptual models of cell migration that comprise the Central Dogma. For example, the adhesion-based model holds that cell speed is well predicted by a biphasic relationship with adhesion, and that this relationship can be explained by adhesion-promoted forces pushing the front and resisting motion at the rear of the cell according to a ''molecular clutch'' model of focal adhesion-actin-myosin dynamics. 6,[42][43][44][45][46][47][48][49][50][51] The actin-based model holds that speed is predicted by actin flow and its effect on polarity signals. 40 None of the well-accepted conceptual models indicate a role for matrix remodeling, likely because they originate primarily from experiments on 2D substrates or permissive 3D environments, and they do not address how migration subprocesses integrate to produce heterogeneous behaviors. Indeed, prior works seeking to understand the coupling between speed and turning (persistence) of immune cells moving in 3D in vivo 52,53 have not considered the role of ECM remodeling. Nonetheless, there is ample evidence that matrix remodeling is required for their in vivo migratory functions. For example, T cells require MMPs to penetrate infected tissue and contribute to significant tissue remodeling. 54 An important outcome of our models is the ability to reliably predict cell trajectories of multiple cell types from measurements of cell-ECM interactions. The reverse is also true, in that our models enable the prediction of cell-ECM interaction state from cell trajectories. This lends insight into how the subprocesses of migration can produce different regimes of cell speed and persistence coupling that generate migration trajectories with different extents of anisotropy. Establishing this link is important toward building multiscale models of cell migration and generating hypotheses for further mechanistic studies.

The degree of migration anisotropy is enhanced with S and P uncoupling

This work presents a conceptual model that integrates two migration frameworks: (1) speed and persistence coupling and (2) PRW models. The integration of these models revealed that speed and persistence are not universally coupled, and that the extent of coupling is directly related to the randomness of migration. Cells that follow the original UCSP law could only be accurately predicted within the PRW framework. The remaining cells were found to be either globally uncoupled or uncoupled only in the non-primary axis of migration, which can be modeled by the APRW framework. These classifications were necessary to link migration behavior with cell-ECM interaction measurements and achieve high predictive power. Integrating these two paradigms of cell migration modeling represents an (C) Graph of PRW parameters cell speed (S) vs. persistence time (P). Line indicates regression of points in the circled region to the universal coupling between speed and persistence (UCSP) equation. S and P are coupled for the majority of the cells, though a large population does not follow this trend. (D) PLSR was performed using every possible combination of measurements to fit cell-ECM interactions to persistence. The fit (R 2 ) vs. predictive ability (Q 2 ) of all possible regressions for persistence (P), and the optimal predictive model is shown by the red square. (E) The best predictive model for P accounts for more than 52% of the variance using just the first two PLS components. (F) Plotting the PLS score of each cell along these two components shows separation of persistence and the cell-ECM interaction measurement component loading projections. Loadings for the two beads displacement measurements overlap. (G) Observed vs. fitted values of z-normalized persistence show overall good agreement. (H-M) The same PLSR approach was applied for cell speed. (H) The fit (R 2 ) vs. predictive ability (Q 2 ) of all possible regressions of cell speed (S), and the optimal predictive model is shown by the red square. (I) This model can account for more than 52% of the total variance in cell speed using just the first two PLS components. (J) Plotting the PLS score of each cell along these two components shows separation of cell speed and the cell-ECM interaction measurement component loading projections. (K) Observed vs. fitted values of z-normalized cell speed show overall good agreement. (L) Trajectories simulated using the PRW from S and P values predicted by the PLSR model using cell-ECM interactions or extracted from the MSD fit. (M) recreate the heterogeneity and migration trends observed experimentally. Scale bars, 5 mm. See also Figures S4 and S5. (legend continued on next page) ll important step toward building comprehensive models of cell migration that capture heterogeneity.

Loss of adhesion uncouples S and P

We observed that many cells display coupling between speed and persistence, which has previously been described as the UCSP law and has been reproduced in multiple experimental 55,56 and theoretical models. 57,58 However, a population of cells, mostly within the ITGB1-inhibited condition, did not follow the UCSP and instead displayed a negative relationship between speed and persistence. These cells are more persistent at lower-than-expected cell speeds, implicating adhesion as a key mediator of speed and persistence coupling. Increased persistence may be consistent with the model of integrin function in which inhibition of b1 integrin activity feeds back on MT1-MMP surface levels and localization. Previous population-level studies have shown that on average, inhibition of b1 integrin activity tends to induce its association with MT1-MMP and reduce recycling, leading to accumulation of both integrin and MT1-MMP on the cell surface. 59 However, in our single-cell experiments, distinct subpopulations within the ITGB1-inhibited condition may represent cells with different levels of ITGB1-inhibition response or different initial levels of ITGB1 and MT1-MMP. The extent to which the level and activity of ITGB1 and MT1-MMP are balanced on the surface of the cell could explain the distinct coordination modes of protrusions and matrix remodeling we observed (Figures 5M-5O).

Blocking adhesion enhanced cell speed on average, indicating that the cells were originally slower migrating because they were in a higher adhesion state (Figures 1E-1H and S8). This reinforces the established biphasic relationship between adhesion and speed and echoes observations that cells having less ECM adhesion are more migratory. 60,61 In 2D, the force the cell exerts via actomyosin contraction at sites of focal adhesion determines adhesion plaque growth, maturation, and stability. 62,63 In 3D, it has been proposed that contractility must be locally balanced with ECM stiffness to stabilize adhesions. 64 When we diminished the ability of integrins to bind to the ECM, we would expect the stability of adhesions to be compromised, resulting in faster turnover and a less stable force balance between the contractile machinery and the ECM. It is interesting that in this condition in which we might expect multiple highly localized force balance instabilities, we see higher persistence of migration. Perhaps with less competition from other adhesion sites, the fewer stable adhesions drive the polarization of the cell.

The degree of S and P uncoupling is associated with the imbalance of protrusion rate and matrix remodeling Coupling between speed and persistence can be global or broken into anisotropic primary and non-primary directions of migration. Applying the PRW model distinguishes globally coupled cells vs. those that are not. Within the PRW coupled cell population, we found a strong positive correlation between protrusion rate and instantaneous cell displacement, which is consistent with the UCSP model that actin flow stabilizes cell speed. However, protrusion rate was also significantly positively correlated with DQ signal. This could suggest that protrusions physically enhanced degradation, which is supported by the fact that inhibiting actin polymerization or contractility lowered DQ signal (Figure 2H) but MMPi treatment did not significantly impact DQ signal or protrusions (Figures 2I-2K). This finding is consistent with other reports, such as how latrunculin B decreases the ability of MDA-MB-231 cells to bundle collagen I around the cell, more so than marimastat 65 and that the mechanical plasticity of the matrix can facilitate protease-independent migration. 66 However, differences in trafficking of MMPs to the membrane cannot be ruled out.

In contrast to the globally coupled cells, the uncoupled population had a negative correlation between protrusion rate and DQ. This suggests that speed and persistent coupling depends upon balanced activity between protrusion rate and matrix remodeling. High levels of one without the other leads to uncoupled migration, which leads to anisotropic migration.

When globally uncoupled cells are further divided within the APRW framework into S and P coupled or uncoupled in the primary direction, we found that the uncoupled cells drove the negative relationship between protrusion rate and matrix remodeling, whereas the APRW coupled cells had no significant correlation between these cell-ECM interactions. Therefore, there appears to be a progressive transformation of the relationship from strongly positive, to no correlation, to strongly negative as cells are globally coupled, anisotropically coupled, or uncoupled. This implicates matrix remodeling as playing a central role in determining S and P coupling. This, paired with the necessity of the DQ measurement for each of the regression models (Figures 4F and 4G), strengthens the argument that matrix remodeling is a key determinant of cell migration in confining 3D conditions.

Matrix remodeling is a necessary input to every predictive model we discovered, and its coordination with protrusive activity of cells helps differentiate between S and P coupling behaviors. The overarching conceptual model we propose captures these essential motility subprocess coordination modes and could be useful as a basis to model processes in which cell migration trajectory is important, such as cancer invasion or immune cell homing. Future studies to elucidate the molecular mechanisms underlying each coordination mode will enable the next generation of physical models of migration.

Mapping speed and persistence coupling to modes of migration The 2.5 mg/mL collagen I matrices used in this study have an average pore size of 2 mm 2 , although the 99th percentile of sizes Table S1: Coefficients from each PLSR model, related to Figure 4 Dependent variable Y is predicted by a linear combination of the z-normalized cell-ECM interaction measurements such that Y = X*B. Cell Population Y Const. Bead disp Bead speed Pct Movers Inst. Cell Speed DQ Protrusion rate Protrusion length Protrusion duration Bead-cell speed ratio All Total Disp 24.02 2.928 1.920 3.654 2.372 1.405 -2.631 PRW Linear P 0.906 0.063 0.106 0.066 0.112 0.093 -0.120 S 0.089 0.013 0.017 0.004 -0.023 APRW Linear P1 1.372 0.138 0.100 0.353 -0.277 S1 0.167 -0.041 -0.009 0.049 0.049 0.016 P2 27.23 -83.97 -14.87 41.609 39.05 S2 0.092 -0.028 -0.035 0.035 0.002 0.007 APRW Outliers P1 100.3 97.92 -1150 -629.49 -343.5 615.753 S1 0.089 -0.072 0.075 -0.253 -0.083 -0.134 P2 58.51 212.4 -163.3 -194.014 -132.37 -39.29 S2 0.068 0.006 -0.059 -0.068 0.043 -0.047