A microfluidic platform to investigate the role of mechanical constraints on tissue reorganization

During embryo morphogenesis or wound healing, biological tissues are sculpted by active processes, such as cell growth and division, cell migration or cell differentiation (Heisenberg and Bellaïche, 2013). These processes generate mechanical stresses within the tissue that can trigger tissue-scale flows, cell deformation and cell rearrangements. The nature of the response of the tissue to mechanical stresses depends on its rheological properties at the scale of a ‘cell group' or ‘representative volume element’ or ‘cohort of cells’ (Khalilgharibi et al., 2016). In turn, these cell-group-scale properties arise from the cell-scale interplay between cell viscoelasticity and cell-cell interactions, such as adhesion.

To understand fully and measure these properties, 3D in vitro tissues, such as cancer cell spheroids or embryonic organoids (Takebe and Wells, 2019), play an increasingly important role in biophysical and medical studies as they are easy to handle and their properties are relevant to further understanding in vivo tissues.

To study this in vitro tissue rheology in the context of developmental biology, several experimental setups have already been used, including cell monolayers stretching (Harris et al., 2012), rounding of cell aggregates (Mombach et al., 2005; Marmottant et al., 2009; Grosser et al., 2021), shearing of cell monolayers or cell aggregates (Fernández et al., 2007; Stirbat, 2013; Guevorkian et al., 2010; Yousafzai et al., 2020 preprint), and compression of cell aggregates (Marmottant et al., 2009; Montel et al., 2012; Mary et al., 2022). These experiments have enabled bulk parameters, such as shear and compression moduli, viscosity, and yield deformation or yield stress, to be measured. To connect these mechanical properties to cellular processes, it is essential to have access at the same time to the cell structure and to its dynamics (e.g. cell shapes or rearrangements) (Hopkins et al., 2022).

Live imaging in vivo (during embryo development or insect metamorphosis, for instance) or in vitro provides the opportunity to track cell movements and divisions, cell shape variations and cell rearrangements, and to quantify their contributions to the cell-group or tissue-scale dynamics, by measuring their contributions to the cell-group or tissue deformation rate, for example (Guirao et al., 2015; Etournay et al., 2015).

Methods combining the application of an external stress, force, deformation or deformation rate controlled by the experimentalist with live imaging at subcellular resolution have proven to be powerful approaches to investigate biological tissues and, more generally, the rheological properties of complex materials. For example, cell monolayers can be stretched by being suspended between rods or on an elastic membrane while imaging the tissue response at the cellular level (Wyatt et al., 2015; Fouchard et al., 2020). Nanometric magnetic beads naturally endocytosed by cells have also been used to apply a global force on cell aggregates (Du et al., 2017; Mary et al., 2022).

A fruitful approach for 3D tissues is to embed ferro-magnetic droplets within the tissue, which can be actively deformed by an externally applied magnetic field while imaging the tissue response at the cell scale. This method has enabled the measurement of spatial gradients of yield stress leading to a liquid-solid transition in the zebrafish mesoderm, and also gives access to the visco-elastic response of cells at short time scales (Mongera et al., 2018). This approach provides a powerful method to probe local tissue rheology for small deformations. To measure effective properties for large deformations, pipette aspiration is used broadly on 3D tissues, such as cancer spheroids (Guevorkian et al., 2010), embryos (e.g. mouse, Xenopus; David et al., 2014), chicken (Regev et al., 2022) and explants. However, it is not easy to image the tissue microstructure and its dynamic response in such aspiration geometries. Cell-cell rearrangements are of high importance in developmental biology and are rarely observed in experiments in 3D. There was a need to design a system whereby they could be properly quantified in order to understand which biophysical principles govern their dynamics.

Here, we present a method to study cell-group scale mechanical properties in 3D tissues and to relate these properties to cell-scale characteristics, such as cell shape and heterogeneities in cell protein and gene expression. This approach is inspired from studies on non-living amorphous materials (Cheddadi et al., 2011; Golovkova et al., 2020). Our method combines microfluidic aspirations with two-photon imaging of cell shapes at subcellular resolution. It enables the spatiotemporal patterns of deformations generated within the 3D tissue, the visco-elasto-plastic properties of cell groups, and the dynamics of cell rearrangements to be determined and correlated.

