Research ArticleBIOPHYSICS

The interplay between matrix deformation and the coordination of turning events governs directed neutrophil migration in 3D matrices

See allHide authors and affiliations

Science Advances  14 Jul 2021:
Vol. 7, no. 29, eabf3882
DOI: 10.1126/sciadv.abf3882


Neutrophils migrating through extravascular spaces must negotiate narrow matrix pores without losing directional movement. We investigated how chemotaxing neutrophils probe matrices and adjust their migration to collagen concentration ([col]) changes by tracking 20,000 cell trajectories and quantifying cell-generated 3D matrix deformations. In low-[col] matrices, neutrophils exerted large deformations and followed straight trajectories. As [col] increased, matrix deformations decreased, and neutrophils turned often to circumvent rather than remodel matrix pores. Inhibiting protrusive or contractile forces shifted this transition to lower [col], implying that mechanics play a crucial role in defining migratory strategies. To balance frequent turning and directional bias, neutrophils used matrix obstacles as pivoting points to steer toward the chemoattractant. The Actin Related Protein 2/3 complex coordinated successive turns, thus controlling deviations from chemotactic paths. These results offer an improved understanding of the mechanisms and molecular regulators used by neutrophils during chemotaxis in restrictive 3D environments.


Three-dimensional (3D) cell migration is an important process during embryonic development, angiogenesis, cancer cell metastasis, immune surveillance, and inflammatory responses (15). During inflammatory responses, neutrophils, which are a subtype of white blood cells (leukocytes) and the first responders of the innate immune system, travel along chemotactic gradients initiated at or near sites of injury or infection (6). This journey involves neutrophil activation, rolling and migration within postcapillary blood vessels, transmigration across vessel walls, and 3D migration through the extravascular space guided by chemokine gradients (5). During these processes, neutrophils must negotiate a wide distribution of pore sizes and densities to translocate (7). Migration in the 3D extravascular space is particularly complex, as it requires navigating through crowded spaces composed of other cells, glycoproteins, and extracellular matrix (ECM) proteins (8, 9).

To explain how neutrophils migrate through 3D extravascular spaces, previous studies have pointed to two main maneuvers: (i) frontal protrusions and (ii) rear contractions (8, 10). Frontal protrusions are typically in the form of lamellipodia, which are an arrangement of rosettes formed by multiple, interleaved thin lamella (11). These structures form as a result of branched polymerization of actin at the front of polarized cells, which is, in turn, mediated by a set of actin binding proteins known as the Arp2/3 complex (12, 13). The Arp2/3 complex plays an important role in pathfinding during 3D leukocyte migration (11). Specifically, the protein complex facilitates lamellipodial formation for environmental probing and helps cells find neighboring pores large enough to squeeze through. Rear contractions are mediated by nonmuscle myosin-II, a motor protein that is important for cytokinesis, organelle transport, and cell migration (14). Non-muscle myosin-II has long been appreciated for its role in cell migration and has recently been found to play an important role in 3D leukocyte migration in restrictive microenvironments (15).

Efficient pathfinding depends on the ability of motile cells to change directions before committing to a path (16). After deciding on a path, cells must then translocate their bodies, including their bulky nucleus, through pores while still maintaining their compass aligned with the chemotactic gradient (17). Neutrophils are assisted during confined migration by their lobulated nuclei, which results in less steric hindrance than what would be experienced by other cell types with rounder nuclei (18). In addition, to perform 3D migration along a given path, neutrophils combine traction forces, which help them squeeze through and enlarge matrix pores, with turning to circumvent impassable pores (17, 18). These two mechanisms are challenged by increasing collagen concentration [col], which would decrease both pore size and matrix deformability. However, neutrophils are able to exhibit directed migration in a wide variety of extracellular environments, including highly restrictive ones (19).

Separately, previous studies have examined the effect of collagen density on 3D cell migration, providing insight into the relationship between changes in cell kinematics and surrounding extracellular matrices (15, 2024). However, many of these studies only considered randomly migrating cells and analyzed their motion at different time scales; this makes it difficult to separate cell speed from trajectory shape. Moreover, many of these studies focused on cancer cells, which are able to secrete matrix metalloproteinases to degrade the ECM material. Matrix degradation reduces the need for cells to bypass impenetrable sections along a migratory path (e.g., narrow pores) and therefore enables relatively straight trajectories (23). Findings from these cancer cell–based studies may therefore not be fully applicable to contexts where degradation is reduced or not possible. Last, most previous works did not quantify 3D matrix deformations induced by traction forces. The present manuscript addresses these gaps in knowledge by analyzing the trajectories of chemotactic neutrophils at different length scales, with a special emphasis on how these cells combine physical forces and turning events to maintain directional motion toward the chemoattractant source.


Automated analysis of neutrophil chemotaxis in 3D matrices

We assembled a 3D migration chamber using a well-established design (see Materials and Methods and Fig. 1A) (25). We used type I collagen in our chamber because collagen is the ECM protein most commonly found in mammals and, specifically, type I collagen is the most common type of collagen (14, 26). We used neutrophil-like dHL-60 cells differentiated from the human promyelocytic leukemia cell line [HL-60, American Type Culture Collection (ATCC)]. We used fluorescent microbeads attached to the collagen fibers to track matrix deformation (Fig. 1A, inset). We chose N-formyl-Met-Leu-Phe (fMLP), a small peptide secreted by bacteria after infection (27), as a chemoattractant.

Fig. 1 Experimental setup and migratory statistics.

(A) 3D migration chamber schematic with dHL-60 cells (gray) and fluorescent microbeads (green) embedded in a collagen gel (red). Imaging fluorescently labeled collagen showed microbeads localize to collagen fibers. (B) Image sequence of chemotaxing cell and tracked trajectory (green). (C) Probability density maps of cell coordinates (gray) and representative cell trajectories in the absence (top) and presence (bottom) of fMLP gradient (collagen, 0.75 mg/ml). Dashed circles represent unbiased random motion. In (B) and (C), time progression along each trajectory is represented by increasing color saturation, and the fMLP gradient’s direction is indicated by purple triangle. (D) Confocal reflection images of matrices with collagen concentrations [col] ranging from 0.25 to 2 mg/ml. (E) Left: Probability distributions of matrix pore radii ξ from bubble analysis. Center: Mean pore size 〈ξ〉 versus [col]. Right: Fraction of available pores with ξ ≥ 5 μm versus [col]. (F) Mean cell speeds versus matrix density (expressed as [col] or 〈ξ〉) in the presence (left) and absence (right) of fMLP gradient. (G) Chemotactic index (CI) versus matrix density. Error bars represent 95% confidence interval. All pairwise conditions showed statistically significant differences besides those labeled n.s. (not significant) in (F). *P < 0.05 according to the Mann-Whitney U test.

To study the kinematics of directed 3D neutrophil migration with high throughput, we developed a fully automated algorithm to track cell trajectories from bright-field microscopy time-lapse sequences (Fig. 1B). Overall, we tracked 21,150 cell trajectories under different conditions (see table S1 for details). We confirmed that our system enabled directed neutrophil migration by performing experiments with and without the addition of fMLP. Under both conditions, we mapped the probability density p(x, y) of the cell coordinates in the directions parallel (y) and perpendicular (x) to the chemoattractant gradient (Fig. 1C). For illustrative purposes, we plotted the six longest trajectories, which were long enough to be representative of the qualitative features of cell paths in our experiments. In the absence of fMLP, the distribution of cell trajectories was symmetric, leading to circular contours in p(x, y), as expected for random migration. In contrast, the chemoattractant gradient created by adding fMLP resulted in cell trajectories preferentially directed toward the chemoattractant source, with a p(x, y) being markedly asymmetric and shifted toward the positive y direction. These data confirm that our system recapitulated 3D neutrophil chemotaxis and that our cell-tracking algorithm was able to capture this process from bright-field microscopy, without the need for fluorescent labeling.

Collagen density affects the cell speed but not the chemotaxis of 3D neutrophils

In physiological contexts, neutrophils must travel across diverse microenvironments. To study the effect of collagen density ([col]) on the ability of neutrophils to undergo directed migration, we fabricated gels with 0.25 mg/ml ≤ [col] ≤2.0 mg/ml. Varying [col] has been previously shown to affect the matrix microenvironment, most notably its pore size (20, 28, 29). Thus, we visualized the structure of our collagen gels (Fig. 1D) and quantified their distribution of pore radii (ξ) (Fig. 1E). This distribution became narrower and shifted to lower values of ξ with increasing [col] (Fig. 1E, left). In particular, 〈ξ〉 decreased monotonically from 4.4 to 2.3 μm when [col] was increased from 0.25 to 2.0 mg/ml (Fig. 1E, middle). We took ξ ~ 5 μm as an approximate value for the radius of dHL-60 cell nuclei (30) and calculated p(ξ ≥ 5 μm) to assess the availability of pores through which cells could easily squeeze through without substantial nuclear deformation. This calculation showed that the fraction of pores larger than ξ ~ 5 μm was more than 40% in matrices (0.25 mg/ml), but it decreased to less than 10% as the [col] increased to 2 mg/ml (Fig. 1E, right). Overall, these data highlight that the different collagen concentrations considered in our experiments allowed us to study 3D neutrophil-like migration through microenvironments that ranged from nonrestrictive to highly restrictive.