Aspiration geometry, i.e. forcing a 3D cell aggregate to flow through a constriction, enables a heterogeneous flow to be applied. As a consequence, each individual experiment provides a self-sufficient data set with a large enough variation range of cell-group deformations and deformation rates. Moreover, this geometry enables the induction of large deformations and the modification of tissue boundary conditions to probe tissue elasto-capillary behavior.

We pre-confine the aggregate and use a capillary with straight sides (rather than a cylindrical pipette). The tissue is physically constrained in the vertical dimension corresponding to the optical axis. This forces it to flow within the plane of observation (Fig. 1, Movie 1). Such quasi-2D flow facilitates imaging of cell shapes and rearrangements.

In principle, the flow could have a 3D structure, namely depend on the direction z, perpendicular to the device plane. This would be the case, for instance, for an aggregate with a core-shell structure, or any heterogeneous tissue. In that case, our set-up would be suitable to perform 3D image analysis, provided the flow was slow enough. We present experiments performed on homogeneous 3D aggregates of cells that do not have core-shell structures (Fig. S1, Movie 2), and we confirm that we can neglect the flow variation along the z direction (Movie 3).

We first focused on cell-group scale aspects of the flow properties. The use of brightfield or phase-contrast microscopy provides valuable information and enables the velocity field and deformation rate to be mapped in the aspired aggregate. The velocity profile inside the channel was consistent with a plug-flow (Fig. 1D, inset) and we integrated it over the channel cross-section to measure the aggregate flow rate. Over time, we observed three distinct phases of the aspiration (Fig. 1E): (1) a rapid elastic aspiration phase that corresponds to the rapid entrance through elastic cell deformation of the aggregate into the channel; (2) a stationary phase with a fairly constant flow rate, corresponding to an effective viscous behavior of the tissue; and (3) a final phase of acceleration when the non-aspired aggregate part size is similar to the channel size, which indicates that the dissipation on the channel walls is negligible compared with the viscous dissipation inside the aggregate.

Taking advantage of the pseudo-2D flow structure, we quantified the in-plane spatial variation of the flow and deformation rate within the aggregate (Fig. 2A,B). Models (Guevorkian et al., 2010) predict that in a pipette experiment, dissipation occurs only in a typical volume of effective radius R_channel. In our setup, we were able to obtain a spatially resolved measurement of the deformation rate: we measured the velocity gradient, symmetrized it to determine the cell-group deformation rate, and plotted its anisotropic part. In our geometry, the typical size of the dissipation zone was similar to R_channel, but we observed an extended zone along the walls of the chamber. The isotropic part of the cell group deformation rate (Fig. 2C) revealed compression against the walls before the channel entrance and expansion after the entrance of the channel, which illustrates tissue compressibility.

To investigate the link between the material properties and the cell dynamics, we imaged the tissue at the cell scale to measure cell deformation. Here, we used the same formalism that we used previously on a 2D system (MDCK monolayers; Tlili et al., 2020). Cell divisions can be neglected when they are inhibited, or when the timescale of the experiment T_exp is much smaller than the timescale of cell cycle T_div; this is the case here, where T_exp ≈ 1 h and T_div ≈ 12 h. The cell group (total) deformation rate is the sum of the cell deformation rate coarse-grained at the cell-group scale and plasticity deformation rate as a result of the rearrangement rate. The cell deformation pattern is determined by the balance between: (1) the deformation induced by the constriction, i.e. the source term ⁠; (2) the relaxation of cell deformation due to cell rearrangements (⁠⁠); and (3) the elasticity transport terms (i.e. the fact that non-deformed cells are constantly advected towards the channel). We wanted to estimate the effective viscoelastic relaxation time related to cell rearrangements τ_r (Fig. 3), which can be defined as ⁠, where the superscript ‘dev’ indicates the deviatoric (i.e. anisotropic) part of a tensor.

In Fig. 4C, we compared the experimental measurement of deformation evolution, 〈ε_cell〉_box (red curve), with various simulations: we varied τ_r between 10 min (light gray) and infinity (black), the latter being a purely elastic limit (Fig. S3). The simulated deformation amplitude agreed reasonably with the experiment when we used τ_r ≈ 20 min. The Fourier determination of deformation was sufficiently discriminant and reproducible to determine correctly the amplitude order of magnitude and we could reasonably exclude values of τ_r outside of the interval [15 min, 25 min] (Fig. S4).

Thus, measurements coarse-grained at the cell-group scale yielded an estimate of an effective cell deformation relaxation time. Because divisions play a negligible role on the experiment timescale, the effective cell deformation relaxation time is necessarily associated with the dynamics of cell rearrangements, which were numerous in particular near the constriction entrance (Fig. S4D, Fig. S5, Movie 5).

Note that the determination of τ_r in Fig. 4C is not sensitive to the shape of the deformation evolution. The peaks in deformation are not systematically observed and their origin is unknown; they are likely due to a stick-slip friction on the glass coverslip (rather than to cascades of rearrangements, which we did not detect). Determining accurately the shape of the deformation evolution is beyond the scope of the present study. It would require the time resolution and signal-to-noise ratio to be improved simultaneously beyond the possibilities of the current Fourier method with our current image quality. This would be possible in principle using cell contour segmentation complemented with detailed tensorial analysis.

Using the same data, but at cell scale, we have access to the cell rearrangements dynamics at the scale of a single event (Fig. 5A) by manually identifying and segmenting cell-cell contacts that are disappearing and newly created.

As the deformation rate and cell deformation fields are heterogeneous, we can trigger in the same experiment: (1) ‘quasi-static’ rearrangements in deformed regions of the tissue away from the channel, where both a disappearing link between cells losing contact and a newly appearing link between cells gaining a new contact can be clearly identified; and (2) ‘sliding’ rearrangements at the vicinity of the channel, where cells are highly stretched, flow rapidly and seem to slide across each other. In the first case, where rearrangements can be clearly identified and well-defined, we manually measured the disappearing and appearing junction lengths versus time (Fig. 5B). We defined velocities of appearance and disappearance before and after the rearrangement by taking the tangent of the junction length curves before and after the four-fold vertex stage. To establish a robust definition, we defined the velocity of disappearance V⁻ = L⁻ (t = − T)/T and of creation V⁺ = L⁺ (t = + T)/T at T = 6 min. Their values were unchanged at T = 10 min, as the length evolution is locally linear around the four-fold vertex stage. We often observed an asymmetry in V⁺ and V⁻ (Fig. 5C).

To relate these quantities with the local cell deformation during the rearrangement, we projected the tensor ε_cell in the direction of the junction that disappeared/appeared, yielding where θ⁻ is the angle between the disappearing link and the principal direction (i.e. the eigenvector corresponding to the positive eigenvalue) of ε_cell; and where θ⁺ is the angle between the appearing link and the principal direction of ε_cell (Fig. 5D). We observed that the disappearing junction is compressed, and that the newly created junction is elongated, which is consistent with the fact that observed rearrangements are triggered by the cell-group-scale flow imposed by the aspiration.

We plotted V⁺ and V⁻ versus and ⁠, which were averaged at the same time, T = 6 min, as the velocities in Fig. 5D. For − 0.3 < ε_cell < 0.3, the velocities V⁻ and V⁺ were around the mean value ≈ 0.5 μm/min. For ⁠, the mean value of V⁺ increased suddenly up to ≈ 2.5 μm/min (see also Movie 5). Examples of such highly deformed rearranging cells are shown in Figs S4D and 5A; these occurred mainly near the constriction entrance. They were usually compatible with the tissue confluence and integrity. High aspiration pressures and velocity sometimes resulted in aggregate fractures (Movie 6): we did not analyze these experiments.