Next, we were interested in how these different matrix microenvironments affected directed 3D neutrophil migration. We found that cells migrated with mean speeds that decreased from 〈v〉 = 8.4 μm/min to 〈v〉 = 5.7 μm/min as [col] was increased and 〈ξ〉 decreased (Fig. 1F, left). A similar dependence of 〈v〉 on 〈ξ〉 was observed for randomly migrating cells in chemotaxis devices without fMLP (Fig. 1F, right). Despite migrating more slowly as [col] was increased, cells were still able to undergo directed migration toward an added fMLP source in all matrices investigated (Fig. 1G, left). This is in contrast to migration in chambers without fMLP (Fig. 1G, right). Overall, these results indicate that decreasing matrix pore size and availability has a marked effect on the migration speed of wild-type neutrophil-like cells, but not on their ability to perform chemotaxis.

3D neutrophil migration requires turning but not mechanical forces to negotiate restrictive matrices

In addition to reducing pore size, increasing [col] from 0.25 to 2 mg/ml is known to cause significant matrix stiffening (i.e., bulk Young’s modulus raised from 440 to 5200 Pa) (29). Thus, we were interested in finding whether changes in [col] affected the mechanical interactions between neutrophils and their environment. To this end, we imaged fluorescently labeled, chemotaxing cells in our custom-made devices along with fluorescent microbeads that were attached to matrix fibers. We used Airyscan confocal microscopy to image both the cells and microbeads in collagen gels (0.25 and 2 mg/ml) (Fig. 2A, left). Matrix deformations were evident in gels (0.25 mg/ml) by comparing the location of beads at any given time point (green) with that at the previous time point (orange). The 3D incremental matrix deformations u=(u,v,w) were quantified using Particle Image Velocimetry (PIV) (Fig. 2A). Our measurements showed clear deformations in the matrices (0.25 mg/ml), with a magnitude (~3 μm) comparable to the cell nucleus size. These data suggest that, in addition to navigating the matrix searching for large pores, cells in low-[col] matrices are able to significantly deform and enlarge matrix pores (Fig. 2B). In contrast, cells in high-[col] matrices did not exert forces large enough to cause appreciable matrix deformations (Fig. 2C). These observations were confirmed by our statistical analysis of 3D matrix deformation moments generated by neutrophils for [col] = 0.25, 0.75, and 2 mg/ml. By extension of 2D definition of Butler et al. (31), the matrix deformation moments were defined as tr(M), where the elements of the M tensor are the first-order moments of vector u in a sphere of 200 μm surrounding the cell (see Eq. 2 in Materials and Methods). The statistics of the absolute value ∣tr(M)∣ indicated a 2.5-fold decrease in matrix deformations between [col] = 0.25 and 0.75 mg/ml and a 10-fold decrease between [col] = 0.25 and 2 mg/ml (Fig. 2D, left). Furthermore, we found the sign of tr(M) to be preferentially negative regardless of [col], indicating that cells preferentially pulled the matrix toward them (fig. S1A). Projecting the M tensor in the directions parallel (∥) and perpendicular (⊥) to the fMLP gradient (Fig. 2D, right) revealed that cell contractions were significantly larger in ∥ the direction than in the ⊥ direction for [col] = 0.25 and 0.75 mg/ml; however, these contractions became isotropic in the more restrictive [col] = 2 mg/ml of matrices.

Fig. 2 Mechanical remodeling of surrounding matrix for low- and high-density collagen matrices.

(A) Top: Reconstruction of a fluorescently labeled cell (gray) and fluorescent microbeads at t = (i − 1)∆t (orange) and t = it (green) during chemotactic migration in collagen gel (0.75 mg/ml); ∆t = 30 s. Bottom: Absolute magnitude of incremental deformations computed by 3D particle image velocimetry. (B) Representative time-lapse reconstructions of cells (gray), beads (green), and incremental 3D matrix deformation (magenta-blue color map) during migration in a low [col] matrix (0.25 mg/ml). (C) Same as (B) for high [col] matrix (2 mg/ml). (D) Left: Mean magnitude of 3D total matrix deformation moment, ∣tr(M)∣, versus collagen density. Right: Projections of the deformation moment in the directions parallel (∣M∣, colored bars) and perpendicular (∣M∣, hollow bars) to the fMLP gradient. *P < 0.05 and **P < 0.01 according to the Mann-Whitney U Test.

Notably, the representative cell in the collagen matrix (2.0 mg/ml) initially moved toward the chemoattractant source but seemed unable to continue moving in that direction at t = 1 min (Fig. 2C). Subsequently, it performed a sharp, ~π/2-radian turn and continued to migrate in the direction perpendicular to the fMLP gradient for the next minute. In comparison, the representative cell in the matrix (0.25 mg/ml) was able to sustain directed motion toward the fMLP source for t = 2.5 min. These observations support the idea that as [col] increases, these cells must adapt their migration strategy to circumvent the more frequent impassable, nondeformable sections of the matrix.

Neutrophils migrating in 3D matrices retain overall chemotactic migration despite turning more frequently as matrix pore size decreases

The results presented thus far led us to hypothesize that turning events are a major factor that enabled the directed migration of neutrophils in restrictive 3D matrices at the expense of decreasing speed. To examine this hypothesis, we aimed to statistically characterize the cell trajectories in matrices with different collagen concentrations. We first examined the probability distribution p(x, y) of cell coordinates along their trajectories for different collagen concentrations (Fig. 3A). In all the matrices investigated, the trajectories were biased toward the chemoattractant source, and the maps of p(x, y) had eccentric, elongated contours in the y direction. However, the trajectories were shorter, and the p(x, y) contours extended less into the y > 0 plane as [col] increased, in agreement with the decreasing cell speed (Fig. 1F). Qualitative inspection of cell trajectories suggested that cells perform sharper turns with increasing [col]. Visually, this trend was most evident at several points along the longest trajectories (see arrows in Fig. 3A). To characterize turning events statistically, we computed the autocorrelation function of the vector tangent to the trajectory in intrinsic coordinates, Ctt=t(s)t(s+Δs) (see diagram in Fig. 3B). This function was plotted in Fig. 3C for chemotaxing and randomly moving cells in different matrices, showing that Ctt decayed with increasing Δs from its value of Ctt(0) = 1 at the origin as the cells added turns to their trajectories. For chemotaxing cells, this decrease was more pronounced in high-[col] matrices than in low-[col] ones. By contrast, randomly moving cells showed a weak dependence of Ctt with [col].

Fig. 3 Persistence in cell trajectories for varying collagen densities.

(A) Probability density maps of cell coordinates (gray contours) and representative cell trajectories for nonrestrictive ([col] = 0.25 mg/ml), intermediate ([col] = 0.75 mg/ml), and restrictive ([col] = 2 mg/ml) matrices. Arrows indicate representative turns. Time progression along each trajectory is represented by increasing color saturation. Dashed circles represent unbiased random motion. (B) Schematic of a cell trajectory, indicating a segment of length ∆s, along which the tangent vector t changes its orientation θ with respect to the gradient by Δθ. (C) Tangent vector autocorrelations versus distance traveled, Ctt(∆s), for no fMLP (dashed blue) and fMLP (red solid) conditions, in matrices of varying density. Dashed line marks the cutoff Ctt = 0.05 used to determine the trajectories persistence (Sp). (D) Mean squared displacements (MSDs) versus distance separation along cell trajectories, r2(∆s). Dashed and solid lines represent r2 ~ ∆s2 and r2 ~ ∆s, respectively. (E) Sp of cell trajectories versus matrix mean pore size, 〈ξ〉, for no fMLP (blue squares) and fMLP (red circles) conditions. Error bars Sp are superimposed on top of data points and represent 95% confidence intervals of the mean. *P < 0.05 and **P < 0.01 according to the Mann-Whitney U test.

Next, we determined the mean squared displacements (MSDs) of the cells’ trajectories. The MSD versus time, i.e., r(t)=x(t+t)x(t)2, have been widely used to quantify 2D and 3D cell migration by exploiting the analogy between these processes and a random walk (21, 23, 3234). In particular, the logarithmic slope of the MSD plot is often used to assess whether migration is random or directed. This time-domain analysis has proven useful to study rate-limiting processes governing cell migration. However, considering this study’s focus was analyzing the shape of cells trajectories affected by physical confinement in 3D matrices, we considered the MSDs as a function of distance traveled, i.e., in intrinsic coordinates r(s)=x(s+s)x(s)2. Intrinsic coordinates are widely used in differential geometry and allow for studying the geometry of cell trajectories independent of cell speed. This consideration is relevant because cell speed can depend on [col] and varies significantly along each cell’s trajectory as the cell turns to avoid impassable sections of the matrix (shown below). The rate of increase in the MSDs with ∆s depends on whether cells maintain directed translocation over a distance equal to ∆s. Note that since the MSDs can be expressed asr(s)=0Δs0ΔsCtt(Δs)dΔsdΔs(1)the MSDs are tightly linked to the autocorrelation of the direction of cell migration. The MSDs of chemotaxing and randomly moving cells were plotted in Fig. 3D for various collagen concentrations. For short travel distances, r(∆s) ~ Δs2 under all conditions tested, a hallmark of directed motion that is consistent with the Ctt ≈ 1 values that are shown in Fig. 3C. Subsequently, the MSDs tapered toward r(∆s)~Δs as Δs increased, and the cells began to add changes of direction (turning events) along their paths as Ctt progressively dropped. The MSD and Ctt plots in the time domain are represented in fig. S2 (A and B) for reference.