We also looked at the lifetime T_vertex of a ‘four-fold vertex stage’ (Fig. 5E). We defined it as the stage when the cell-cell junction was below 1.5 μm, which is small compared with the cell size, and large enough to avoid artefacts due to the resolution of the length determination when the junction is small. We observed that the lifetime of a four-fold vertex at null deformation rate can be large, which is qualitatively different from the physics of liquid foams, in which the four-fold vertex is unstable (Cantat et al., 2013). The lifetime of the four-fold vertex in a biological tissue can be related to the dynamics of adhesion remodeling (Bardet et al., 2013; Guirao and Bellaïche, 2017), i.e. the timescale necessary for non-contacting cells to create a new junction. At a high cell deformation but null deformation rate, we could imagine that the four-fold vertex lifetime is dominated by active fluctuations and adhesion remodeling. Conversely, at both high cell deformation and deformation rate, a new contact between cells can be created without maturing cell-cell adhesion (i.e. new contacts and sliding of cells are triggered by the deformation rate).

The aspired aggregate can either relax in a reservoir at the end of the channel or be pushed back into the entrance reservoir. In both cases, the relaxation remains an actual 2D relaxation as the tissue is still confined in the z direction. The relaxation dynamics qualitatively depend on the time the tissue has spent in the channel, and two relaxation regimes can be distinguished.

In turn, this new elastic deformation can eventually fully relax. This occurs at low shear rate, and long timescales, where cell divisions may also play a role.

To recapitulate the different behaviors observed in the previous sections, we propose a multi-scale rheological diagram of the aggregate (Fig. 8, Table 1). It combines the contributions of intra-cellular rheology (cell deformations), inter-cellular rheology (cell rearrangements) and aggregate surface tension (elasto-capillarity). This complete rheological model could, in principle, apply to all experiments presented here (Figs 4C, 6E and 7E). Ideally, all of them could be fitted with a single set of parameters, at least if elastic and viscous moduli are supposed constant. We, however, emphasize that only the aspiration part of the experiment is here submitted to a quantitative test. The relaxation part could be tested more quantitatively through a finite element simulation taking into account the geometry details.

Our experimental set-up uses a microfluidic constriction to impose a highly heterogeneous flow to a 3D cell aggregate, by controlling both the force applied to the aggregate and the flow geometry. It triggers (1) cell deformations, (2) cell rearrangements and (3) global aggregate geometry changes, and enables the interplay between the three to be studied. Using two-photon microscopy, we simultaneously image cell contours at subcellular resolution and quantify mechanical response at the scales of cell, cell group and whole aggregate.

Although the basic principle is similar to micropipette aggregate aspiration, the present set-up has at least three advantages. First, we can study separately, and in detail, the effect of cell group elasticity, aggregate surface tension, and coupling between both (so-called ‘elasto-capillarity’). Second, we can reach higher cell deformation, or more precisely we obtain a large product of deformation rate by relaxation time (the so-called ‘Weissenberg number’, here of order of one); this is probably relevant in fast morphogenetic changes phases such as gastrulation, whereby tissues encounter both high cell deformations and cell deformation rates. Third, we can map different scales and local differences in mechanical properties; we are thus not restricted to homogeneous aggregates.

3D cell aggregates are model systems of active fluctuating cellular materials. We propose here to measure tissue-scale effective rheological parameters averaged over time and space. Improvements could include changes in constriction walls, which could be rounded or at angles smaller than 90°. Embedding deformable stress sensors inside the tissue (Mongera et al., 2018; Campàs et al., 2014; Dolega et al., 2017; Souchaud et al., 2022) would enable cell deformation to be determined in relation to the local tissue stress to validate further rheological models; beyond the linear models we test here, future work might investigate more realistic power-law rheologies.

Further progress could come from improved time and space resolution, segmentation and tracking, to obtain better statistics; dynamic stress inference (which applies even to moving junctions) (Brodland et al., 2010; Roffay et al., 2021) could infer junction tensions and determine how they are modulated by cell deformations and rearrangements.

How cell membrane trafficking, adhesion proteins and molecular motors dynamics affect rearrangement dynamics and are relocalized during junction remodeling in 3D tissues is still not well understood (Iyer et al., 2019; Nishizawa et al., 2022 preprint). Combining our method with cell lines endowed with fluorescent reporters of such proteins would enable quantification of their dynamics using fluorescence recovery after photobleaching techniques, for example. Other exciting insights could come from blocking the action of some proteins (using knockout cell lines or chemical inhibition of various pathways and effector proteins) in order to quantify the relative contributions of various components of the cellular machinery to tissue mechanics. Studying the effects of drugs such as blebbistatin, which inhibits myosin II, could also bring new interesting pieces of information, as shown by our proof-of-concept experiments (Fig. S7B, Movie 12). In the same way, cell-cell stress propagation can also be modified for instance using an α-catenin-null cell line. Because α-catenin links the cytoskeleton with cadherins, cell-cell junctions are modified. Aggregates still form cadherin-cadherin junctions, but cell shapes are more irregular, there is a shorter correlation range in both shape and velocity, and rearrangements are dramatically impeded (Movie 13).

Cell rearrangement dynamics have been studied in in vivo systems, such as Drosophila (Guirao et al., 2015; Etournay et al., 2015; Bardet et al., 2013; Guirao and Bellaïche, 2017; Rauzi et al., 2008) or Xenopus embryos (David et al., 2014), in vitro tissues (Tlili et al., 2020) and numerically for in silico tissues (Bi et al., 2014). However, in the context of embryo development or collective cell migration in vitro, it is difficult to distinguish the rearrangements that build or relax deformation (or stress) (Guirao and Bellaïche, 2017), i.e. those for which cells act as motors, versus those for which cells resist a force imposed from outside. Here, we present an experimental platform suitable for quantitative characterization of the resistant cell rearrangements triggered by an imposed flow.

This microfluidic platform could also be a useful, simple tool for perturbing mechanosensitive pathways by applying shear stresses and to quantify their effect on biomechanical networks (using fluorescent reporters or doing RNA sequencing afterwards, as it is easy to collect the aspired aggregates after the experiment). For this kind of application, the measurement reproducibility (Fig. S4, Movies 14 and 15) enables development of parallel and/or sequential constriction set-ups, as has already been done for single-cell experiments (Brückner et al., 2019).

It has been recently shown that gradients of expression of cadherins are associated with embryonic organoid symmetry breaking and elongation (Hashmi et al., 2022). Whether these gradients are associated with gradients in rheological properties remains unknown in such systems. An advantage of our method is that it relies on the observation and extraction of local mechanical responses in 3D: it could be used to identify spatial heterogeneities in cell-group mechanical responses in systems such as embryonic organoids and correlate them with local protein or gene expression.

This set-up will be very useful for understanding the mechanisms at stake in the response of 3D tissues to mechanical stresses and in particular to obtain new insights into the biophysics of cell-cell rearrangement, which plays a central role in development. We envisage that it could be used in many different systems, such as Hydra (from which aggregates can easily be formed), mouse embryos, organoids, and 3D bulk tissues of Xenopus, and could yield new information on key developmental processes.

For these experiments, we used a mouse embryonic carcinoma F9 cell line (Nagafuchi et al., 1987) (American Type Culture Collection, CRL-1720). Aggregates were obtained by resuspending single cells after trypsinization in a dish and preventing them from settling down using an orbital shaker. After 2 days, we obtained aggregates with a diameter around 250 µm.