The rate of decay of Ctt governs the transition between the linear and quadratic regimes of the MSDs according to Eq. 1 and is related to the persistence length of the trajectories, Sp, (35). Thus, we estimated Sp using an exponential decay model for Ctt (see Material and Methods for details) and plotted the results in Fig. 3E. These data show that the persistence of the paths of neutrophils undergoing directed migration through 3D matrices plateaued at Sp ≈ 20 μm for matrices with 〈ξ〉 ≳ 3 μm and then fell sharply to reach Sp ≈ 7 μm for 〈ξ〉 < 3 μm . In contrast, randomly moving neutrophils had trajectories with Sp ≈ 5 μm regardless of 〈ξ〉. The persistence times obtained from Ctt in the time domain follow a similar trend with 〈ξ〉 (fig. S2C). Together with the matrix deformation measurements presented in Fig. 2, these results suggest that neutrophils undergoing chemotaxis perform more frequent turns in restrictive matrices with narrower pores that cannot be easily deformed. In our most restrictive matrices, these turns became almost as frequent as those performed by randomly moving cells.

Previous studies have uncovered key features of 3D random cancer cell migration by characterizing the probability distribution of turning angles in cell paths (23). Thus, we determined the distribution of Δθ versus Δs for the chemotactic neutrophil trajectories tracked in our 3D experiments, i.e., p(Δθ|Δs). Figure 4 (A and B) displays contour maps of p(Δθ|Δs) in collagen matrices (0.5 and 2 mg/ml). They show that Δθ is narrowly distributed around Δθ = 0 for small Δs, implying that the neutrophils followed relatively straight paths for short spatial separations. With increasing Δs, the distributions become progressively more uniform as the neutrophils implement changes of direction to navigate the matrix and their paths lose directional persistence. As expected from the Ctt and MSD data described above, the loss of persistence took place at shorter distances with increasing [col]. However, note that, while the distribution of Δθ became quasi-uniform at Δs ~ 20 μm for randomly moving neutrophils (Fig. 4, A and B, insets), chemotaxing neutrophils maintained long-range directionality, as denoted by p(Δθ|Δs > 100 μm) that still maintains a peak at Δθ = 0. This persistence was also observed in the most restrictive matrices, despite the cells turning more often in these matrices.

Fig. 4 Statistics of turning angles and cell speeds versus matrix density.

(A and B) Probability distribution maps of turning angles ∆θ conditional to distance separation ∆s, p(∆θ∣∆s). (C and D) Probability distribution maps of ∆θ conditional to the orientation θ with respect to the gradient, for distance separation equal to persistence length, p(∆θ∣θ, ∆s = Sp). (E and F) Maps of p(∆θ∣θ, ∆s = 5Sp). (G and H) Joint probability distribution maps of turning angle and cell speed, for ∆s = Sp/2,i.e., p(∆θ, v, ∆s = Sp/2). The maps in (A), (C), (E), and (G) correspond to [col] = 0.25 mg/ml, whereas those in (B), (D), (F), and (H) correspond to 2 mg/ml. In each panel, the main probability map comes from fMLP conditions, while the inset comes from randomly moving cells in the absence of fMLP. The dotted and solid lines in (C) to (F) indicate θ = ∆θ and θ = ∆θm(θ), respectively, where ∆θm is the mode of the ∆θ distribution. The solid and dotted lines in (G) and (H) represent the mean cell speeds of chemotaxing and randomly moving cells versus turning angle for s=Sp2, i.e., v(∆θ, ∆s = Sp/2). The solid circles in (G) and (H) indicate the total mean speed 〈v〉 for each matrix. In (H), the dashed red line represents the mean speed of chemotaxing cells in matrices (0.5 mg/ml).

Chemotaxing neutrophils integrate short-range displacements dictated by diverse 3D microenvironments to achieve overall chemotactic migration

The data presented above suggest that the persistence length of neutrophils migrating in 3D is dictated by small-scale features of the matrix such as pore size. To investigate how these cells balance the short-range persistence imposed by their microenvironment with the long-range directional bias necessary for chemotaxis, we calculated the probability distribution of their turn angles conditional to their prior orientation with respect to the gradient and distance separation, i.e., p(Δθ|θ,Δs) where θ=cos1[t(0,1)] (see diagram in Fig. 1F). Each horizontal section of these “topographical” probability maps is the probability of turning an angle Δθ of the cells that are moving at an angle θ with respect to the chemoattractant (see fig. S3 for further details about this representation). Figure 4 (C and D) shows maps of p(Δθ|θ,Δs) for a fixed distance separation Δs = Sp for [col] = 0.5 and 2 mg/ml. Because of the definition of Sp, this distance separation corresponds roughly to one turn on average. In both matrices, p(Δθ|θ,Δs = Sp/2) is independent of θ, indicating that neutrophils moving both toward and away from the gradient have similar turning angle distributions for distances comparable to the persistence length. This is, by no means, a trivial result; it would have been tempting to speculate that cells moving along the gradient (θ = 0) would turn less sharply than those moving perpendicular to the gradient (θ = π/2). This is the case in 2D chemotaxis (36) but not in our 3D experiments. Furthermore, the angular distributions of chemotaxing neutrophils at ∆sSp were similar to those for randomly moving ones at this distance separation [see insets in Fig. 4 (C and D) and fig. S3A].

Next, we examined the conditional Δθ distributions of turning angles for longer segments of trajectory, i.e., Δs = 5Sp, which likely concatenate multiple turns (Fig. 4, E and F, and fig. S3B). Regardless of [col], these distributions peaked near Δθ = θ, in evident contrast with the distributions obtained from randomly moving cells (Fig. 4, E and F, insets), which were roughly uniform and independent of θ. Overall, our data are consistent with a model in which 3D migrating neutrophils implement turns in their trajectories to circumvent impassable sections of the matrix. These short-range turns are dictated by the local microenvironment and result in misalignments with respect to the gradient that are not immediately corrected. Instead, it takes several consecutive turns for neutrophils to realign themselves with the gradient in both nonrestrictive and restrictive matrices. We examine this realignment process in more detail later in this study.

Neutrophils migrating in 3D matrices decrease their speed when turning

To study whether turning events and cell migration speed are interdependent, we generated probability distributions of cell velocity and Δθ as a function of the distance traveled in our different matrices, i.e., p(Δθ,vs). We observed a clear dependence of v on Δθ when we considered a segment with approximately one turn by choosing Δs comparable to Sp. To illustrate this dependence, we plotted 2D probability density maps p(Δθ,vs) for fixed Δs = Sp/2 (Fig. 4, G and H). These data revealed that neutrophil speeds were reduced significantly while turning and that this reduction was greater for sharper turns (i.e., higher Δθ). When [col] increased from 0.5 to 2 mg/ml, cell speed decreased further (compare red dashed and black solid lines in Fig. 4H), and this decrease was more pronounced when the cells were turning. Given that cells turned more often as [col] increased, the difference between the speed of straight motion [v(Δθ = 0)] and the mean cell speed (solid circles in Fig. 4, G and H) was larger in high-[col] matrices than in low-[col] ones.

Since the length of a curved trajectory is longer than its secant line, tracking migrating cells with a finite-time resolution could underestimate their turning velocity and create a mistaken appearance that cell speed decreases with turning angle. To confirm that this artifact was not present in our data at ∆t = 30 s, we performed additional experiments on cells migrating through matrices (0.5 mg/ml) under fMLP gradient and tracked these cells with a much finer temporal resolution of ∆t = 1 s. Given that Sp ≈ 20 μm and v ≈ 7 μm/ min under these conditions, each turning event was resolved by >170 data points on average at ∆t = 1 s, which should be enough to avoid artifacts. The p(Δθ, vs = Sp/2) map from the ∆t = 1-s experiments was very similar to that obtained with ∆t = 30 s (fig. S4A). To further examine whether the observed dependence of cell migration speed on turning angle was a consequence of cell-to-cell heterogeneity or if each single cell slowed down when turning, we plotted v(s) versus Δθ(s) or single cells from the ∆t = 1-s experiments (fig. S4B). These data clearly show that each cell migrates slower when turning.

Short-range persistence impedes overall chemotactic migration in restrictive matrices

After finding that neutrophils increase their frequency of turning in restrictive matrices, thereby decreasing their migratory persistence, we examined how perturbations leading to increased persistence would affect 3D cell migration in different microenvironments. To this end, we treated cells with the Actin Related Protein 2/3 complex (Arp 2/3) inhibitor ck666, which restricts the ability of the Arp2/3 complex to polymerize branched actin filament networks and form lamellipodial protrusions (37, 38). Previous studies have postulated that Arp2/3-mediated lamellipodia formation is crucial for exploring the local microenvironment during 3D migration and finding the best paths (pores) to migrate through (11, 16). These studies have also shown that the inhibition of the Arp2/3 complex leads to cell trajectories that are more persistent in space.

We first examined the probability distributions of cell positions after ck666 treatment and found a consistent trend of decreased traveled distances and straighter trajectories compared to untreated cells, even for high [col] (compare Figs. 5A and 3A). To ensure the appropriateness of our untreated controls with ck666-treated cells, we performed vehicle control experiments and verified that the ck666 vehicle dimethyl sulfoxide (DMSO) did not have a significant effect on dHL-60 cell kinematics (fig. S5). Consistently, the turning angles of ck666-treated cells were narrowly distributed near Δθ = 0 for Δs < 100 μm regardless of [col] (compare Figs. 5A and 4, A and B). Relative to untreated cells, the MSDs versus ∆s of the ck666-treated cells remained relatively close to r(∆s) ~ Δs2 (Fig. 5B), and Ctt decayed more slowly for these cells than for untreated cells (Fig. 5C). Consequently, the persistence distance of the migrating ck666-treated cells was generally higher than that of untreated cells (Fig. 5D). Similar to untreated cells, we observed a transition toward shorter Sp with decreasing ξ. Furthermore, this transition seemed to occur for larger ξ than in untreated cells, suggesting that the migration of Arp2/3-inhibited cells was disrupted by less restrictive matrices than the migration of untreated cells. However, ck666-treated cells were unable to reduce their persistence as much as untreated cells to navigate matrices with small pore sizes. For instance, Sp(ξ = 2.3 μm) was twice as long for ck666-treated than for untreated cells (14 μm versus 7 μm; see Fig. 5D).

Fig. 5 Migration persistence for Arp2/3-inhibited cells.

(A) Probability density maps of cell coordinates and representative cell trajectories together with maps of p(∆θ∣∆s) for nonrestrictive ([col] = 0.25 mg/ml), intermediate ([col] = 0.75 mg/ml), and restrictive ([col] = 2 mg/ml) matrices. The data come from chemotaxing cells treated with the Arp2/3 inhibitor ck666. Time progression along each trajectory is represented by changes from white to saturated colors. Dashed circles represent unbiased random motion. (B) MSDs of ck666-treated cells in matrices of varying density, r2(∆s). Dashed and solid lines represent r2 ~ ∆s2 and r2 ~ ∆s, respectively. (C) Ctt(∆s) of ck666-treated cells in matrices of varying collagen density. The gray dashed line marks the cutoff Ctt = 0.05 used to determine Sp. In (B) and (C), dashed dark red lines come from untreated cells in matrices (2 mg/ml). (D) Persistence lengths of cell trajectories, Sp, as a function of matrix mean pore size, 〈ξ〉, for ck666-treated cells (purple triangles) and untreated cells (red circles). Error bars for persistence length calculations are superimposed on top of data points and represent 95% confidence intervals of the mean.

On average, ck666 treatment caused a decrease in neutrophil speed, and this decrease was more accentuated for high [col] (Fig. 6A). A more detailed analysis of neutrophil trajectories revealed that, similar to untreated neutrophils, the speed of ck666-treated cells dropped mildly with Δθ for [col] = 0.25 mg/ml (Fig. 6B, left) and sharply for [col] = 2.0 mg/ml (Fig. 6B, right). Notably, the largest differences in speed between untreated cells and ck666-treated cells occurred at low turning angles in both matrices (compare black lines and red dashed lines in Fig. 6B), implying that Arp2/3 inhibition not only impaired pathfinding but also slowed down straight motion.

Fig. 6 Migration statistics for Arp2/3-inhibited cells.

(A) Mean speeds <v> of ck666-treated cells versus matrix density. The red dashed line represents <v> from untreated cells. The error bars represent 95% confidence interval. (B) Maps of p(∆θ, v, ∆s = Sp/2) for ck666-treated cells in nonrestrictive ([col] = 0.25 mg/ml; left) and restrictive ([col] = 2 mg/ml; right) matrices. Solid and dashed lines represent mean cell speed v(∆θ, ∆s = Sp/2) for ck666-treated cells and untreated cells, respectively. The dotted line in the right panel represents v(∆θ, ∆s = Sp/2) for ck666-treated cells in the matrix (0.25 mg/ml). (C) CI of ck666-treated cells versus matrix density. Red dashed line represents CI from untreated cells. The error bars represent 95% confidence interval. (D) Maps of p(∆θ∣θ, ∆s = 5Sp) for ck666-treated cells in sparse (0.25 mg/ml; left) and dense (2 mg/ml; right) matrices. The dotted and solid lines indicate respectively θ = ∆θ and θ = ∆θm(θ). The dashed red lines correspond to ∆θm(θ) from untreated cells. All pairwise conditions showed statistically significant differences besides those labeled n.s. in (A). Otherwise, *P < 0.05 according to the Mann-Whitney U test.

Next, we asked how the increased short-range persistence (i.e., higher Sp) caused by Arp2/3 inhibition affected chemotactic efficiency, as measured by chemotactic index (CI). Our results (Fig. 6C) indicated that ck666-treated cells had a higher CI than untreated cells for low [col]. However, the CI of ck666-treated cells decreased with [col], falling below the CI of untreated cells in our three highest-[col] matrices. In an effort to understand these trends, we analyzed the p(Δθ|θ,Δs = 5 Sp) for ck666-treated cells in the matrices (0.5 and 2 mg/ml) (Fig. 6D). Our analysis showed that, despite their decreased rate of turning per unit length of their trajectory (i.e., increased Sp), the Δθ distributions peaked near Δθ = θ in the matrix (0.25 mg/ml) (Fig. 6D, left). Thus, ck666-treated cells were able to undergo sufficient turning to maintain directional bias toward the chemoattractant source in nonrestrictive matrices, similar to untreated cells. However, the Δθ distributions from ck666-treated cells were notably different in restrictive, high-[col] matrices (Fig. 6D, right). In particular, they became wider with an increasing misalignment angle θ but peaked near Δθ = 0 regardless of θ. Data from other values of Δs yielded similar results. These findings suggest that the ck666-treated cells underwent decreased biased turning and were not able to correct their misalignment with respect to the chemotaxis axis in restrictive matrices.

Myosin-II contractility is not crucial for turning but facilitates fast 3D neutrophil migration

Previous studies have shown that myosin-II–mediated rear contractility is an important mechanism in 3D cell migration (39). Inhibition of myosin-II significantly decreases the ability of neutrophils to squeeze through narrow constrictions (15). Given that our data for both untreated and ck666-treated cells suggested an important role for turning in maintaining the long-range directional bias required for chemotaxis, we examined whether myosin-II contractility was necessary for turning. To investigate this, we treated cells with the nonmuscle myosin-II inhibitor blebbistatin, which disrupts myosin-II–mediated cell contractility (40).

The distributions of cell positions after blebbistatin treatment (Fig. 7A) suggested decreased traveled distances when compared to untreated cells (Fig. 3A) for all collagen concentrations; however, a comparison of the six longest cell trajectories did not reveal noticeable differences in turning events. As with untreated cells, the turning angles for blebbistatin-treated cells were narrowly distributed around Δθ = 0 for small Δs and progressively became more uniform with increasing Δs (Fig. 7A). The MSDs versus ∆s results for blebbistatin-treated cells were also similar to those for untreated cells in that r(∆s) ~ Δs2 for small Δs and tapered to r(∆s)~Δs as Δs increased (Figs. 7B and 3D, respectively). However, careful analysis of the persistence of blebbistatin-treated cells from their Ctt(∆s) (Fig. 7C) revealed that decreased contractility made cells more sensitive to their microenvironment, as the transition to an increased rate of turning occurred for larger ξ and lower [col] than in untreated cells (Fig. 7D).

Fig. 7 Migration persistence for myosin-II–inhibited cells.

(A) Probability density maps of cell coordinates and representative cell trajectories together with maps of p(∆θ∣∆s) for nonrestrictive ([col] = 0.25 mg/ml), intermediate ([col] = 0.75 mg/ml), and restrictive ([col] = 2 mg/ml) matrices. The data come from chemotaxing cells treated with the myosin-II inhibitor blebbistatin. Time progression along each trajectory is represented by changes from white to saturated colors. Dashed circles represent unbiased random motion. (B) MSDs of blebbistatin-treated cells in matrices of varying collagen density, r2(∆s). The dashed and solid lines represent r2 ~ ∆s2 and r2 ~ ∆s, respectively. (C) Ctt(∆s) of blebbistatin-treated cells in matrices of varying density. The gray dashed line marks the cutoff Ctt = 0.05 used to determine Sp. In (B) and (C), dashed dark red lines come from untreated cells in the matrix (2 mg/ml). (D) Persistence lengths of cell trajectories, Sp, as a function of matrix mean pore size, 〈ξ〉, for blebbistatin-treated cells (green diamonds) and untreated cells (red circles). Error bars for persistence length calculations are superimposed on top of data points and represent 95% confidence intervals of the mean.