The experiment was performed on a microfluidic chip (Fig. 1A-C). The PDMS mold was obtained using standard lithography techniques. For the experiments presented here, the channel width was 50 µm and its height 100 µm. After bonding the PDMS chip to a glass coverslip for 60 s using a Plasma Harrick Plasma cleaner on both glass and PDMS, the channel was incubated for 30 min with a 1% Pluronic F-127 solution (Sigma-Aldrich) in order to avoid adhesion of the aggregate to the wall while flowing. It was then rinsed with Leibovitz's L-15 Medium without Phenol Red (Gibco) complemented with 10% fetal bovine serum, which was further used for the experiment. For the introduction of the aggregate, we used a 200 µl cone, and the outlet was then linked to a fluid reservoir the height of which could be adjusted in order to generate a difference in pressure controlled by the difference in height between the sample and the reservoir.

The aggregate was first dropped into a cone pipette. Then, it was sedimented into a 100-µm-high reservoir and deformed into a pancake shape (Fig. 1A-C). After waiting at least 10 min for stress relaxation, the aggregate was pushed towards the channel entrance. Because of the applied pressure difference, the aggregate slid into the reservoir towards the channel entrance and started to deform inside the channel (Fig. 1D).

The aggregates were imaged on a two-photon inverted microscope using an oil immersion 40× objective (HCX PL APO CS 40.0×1.25 OIL UV) equipped with MaiTai Spectra Physics laser (690-1060 nm) used at 800 nm, and a lowpass filter below 600 nm, under a temperature control environment set-up at 37°C. To visualize the cell contours, we used a solution of sulforhodamine B (Marmottant et al., 2009). We diluted 16 µl of a stock solution (1 mg/ml) in 1 ml of Leibovitz complete medium. With this method, bleaching is very low as the fluorophore rapidly diffuses from the medium to the extracellular space.

We measured the two-dimensional velocity field using a custom-made MATLAB optic flow code based on the Kanade–Lucas–Tomasi (KLT) algorithm (Lucas and Kanade, 1981) with a level 2 pyramid. This coarse-grained velocity field can be averaged in time, yielding components v_i(x_j) where i, j = 1 or 2. Using finite differences, we obtained in each box the deformation rate as the symmetric part of the velocity gradient, ⁠, with components (∂_iv_j+∂_jv_i)/2. This symmetric tensor can be diagonalized; we graphically represent in each box the tensors by bars, the direction and length of which correspond respectively to eigenvectors and eigenvalues, with blue bars for compression and red bars for elongation (Fig. S2).

To extract the coarse-grained cell anisotropy, we used Fourier transform (for details and validation, see Durande et al., 2019). In short, the Fourier transform of an image of a connective tissue estimates the cell size, shape and anisotropy averaged over a cell group. This method provides an efficient measurement of the deformation coarse-grained at cell-group scale, without having to recognize and segment each individual cell contour. We used the same grid and box size as for the velocity (box of 256 pixels, overlap 0.75). To improve the signal-to-noise ratio, the Fourier transform norm was smoothed in time by using a sliding average over three consecutive time points.

We defined the cell shape deformation tensor ε_cell with respect to a rest state which we assumed to be isotropic, with two equal eigenvalues L₀ = (L_maxL_min)^1/2. Hence, ε_cell has the same eigenvectors as the Fourier pattern, and two eigenvalues: and ⁠. The deformation tensor ε_cell is equivalent to other definitions of the deformation (Graner et al., 2008; Durande et al., 2019) within a linear approximation; in addition, ε_cell has the advantage of having well-established transport equations (Tlili et al., 2015).

We represent the deviator Tr(ε_cell)I/2, where I is the identity tensor; is the anisotropic part of ε_cell and only depends on the ratio ⁠, which is represented by ellipses the color of which indicates the value of the deformation.

We thank the Nanoptec platform at ILM and especially Christophe Moulin, and the ImagoSeine core facility. We thank Akira Nagafuchi for his generous gift of F9 cells, David Charalampous for developing the pressure control device, Alain Richert for help on tissue culture, Rémy Fulcrand for help on wafer fabrication in clean room, Pascal Hersen for letting us use his wet lab for microfluidic device preparation.

The peer review history is available online at https://journals.biologists.com/dev/lookup/doi/10.1242/dev.200774.reviewer-comments.pdf.