Blebbistatin treatment resulted in decreased average neutrophil speeds that continued to decrease with [col] (Fig. 8A). This drop in speed was nearly constant across all turning angles (Fig. 8B). The CI of blebbistatin-treated cells did not vary with [col] and was slightly higher than that for untreated cells (Fig. 8C). Consistent with the CI data, the long-range (i.e., Δs = 5Sp) Δθ distributions conditional to θ did not vary qualitatively after blebbistatin treatment in both nonrestrictive and restrictive matrices (Fig. 8D). Together, these results suggest that myosin-II contractility is important for fast 3D neutrophil migration and affects the transition to decreased short-range persistence as [col] increases. However, neutrophils are able to navigate restrictive matrices by coordinating consecutive turning events to maintain long-range directional bias with or without fully active myosin-II contractility.

Fig. 8 Migration statistics for myosin-II–inhibited cells.

(A) Mean speeds <v> of blebbistatin-treated cells versus matrix density. The red and purple dashed lines represent <v> from untreated and ck666-treated cells, respectively. Error bars represent 95% confidence interval. (B) Maps of p(∆θ, v, ∆s = Sp/2) for blebbistatin-treated cells in sparse (0.25 mg/ml; left) and dense (2 mg/ml; right) matrices. The solid and dashed lines represent the mean cell speed v(∆θ, ∆s = Sp/2) for blebbistatin-treated cells and untreated cells, respectively. The dotted line in the right panel represents v(∆θ, ∆s = Sp/2) for blebbistatin-treated cells in the matrix (0.25 mg/ml). (C) CI of blebbistatin-treated cells versus matrix density. The red and purple dashed lines represent CI from untreated cells. Error bars represent 95% confidence interval. (D) Maps of p(∆θ ∣θ, ∆s = 5Sp) for blebbistatin-treated cells in sparse (0.25 mg/ml; left) and dense (2 mg/ml; right) matrices. The dotted and solid lines indicate respectively θ = ∆θ and θ = ∆θm(θ). The dashed red lines correspond to ∆θm(θ) from untreated cells. All pairwise conditions showed statistically significant differences besides those labeled n.s. in (A). Otherwise, *P < 0.05 according to the Mann-Whitney U test.

Directional bias depends on Arp2/3 and requires on average three consecutive turning events, regardless of matrix density and myosin-II contractility

Neutrophils chemotaxing on 2D surfaces have been shown to immediately correct their misalignments with respect to a chemotactic gradient in one turn (36). However, our data suggest that, in 3D matrices, these cells need a longer distance to achieve the same objective (Fig. 3). We denoted this distance Sθ (see diagram in Fig. 9A) and estimated it in two steps. First, we obtained p(Δθ,vs) as described above and calculated the difference between the most likely turning angle in this distribution, ∆θm(θ, ∆s), and a full directional correction [i.e., ∆θm(θ, ∆s) = θ] for each value of ∆sdθm(s)=π20πθm(θ,s)θdθ(2)

Fig. 9 Cell realignment with the chemoattractant gradient.

(A) Top: Schematic indicating the bias length scale, Sθ, defined as the distance that cells need to correct an angular misalignment θ with respect to the chemoattractant gradient. Bottom: Metric d∆θm (equation in figure) that quantifies how close cells are to a total alignment correction ∆θ = θ as a function of distance separation ∆s. The data come from untreated chemotaxing cells in matrices of different densities. The gray dashed line indicates the threshold used to determine Sθ. (B) Left: Sθ versus 〈ξ〉 for different treatments and matrix densities. Error bars for bias length scale calculations are superimposed on top of data points and represent 95% confidence intervals of the mean. Right: Sθ versus Spfor different treatments and matrix densities. The dashed dark and light lines represent Sθ = 2Sp and Sθ = 4Sp respectively. Untreated cells, ck666-treated cells, and blebbistatin-treated cells are represented by red circles, purple triangles, and green diamonds respectively, and their color indicates collagen density (0.25, 0.5, 0.75, 1, and 2 mg/ml from light to dark).

It is straightforward to show that this definition yields d∆θm = 1/2 for a narrow (e.g., Dirac delta) distribution with ∆θm = 0. Consistently, we found that d∆θm (∆s) decreases with ∆s from its value d∆θm (0) ≈ 0.5 at the origin (Fig. 9A, bottom). The second step to determine Sθ was to use the intersection between the d∆θm (∆s) curve and the threshold value of 1/4 to estimate Sθ for each [col], both for untreated and pharmacologically manipulated cells. We chose the 1/4 threshold because d∆θm = 1/4 for a uniform distribution with ∆θm = π/2.

The results of this analysis (Fig. 9B, left) indicated that untreated cells and blebbistatin-treated cells were able to correct misalignments with respect to the gradient over a distance Sθ ≈ 50 μm in nonrestrictive, low-[col] matrices with 〈ξ〉 > 3 μm. However, this distance decreased to Sθ ≈ 30 μm in restrictive, high-[col] matrices with 〈ξ〉 < 3 μm. To investigate this seemingly counterintuitive result that cells can restore their chemotactic compass over a shorter distance in more restrictive matrices, we plotted Sθ versus the persistence distance Sp (Fig. 9B, middle). We found that Sθ spanned a relatively narrow range when scaled with Sp [Sθ ≈ (2 − 4) Sp]. This result implies that it takes neutrophils between two and four turns on average to realign their motion with the chemoattractant gradient. Furthermore, since cells turn more often in more restrictive matrices, they need a shorter physical distance to realign. Our data suggest that [col] or myosin-II contractility did not seem to affect this result, but ck666 inhibition of the Arp2/3 complex markedly increased Sθ (Fig. 9B, middle). Thus, ck666-treated cells that performed turns to avoid obstacles in their restrictive microenvironments had a reduced ability to control the direction of turning to remain aligned with the chemotaxis gradient.

Given that Sp and Sθ can be loosely interpreted as the frequencies per unit distance of random and biased turns along a cell’s trajectory, we hypothesized that the chemotactic efficiency of the cells would depend on the ratio Sp/Sθ. To test this hypothesis, we plotted the average of cosθ versus Sp/Sθ for all our experimental conditions (Fig. 9B, right). These data confirm that cell motion remained more aligned with the chemoattractant gradient for larger values of Sp/Sθ. In addition, we built mathematical models of biased random walks with persistence distance Sp and bias length scale Sθ arising from an angular controllerθ·(s)=1Sθsinθ(s)+Γ(s)(3)where Γ(s) is a white noise source with Γ(s)Γ(s)=2Spδ(ss)[see the Supplementary Materials and (41)]. The simplest version of these models, which assumes constant cell speed independent of θ, is particularly interesting because it has an exact analytical solution for the average of cosθ that does not have any adjustable parameter, i.e., cosθ=I1( Sp/Sθ)I0( Sp/Sθ), where I1 and I0 are modified Bessel functions of the first kind. This model prediction agrees well with our experimental data (Fig. 9B, right).

Cell-generated mechanical forces affect the rate of turning and speed of neutrophils migrating in 3D, but not their chemotactic efficiency

We examined to what extent neutrophils relied on mechanical forces to migrate in 3D environments by quantifying the matrix deformations caused by untreated, ck666-treated, and blebbistatin-treated cells. In our matrices (0.25 mg/ml), both treatments caused a notable decrease in 3D matrix deformations (Fig. 10, A to C) that was statistically significant in the case of blebbistatin treatment as measured by the 3D matrix deformation moment ∣tr(M)∣ (Fig. 10D), with the matrix deformation first-order tensor M defined asM=RxuT dx(4)

Fig. 10 Relationship between cell-generated matrix deformations and cell migration parameters.

(A to C) Representative reconstruction of cells (gray), beads (green), and incremental 3D matrix deformation (magenta-blue color map) during migration in a sparse (0.25 mg/ml) collagen matrix for no treatment (A), blebbistatin treatment (B), and ck666 treatment (C). (D to F) Mean magnitude of 3D contractile moment, ∣tr(M)∣, for different treatments in sparse (0.25 mg/ml) (D), intermediate (0.75 mg/ml) (E), and dense (2 mg/ml) (F) matrices. Error bars represent 95% confidence intervals of the mean. *P < 0.05 and **P < 0.01 according to the Mann-Whitney U test. (G) Mean cell speed <v> versus ∣tr(M)∣ for different treatments and matrix densities. (H) Persistence length of cell trajectories, Sp, versus ∣tr(M)∣ for different treatments and matrix densities. (I) Chemotaxis index CI versus ∣tr(M)∣ for different treatments and matrix densities. *P < 0.05 and **P < 0.01 according to the Mann-Whitney U test.

However, the treated cells were still able to appreciably deform their microenvironment. Given that myosin-II inhibition caused a marked decrease in Sp for intermediate [col] (i.e., 0.75 mg/ml; see Fig. 7D), we measured the matrix deformation moments created by blebbistatin-treated cells in these matrices (Fig. 10E) and found that ∣tr(M)∣ was significantly reduced after treatment. The matrix deformation moments in our most restrictive matrices were about 10 times lower than in nonrestrictive matrices and also decreased after ck666 or blebbistatin treatment (Fig. 10F). Notably, these treatments also caused the distribution of the signs of sign(tr(M)) to be symmetric, implying that they reduced matrix contraction toward the cell (fig. S1).

Since blebbistatin- and ck666-treated cells had lower migration speeds than untreated cells, we examined the relationship between ∣tr(M)∣ and 〈v〉 (Fig. 10G). This analysis revealed that ∣tr(M)∣ and 〈v〉 collapsed on a single curve regardless of treatment and [col], suggesting that the ability of 3D migrating neutrophils to mechanically remodel their environment is a key factor governing their speed. Given that cell speed decreased when cells were turning (Figs. 4, G and H, 6B, and 8B), we plotted Sp versus ∣tr(M)∣ (Fig. 10H) and found that, for both treatments and control conditions, Sp decreased when ∣tr(M)∣ fell below ~5 × 103 μm2 regardless of [col]. Thus, neutrophils may transition to a migratory strategy with an increased rate of turning when they become unable to sufficiently deform their surrounding matrix. Last, we plotted the CI versus ∣tr(M)∣ (Fig. 10I) and found that this index did not depend on ∣tr(M)∣ except when Arp2/3 was inhibited.


The present study investigated the ability of neutrophils to migrate up a chemoattractant gradient in 3D environments with different degrees of complexity. We hypothesized that turning events would be a key factor for migration toward a chemoattractant source by allowing cells to circumvent obstacles in the matrix (e.g., narrow pores) and quickly resume their motion toward the chemoattractant. We tested these hypotheses using a custom-made 3D collagen matrix chemotaxis chamber, label-free automated cell tracking, and statistical analyses to quantify the migratory trajectories of neutrophil-like HL-60 cells and the matrix deformation these cells exert during migration.

Previous studies used time-dependent kinematic analyses of leukocyte trajectories in 3D matrices to find that cell speed and persistence decrease with decreasing matrix pore size (11, 20, 22). Although more scarce, there are previous reports of 3D matrix deformations caused by migrating leukocytes (42). In the present study, we analyzed cell trajectories based on distance traveled rather than time elapsed during migration. This interpretation allowed us to study the geometry of cell trajectories and decouple it from the speed of cell migration, which varies along the trajectories as the cells navigate different matrix sections. Furthermore, our label-free automated cell tracking made it possible to analyze more than 20,000 migrating cells, yielding highly detailed, multiparametric statistics of the migratory process. Together with our 3D measurements of matrix deformation, these data have revealed mechanistic interplays between cell contractility and migratory parameters that had not been quantified before.

We first characterized how increasing matrix collagen concentration ([col]) caused a decrease in matrix pore (ξ) size and heterogeneity, leading to a stark reduction in available pores larger than the cell nucleus in our most restrictive matrices. We found that the decreasing ξ slowed cells down but did not diminish their ability to migrate directionally up the chemoattractant gradient, consistent with previous studies involving other leukocyte subtypes (20, 22). Our measurements of 3D matrix deformations indicated that neutrophils could transiently deform low-[col] matrices using physical forces. In particular, the cell-generated collagen deformations were oriented preferentially toward the cell and along the chemoattractant gradient direction. However, the cells became unable to appreciably deform high-[col] matrices and the small deformations they generated in the collagen were more isotropic. The cells adapted to increased matrix rigidity in high-[col] matrices by turning more often, which caused a marked decrease in the spatial persistence of the direction of motion. This transition was triggered by the cells’ inability to physically remodel their microenvironment, which was confirmed by pharmacological inhibition of contractile forces mediated by myosin-II, or protrusive forces exerted by Arp2/3-mediated dendritic actin polymerization.

The physical properties of different types of tissue, in which the ECM is a significant component, can vary greatly across different locations in the body (43). However, during immune responses, this variation does not impede neutrophils from migrating toward sites that initiate inflammatory responses (8). The extracellular matrices in our study contained average pore sizes well within physiological ranges (44). Therefore, the versatility in migration strategies that we found in vitro has implications for understanding the ability of leukocytes to migrate across a range of environments during innate immune responses (9). This adaptability is crucial for effective chemotaxis because, unlike cancer cells, neutrophils do not heavily rely on proteolytic matrix degradation when facing restrictive 3D environments (20, 45). Cancer cells can transition to proteolytic matrix degradation when presented with rigid microenvironments that they cannot deform by means of physical forces (46). The ensuing permanent matrix remodeling has been shown to cause migratory statistics consistent with an anisotropic random walk, thus ensuring a long-range directional bias (23). In contrast, 3D migrating neutrophils adapt to noncompliant microenvironments by turning more often; thus, long-range directional bias is achieved by a different means.

In our experiments, the short-range statistics of 3D neutrophil migration became progressively disturbed with decreasing ξ, becoming similar to those of randomly moving cells. However, neutrophils conserved their chemotactic efficiency, measured as the ratio of distance traveled toward the chemoattractant to total distance traveled, even in the most restrictive matrices. This result hints at a trade-off between short-range and long-range goals in the directed 3D migration of neutrophils. A key factor in this trade-off is the ability to correct the misalignments with respect to the gradient that arises when circumventing impassable sections of the matrix. Our data show that, in contrast to 2D chemotaxis, neutrophils migrating in 3D matrices need a distance significantly longer than their persistence length to realign their motion with the gradient. While cell migration is a continuous process in which the cell continuously changes direction and it would be oversimplistic to describe it as a sequence of discrete events, these data imply that neutrophil directed motion through 3D matrices involves the coordination of consecutive turning events. We found that a biased random walk model (41) with no adjustable parameters reproduced the evolution of cell motion orientation along the cell trajectories (fig. S6). In the biased random walk model, cells steer toward the chemoattractant source over a bias distance, Sθ,which dictates the chemotactic efficiency in competition with the persistence of cell migration, Sp; i.e., cells with higher values of Sp/Sθ have longer portions of their trajectories aligned with the chemotactic gradient. Our experimental data suggest that 3D migrating neutrophils maintain their chemotactic efficiency as ξ is decreased by varying their bias distance in concert with their persistence so that Sθ ≈ 3Sp. Thus, neutrophils seem to have developed a 3D migration strategy that uses obstacles in the matrix as pivoting points to progressively reorient cell motion toward the chemoattractant source. In this strategy, denser matrices present not only more obstacles but also more pivoting points per unit volume, allowing cells to maintain efficient chemotaxis for a wide range of collagen densities.

We tested our hypothesis about turning events and chemotactic efficiency by inhibiting Arp2/3, a complex shown to play a role in pathfinding during 3D neutrophil migration (11). Neutrophils migrating through straight microchannels have been reported to increase their speed and persistence after Arp2/3 inhibition (47). However, the microchannels in those experiments did not challenge the cells to navigate physical barriers in search of migratory paths, which is different from our system and in many physiologically relevant contexts (44). The neutrophils in our study had lower speeds after Arp2/3 inhibition due to their encounters with narrow matrix pores and the difficulties of either searching for alternate paths or exerting traction forces to remodel the surrounding matrix. Consistent with these results, previous work has shown that the Arp2/3 complex plays an important role in traction force generation for cancer cells (48). Furthermore, our statistical analysis of trajectory orientation suggests a crucial role for the Arp2/3 complex in coordinating consecutive turns to circumvent narrow matrix pores, while keeping long-range directional motion toward the chemoattractant source. Specifically, Arp2/3 inhibited cells experienced a marked decrease in Sp/Sθ, implying that they needed to concatenate a large number of turning events to correct misalignments in their motion. The effects of Arp2/3 inhibition were more prominent in high-[col] matrices, where neutrophils were more likely to encounter narrow, rigid pores. These results are consistent with previous studies involving dendritic cells and neutrophils, in which lamellipodia formation appeared to be crucial for navigational decision-making through complex geometries, especially as the microenvironment becomes more restrictive (16).

Recent studies with dendritic cells have shown that myosin-II becomes increasingly important during 3D neutrophil migration as matrix pore sizes became smaller (15, 49). In particular, the transition of the nucleus through matrix pores has been found to be a rate limiting factor in 3D migration (49). In our experiments, even after myosin-II inhibition, neutrophils revealed their remarkable migratory versatility by increasing the rate of turning to circumvent matrix pores they would have otherwise been able to physically expand or squeeze through via contractile forces. While this adaptation translated into a decrease in cell speed, it did not seem to interfere with the ability of neutrophils to perform chemotaxis.

In conclusion, we analyzed the interplay between matrix deformation and turning events in 3D neutrophil migration through matrices of varying pore sizes. We found that, as [col] increases, neutrophils transition from a migratory regime characterized by transient large matrix deformations and straight trajectories to a second regime characterized by low matrix deformation and frequent turning events. More frequent turning translated into slower cell speed at higher collagen concentrations; however, it did not cause a decrease in chemotactic efficiency because the cells coordinated their more frequent turns to quickly resume directed motion. By pharmacologically inhibiting myosin-II contractility and Arp2/3 mediated pseudo-pod protrusion, we were able to trigger the transition to turning based migration at lower collagen concentrations than in untreated cells, confirming that this transition takes place when cells become unable to significantly deform the matrix. In addition, Arp2/3 inhibition caused a decrease in the frequency of turning events, leading to more persistent trajectories. However, it also caused a loss of coherence among turning events that made the Arp2/3-inhibited cells less efficient at chemotaxis. Overall, this work contributes to an enhanced mechanistic understanding of the role that the Arp 2/3 complex, myosin-II, and surrounding collagen density play in 3D neutrophil chemotaxis, which, in turn, may provide crucial insights into molecules and processes to target for regulating pathological inflammatory responses (50).


Cell culture, differentiation, and cytoplasmic labeling

We cultured and differentiated cells from the human promyelocytic leukemia cell line (HL-60, ATCC) into neutrophil-like cells (dHL-60) as previously described (51). We grew HL-60 cells to approximately 1 × 106 cells/ml and passaged the cells every 2 to 3 days in RPMI 1640 (Thermo Fisher Scientific) supplemented with l-glutamine and 10% fetal bovine serum (FBS; Omega Scientific). We differentiated HL-60 cells by taking 1.5 × 106 cells from each passage and placing them in supplemented RPMI 1640 with the additional 1.3% DMSO (Sigma-Aldrich). Next, we incubated both HL-60 and dHL-60 cells at 37°C and 5% CO2. We then performed experiments with dHL-60 cells 4 days after being cultured in the presence of DMSO. For experiments involving fluorescently labeled cells, we labeled day 4 differentiated dHL-60 cells with CellTracker Deep Red (Thermo Fisher Scientific) before mixing them with collagen solutions. Cells were incubated in a 10 nM solution of CellTracker in RPMI 1640 without FBS for 30 min at 37°C and 5% CO2. Afterward, the cells were centrifuged and resuspended in RPMI 1640 containing 10% FBS.

In vitro 3D directed migration assay

We adapted a version of a previously published custom-built chemotaxis device to study 3D neutrophil chemotaxis (Fig. 1A) (25, 52). First, we treated 25-mm glass coverslips (Thermo Fisher Scientific) on one side with 250 μl of 0.1 M NaOH for 5 min. We then removed the NaOH and rinsed the coverslips with distilled H2O. We let the coverslips dry and added (3-aminopropyl)triethoxysilane (Sigma-Aldrich) to their treated sides afterward. After 30 min, we rinsed the coverslips again with distilled H2O and placed them on Kimwipes for air-drying, with the treated surfaces face up. Next, we punched 12-mm-diameter holes in the center of 50 × 9-mm Falcon petri dishes (BD Falcon). We then attached treated glass coverslips to the bottoms of the petri dishes with their treated sides face up using vacuum grease (Beckman Vacuum Grease Silicone). We attached separate treated glass coverslips to the tops of the holes with their treated sides face down and covering the majority of the holes’ diameters. This sandwiching of coverslips created a pocket that was used later to place collagen gel solutions. Before fabricating collagen gels, we rested the devices on top of a metal block that sat in a container filled with ice.

Collagen gel fabrication and fluorescent labeling

We fabricated collagen gels by neutralizing solutions of rat-tail type-I collagen dissolved in a 0.2-N acetic acid (Corning) with 1-N NaOH in proportions recommended by the manufacturer. We added a suspension of supplemented RPMI 1640 containing 4.7 × 103 dHL-60 cells instead of the recommended 10× phosphate-buffered saline and water. We also mixed a 1.2% solution of fluorescent microspheres (Molecular Probes) with the collagen solution. While mixing the different constituents of the collagen gel, we placed all solutions in a small container filled with ice. After the collagen gel solution was made, we pipetted 113 μl of the solution into the pocket of each device and placed each device in an incubator at 37°C and 5% CO2 for 30 min. We placed the devices in an incubator vertically to allow the gel solution to settle to the lower half of the pocket. We then added 113 μl of supplemented RPMI 1640 into the remnant volume of the pocket and placed the device in the incubator for another 30 min. Last, before imaging, we replaced the supplemented RPMI 1640 with either fresh supplemented RPMI 1640 or a 50 nM solution of fMLP (Sigma-Aldrich).

For experiments involving reflection microscopy imaging, we fabricated collagen gels in the same manner as described above, leaving out the cells and fluorescent microspheres in the supplemented RPMI 1640. We placed 100 μl of the gel solution on top of glass coverslips attached to the bottoms of petri dishes with 12-mm-diameter holes. We did not add a second coverslip on top of the holes. We then put each petri dish in an incubator at 37°C and 5% CO2 for 30 min. Afterward, we filled these dishes with 2 ml of supplemented RPMI 1640.

We fluorescently labeled collagen fibers with 5-carboxytetramethylrhodamine (TAMRA; Thermo Fisher Scientific) by following a previously published protocol (53). We first added 25 mg of TAMRA to a 2.5-ml solution of DMSO and stored at −20°C in an Eppendorf tube. We covered the Eppendorf tube with aluminum foil to protect the solution from light. Next, we made a labeling buffer of 0.25 M NaHCO3 and 0.4 M NaCl and subsequently adjusted the pH to 9.5 using a 10 M solution of NaOH. We then filled a 1-ml syringe with a solution (8.46 mg/ml) of rat-tail type I collagen in a 0.2-N acetic acid solution (Corning) and injected the solution into a 3-ml dialysis cassette. The cassette we used had a 10-kDa molecular weight cutoff. Before retracting the hypodermic needle, we removed air from the cassette by pulling the syringe plunger up. We then performed an overnight dialysis in a beaker containing a 1-liter solution of our labeling buffer at 4°C. One hundred microliters of our TAMRA solution (10 mg/ml) was mixed with 900 μl of our labeling buffer, with both solutions at room temperature. After mixing, we placed the TAMRA solution in a 4°C fridge. The collagen was removed from the cassette using a 2-ml syringe, and 1 ml of the solution was mixed with 1 ml of diluted TAMRA using a pipette. We then placed the resulting solution in a microcentrifuge tube and incubated with rotation overnight at 4°C. At this point, the collagen was labeled with the TAMRA dye and we put this mix in a dialysis cassette and dialyzed it overnight at 4°C against a 1-liter solution of our labeling buffer. Last, we dialyzed the TAMRA-labeled collagen overnight at 4°C against a 0.2% acetic acid solution with a pH of 4. We used the final volume to calculate the concentration of the TAMRA-labeled collagen.

Inhibition drug treatment experiments

We performed Arp2/3 complex and myosin-II inhibition experiments using ck666 (Sigma-Aldrich) and blebbistatin (Abcam). For each treatment condition, we incubated day 4 differentiated cells with either a 100 μM DMSO-diluted solution of ck666 in supplanted RMPI 1640, a 10 μM DMSO-diluted solution of blebbistatin in supplemented RPMI 1640, or DMSO vehicle control solution in supplemented RPMI 1640 at 37°C and 5% CO2 for 30 min. We then centrifuged the treated cells and resuspended the cells in supplemented RPMI 1640. Afterward, we followed the aforementioned protocol for embedding cells in collagen gels, with the addition of ck666, blebbistatin, or both so that the final concentration of each gel in these experiments contained 100 μM ck666 and 10 μM blebbistatin. We added supplanted medium to the chemotaxis chambers after gel formation. Our chemoattractant solutions containing fMLP also contained the drugs at the same concentrations as those used when preparing the collagen gels.


We imaged collagen fibers for pore size analysis using a Leica SP5 microscope in reflection mode with a 40× immersion lens at 2× zoom. We obtained bright-field images for population migration experiments using an enclosed Leica DMI6000 B microscope with a 5× air lens at 37°C and 5% CO2. We retracted the microscope’s magnifier to produce images that made automated cell tracking easier. In each experiment, we acquired four planes at 200-μm spacing in each collagen gel every ∆t = 30 s for 1 hour. We imaged four planes to follow a larger a number of cells and we chose 200-μm spacing to avoid imaging the same cells in more than one plane.

We imaged single cells and matrix deformations on an enclosed Zeiss LSM 880 laser scanning microscope in the fast Airyscan mode, with a 40× water lens and 1.3× zoom. We acquired two 80-μm imaging stacks with 0.7-μm spacing between planes every 30 s. The duration of time-lapse experiments varied. We imaged cells in each experiment until cells of interest left the field of view. Given that dHL-60 cells are relatively fast-moving cells, we switched imaging channels for the labeled cells and microspheres per line scan during acquisition. This was to minimize significant differences in cell positions with respect to matrix deformations during the acquisition of the two stacks for a single time point. The raw Carl Zeiss Image Airyscan output files from the Zeiss LSM 880 microscope underwent Airyscan processing using Zeiss’s imaging software Zen Black. After Airyscan preprocessing, these files were loaded into FIJI (54) and then exported as uncompressed .TIFF files. The .TIFF files were loaded into MATLAB (55) for further quantitative analysis.

Pore size analysis

We used a bubble analysis method to characterize the distribution of pore sizes from confocal reflection images of collagen networks, as previously described (56). Before performing this analysis, we preprocessed the images in FIJI. First, we removed a radial intensity gradient that appeared as an artifact from our imaging system using FIJI’s sliding paraboloid feature with a radius of 50 pixels under the software’s background subtraction option. We then removed isolated pixels in the images using FIJI’s despeckle feature, which applied a median filter to the images. Last, we exported the preprocessed images to MATLAB to perform the bubble analysis.

In our bubble analysis, we binarized each image by setting an intensity threshold equal to the mean image intensity, which resulted in the segmentation of the collagen fibers. We then removed segmented objects containing less than 400 pixels to eliminate noise. Next, we computed the distance to the nearest segmented fiber for all points inside the pores. Afterward, the resulting distance maps were smoothed out with a Gaussian window of size 100 × 100 pixels and 1/e-radius of 50 pixels. We identified pore centers using the local maxima in the distance map, with maximum distance values representing pore sizes. We characterized the pore size distributions in terms of their probability density function (p.d.f.) p(ξ), mean pore size 〈ξ〉, SD, and the probability of finding pores with diameters larger than 5 μm, i.e., p(ξ ≥ 5 μm).

Automated cell tracking

We implemented a custom MATLAB algorithm to automatically track cells in time-lapse sequences of bright-field images, B(x, y, t) (Fig. 1B). First, we processed the initial image of each sequence, B(x, y, 0), in two steps to identify the positions of in-focus cells.

Step 1. Noise removal and segmentation. Since retracting the microscope magnifier produced images where in-focus cells resembled a negative Laplacian of a Gaussian, i.e.h(x,y)=2σ2(xx0)2(yy0)2σ4e((xx0)2+(yy0)22σ2)(5)where (x0, y0) are the coordinates of the cell center and σ is a free parameter related to cell size, we convolved each image with a centered 6 × 6 pixel h(x, y) window. This operation, i.e., B¯(x,y,0)=B(x,y,0)*h(x,y) resulted in noise removal. Next, we computed a histogram of image intensities and used intensities higher than the distribution’s first peak, which corresponded to background, as a global threshold to segment cells in B¯(x,y,0).

Step 2. Calculation of cell positions. We centered a 10 × 10 pixel search window at the centroid of each cell segmented in step 1, cropped, and normalized the image by its maximum intensity. We then fitted h(x, y) to the resulting image using least squares with x0, y0 and σ as fit parameters. This operation yielded the coordinates of the center of each cell, x(0)=(x0,y0).

For subsequent time points, we followed the same two-step procedure described above with a slight modification that allowed us to track cell trajectories: In step 2, the search window for each cell in frame t was centered at the nearest cell centroid determined in the previous frame, i.e., x(tt).

Spurious track removal. We encountered two situations leading to spurious cell tracks. Cells that moved outside of the field of view were lost causing the algorithm to track random points. This resulted in object displacements much larger than the average cell displacements, which were <20 pixels. Thus, objects that had displacements greater than this threshold between any two successive frames were automatically excluded. Cells that moved too close to one another sometime led algorithm to switch their track labels. Given that these events were not frequent and that, in some cases, it was not straightforward to determine the correct labeling, all cell track pairs that intersected each other at the same time instant were automatically excluded.

Stage drift correction. We noted that each experiment contained either dead cells or cells that were nonmotile and tracked the motion of these cells as a surrogate of stage drift. Thus, we computed the net displacements of all tracked cells, d(T)=x(T)x(0), where T = 60 min. We then selected the cells in the bottom 15% of d(T) as nonmotile cells. Their trajectories were fitted to the model xdrift(t)=x(0)+vt, where v represents the drift velocity, and ensemble averaged over all the nonmotile cells, xdrift(t). Since drift was systematically largest within the first 20 min after adding buffer or fMLP, we did not include data from that initial 20-min period in our analysis. Drift was removed with x(t)xdrift(t), where x(t) was the trajectory for a cell in the upper 85% of d(T) the bracket operator that denotes ensemble averaging for all the motile cells.

Analysis of cell trajectories in intrinsic coordinates

We analyzed the geometry of each cell’s trajectory in intrinsic coordinates (i.e., distance traveled and angular orientation). To perform this analysis, we first tracked cell positions versus time as described above. We then computed the distance traveled by the cell versus time, i.e., s(t)=τ=0τ=ts(τ), where s(t)= s(t+t)s(t)=x(t+t)x(t), and used this information to obtain x(s) via interpolation. The vector tangent to the trajectory was obtained as t(s)=dx(s)/ds. The angular orientation of cell motion with respect to the chemoattractant gradient was quantified by the angle θ=cos1[t(0,1)]. The average cell speed was defined from the length s of the cell trajectory over T = 60 min, i.e., v=T10s(T)ds. A cell was considered to be motile if its net displacement d(T) was larger than 30 μm. We computed the chemotaxis index of the trajectories using CI=cos1[d(T)(0,1)d(T)]. The MSDs of the cells were computed as a function of distance lags, Δs, i.e., r(s)=x(s+s)x(s)2, where the bracket operator denotes ensemble averaging for all the motile cells in each experiment and along all values of s along each trajectory. Similarly, we computed the autocorrelation function of the vector tangent to the trajectory, Ctt=t(s)t(s+Δs). We estimated the persistence distance Sp of the trajectories for each experimental condition by assuming that this correlation decays exponentially with distance traveled, i.e., CtteΔsSp. Specifically, we determined the value Δs0 for which the correlation function reached the cutoff of Ctt,0 = 0.05 and estimated the persistence length as Sp ≈ −Δs0/ log Ctt,0. To better understand the changes in direction along the cells’ trajectories, we calculated the changes in orientation of t versus Δs, i.e., θ(Δs)=cos1[t(s)t(s+Δs)], and determined their p.d.f. conditioned to Δs, p(θ∣Δs). We also determined the joint p.d.f. of speed and changes in orientation, p(θ, v∣Δs), as well as the p.d.f. of ∆θ(Δs) conditioned to θ and Δs,p(∆θ∣θ, Δs).

3D matrix incremental deformation fields and matrix deformation moments

We computed the incremental changes in 3D deformation fields inside our collagen gels, u=(u,v,w), using a multigrid PIV approach (57) that was an amended version of a previous in-house 3D PIV algorithm (58). Briefly, we performed a first PIV pass on each confocal stack of fluorescently labeled beads using the stack from the immediately preceding time point as a reference. An interrogation volume of 56 × 56 × 24 pixels with a spacing of 28, 28, and 22 pixels, respectively, was chosen during this pass. We then ran PIV again using interrogation volumes of 28 × 28 × 22 pixels with a spacing of 14, 14, and 11 pixels, respectively. During the second pass, the smaller interrogation windows of the reference stack were shifted by the displacements calculated in the first pass. This resulted in a 3D deformation field that captured large deformations with the relatively large interrogation volumes in the first pass, while also achieving higher spatial resolution from the smaller interrogation volumes used during the second pass (59). We constructed visualizations of the 3D cell shapes and incremental deformation fields by creating Visualization Toolkit files, which we imported into the open source software ParaView (60) for 3D visualization.

To quantify a cell’s ability to deform the matrix, we expanded Butler et al.’s (31) concept of contractile moment to three dimensions and determined the first-order moments of u. Using the 3D cell shape data from confocal stacks of fluorescently labeled cells, we calculated the cell centroids xC and defined a region around each cell, R:(10 μm<xxC <200 μm), to calculate the matrix deformation first-order moment tensor M.

The diagonal elements of M provided the total contraction or extension of the collagen gel around individual cells, in each ordinate direction between two consecutive time points (Δt = 0.5 min ). The net 3D matrix deformation moment displacing the collagen gel was given by tr(M). Note that tr(M) > 0 indicates that the matrix is being deformed away from the cell, whereas tr(M) < 0 indicates that the matrix is being deformed toward the cell. The matrix deformation moment tensor was projected along the directions parallel and perpendicular to the chemoattractant gradient, and the mean absolute values of these projections, ∣M∣ and ∣M∣, were quantified. Last, the time-averaged absolute quantity 〈∣tr(M)∣〉 represented a concise metric for gauging cells’ ability to incrementally deform surrounding ECM fibers.

Statistical analysis

For all bar graphs, statistically significant differences between [col] were tested using a Mann-Whitney U test with Bonferroni corrections to account for multiple comparison testing.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We are grateful to R. A. Firtel (Section of Cell and Developmental Biology) at University of California, San Diego and K. Ley (La Jolla Institute for Immunology) for insightful comments. Funding: J.F. was supported by the Ford Foundation Program through the Ford Foundation Dissertation Fellowship, a National Heart, Lung, and Blood Institute’s T32 Fellowship (HL105373), R01 GM084227, and the University of California San Diego’s School of Medicine Microscopy Core under a National Institute of Neurological Disorders and Stroke P30 core grant (NS047101). Author contributions: J.F., Y.-T.Y., R.M., S.C., J.C.L., and J.C.d.Á. designed the research. J.F., A.S., and C.A. performed the experiments. J.F., A.K., C.A., and J.C.d.Á. analyzed the data. J.F., Y.-T.Y., A.K., S.C., J.C.L., and J.C.d.Á. wrote the paper. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplemental Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article