11institutetext: Dipartimento di Biotecnologie Mediche e Medicina
Traslazionale, Università degli Studi di Milano, I-20090 Segrate, Italy
22institutetext: IFOM-FIRC Institute of Molecular Oncology, I-20139 Milano,
Italy 33institutetext: Dipartimento di Oncologia e Emato-Oncologia, Università
degli Studi di Milano, I-20133, Milano, Italy 44institutetext: University of
Vienna, Faculty of Physics, 1090 Vienna, Austria

# Non-invasive measurement of nuclear relative stiffness from quantitative
analysis of microscopy data

Stefano Villa 11 Andrea Palamidessi 22 Emanuela Frittoli 22 Giorgio Scita 2233
Roberto Cerbino Corresponding author<EMAIL_ADDRESS>Fabio
Giavazzi Corresponding author<EMAIL_ADDRESS>
(Received: date / Revised version: date)

###### Abstract

The connection between the properties of a cell tissue and those of the single
constituent cells remains to be elucidated. At the purely mechanical level,
the degree of rigidity of different cellular components, such as the nucleus
and the cytoplasm, modulates the interplay between the cell inner processes
and the external environment, while simultaneously mediating the mechanical
interactions between neighbouring cells. Being able to quantify the
correlation between single-cell and tissue properties would improve our
mechanobiological understanding of cell tissues. Here we develop a methodology
to quantitavely extract a set of structural and motility parameters from the
analysis of time-lapse movies of nuclei belonging to jammed and flocking cell
monolayers. We then study in detail the correlation between the dynamical
state of the tissue and the deformation of the nuclei. We observe that the
nuclear deformation rate linearly correlates with the local divergence of the
velocity field, which leads to a non-invasive estimate of the elastic modulus
of the nucleus relative to the one of the cytoplasm. We also find that nuclei
belonging to flocking monolayers, subjected to larger mechanical
perturbations, are about two time stiffer than nuclei belonging to dynamically
arrested monolayers, in agreement with atomic force microscopy results. Our
results demonstrate a non-invasive route to the determination of nuclear
relative stiffness for cells in a monolayer.

##  1 Introduction

In all the phases of their life cycle, biological tissues are dynamic
entities: during morphogenesis, cells move and find their optimal location and
shape in the organized spatial distribution of the tissue; in homeostatic
equilibrium, dying cells are continuously extruded and replaced by new ones,
but a variety of processes - such as tissue repair, inflammatory response,
carcinogenesis, and tumor progression - can drive the tissue out-of-
equilibrium. In all these cases, the dynamic state of the tissue changes as a
result of the appropriate physical and chemical stimuli being exchanged at
different scales, from the one of the single cells to the one involving the
entire tissue. However, the details of this process remaining elusive.

With the eyes of a physicist, these changes of state of the tissue can be
thought of as phase transitions that bring the tissue from one state to a
different one. Focusing on cancer, it is possible to consider cancer
progression towards metastasis and invasion as a sequence of steps involving a
transition of the tissue from a solid-like to a fluid-like one; this
fluidization transition is also referred to as unjamming park2015unjamming ;
sadati2013collective ; blauth2021jamming ; malinverno2017endocytic ;
palamidessi2019unjamming ; kim2020unjamming ; oswald2017jamming ;
ilina2020cell . The fluid tissue phenotype is typically accompanied at the
single cell level by a softening of the cell, which may impact the
infiltration of tumor cells alibert2017cancer ; swaminathan2011mechanical ,
and at the same time, by an increased cell motility, which makes it easier for
sick cells to invade surrounding healthy tissues. Since forces external to the
cells propagate through the cytoplasm directly to the nucleus
wang2009mechanotransduction , larger stresses within the tissue resulting from
unjamming can lead to nuclear envelope rupture and consequently lead to DNA
leakage and damage denais2016nuclear ; raab2016escrt ; xia2018nuclear ;
pfeifer2018constricted . In this framework, nuclear mechanical properties and
shape are important regulators of the state of a tissue grosser2021cell ; for
example, non-trivial correlations between nuclear stiffness and tissue
dynamics have been recently observed harada2014nuclear ; parreira2021single .

As a consequence of the complexity of considering cellular tissues as
materials, an exhaustive characterization of their mechanical and rheological
properties is extremely challenging. One of the main reasons is the intricacy
of the cell, which results from the composition of several elements: cell
membrane, cytoplasm, cytoskeleton, nuclear envelope, organelles, nucleoplasm,
chromatine, etc.. Moreover, cell mechanical properties depends on a variety of
factors, including the cell status, the substrate, the treatments and the
applied stress fischer2020effect ; liu2014situ ; fabry2003time ;
guilak2000viscoelastic ; darling2015high ; harris2012characterizing .
Modelization of the cell as an ensemble of few elastic and viscous elements is
often incomplete and unrealistic alibert2017cancer . To further complicate the
picture, different techniques for the measurement of the cell viscoelastic
properties have been deployed (e.g. AFM radmacher1992molecules ;
fischer2020effect ; liu2014situ , micropipette aspiration
mitchison1954mechanical ; rowat2005characterization , optical tweezers
titushkin2006distinct , magnetic twist cytometry fabry2003time ), and they
generally obtain different values for these properties, as they often probe
different cellular components darling2015high ; lim2006mechanical .

Specifically for the cell nucleus, when considering its stiffness within the
force field of a tissue one expects values that differ substantially from the
ones obtained experimentally in controlled settings and, moreover, different
cell lines can have very different mechanical properties, with cells in
similar situations exhibiting nuclei either stiffer or softener than cytoplasm
fischer2020effect . Recently, Parreira et al. parreira2021single  have shown
how the presence of a single stiffer nucleus alters the dynamics of a whole
monolayer, impairing collective migration. Through a refined image
segmentation of nuclei, allowing the collection of several dynamical
parameters, the Authors show that the decrease in net migration velocity is
counterbalanced by an increase of the rotational dynamics of cell nuclei,
accompanied by an increment of nuclear deformation. Interestingly, nuclear
deformation affects all the cells and not only the first neighbors of the
stiffer nucleus, thus suggesting that the presence of the obstacle influence
nuclei far away from it. These results provide an important motivation for the
development of non-invasive in situ tools for quantifying the nuclear
stiffness as a proxy of mechano-biological alterations of the tissue.

In this work, we propose an automated procedure enabling the extraction of
time-resolved information on morphology and dynamics of the nuclei in a
monolayer, and the investigation of the interplay between nuclear shape
fluctuations and collective motility within the cell aggregate. A custom
algorithm is developed to identify, segment and track fluorescent nuclei over
time in time-lapse microscopy image sequences. The data are subsequently
analyzed to investigate correlations between nuclear deformations and
kinematics parameters characterizing the monolayer. Our estimators allow
extracting non-invasively a robust estimate of the relative stiffness between
the nucleus and the cytoplasm. We test our approach with different cell lines
approaching cell jamming, as well as on tissues undergoing an unjamming
transition via flocking induced by the overexpression of the RAB5A protein
malinverno2017endocytic .

##  2 Materials and methods

###  2.1 Sample preparation and imaging

MCF10A cells were a kind gift of J. S. Brugge (Department of Cell Biology,
Harvard Medical School, Boston, USA) and were maintained in Dulbecco’s
Modified Eagle Medium: Nutrient Mixture F-12 (DMEM/F12) medium (Gibco)
supplemented with $5\%$ horse serum, 0.5 mg ml-1 hydrocortisone, 100 ng ml-1
cholera toxin, 10 $\mu$g ml-1 insulin and 20 ng ml-1 EGF. The cell line was
authenticated by cell fingerprinting and tested for mycoplasma contamination.
Cells were grown at 37 °C in humidified atmosphere with $5\%$
$\text{CO}_{\text{2}}$. MCF10A cells were infected with pSLIK-neo-EV (empty
vector control) or pSLIK-neo-RAB5A lentiviruses and selected with the
appropriate antibiotic to obtain inducible, stable cell lines. Constitutive
expression of EGFP-H2B was achieved by retroviral infection of MCF-10A cells
with pBABE-puro-EGFP-H2B vector.

MCF10.DCIS.com cells obtained from Dr John F Marshall (Barts Cancer Institute,
London, UK) and were maintained in Dulbecco’s Modified Eagle Medium: Nutrient
Mixture F-12 (DMEM/F12) medium (Gibco) supplemented with $5\%$ horse serum,
0.5 mg ml-1 hydrocortisone, 10 $\mu$g ml-1 insulin and 20 ng ml-1 EGF.
MCF10.DCIS.com cells were infected with pSLIK-neo-EV (empty vector control) or
pSLIK-neo-RAB5A lentiviruses and selected with the appropriate antibiotic to
obtain stable inducible cell lines. Constitutive expression of mCherry-H2B was
achieved by retroviral infection of MCF10DCIS.com cells with pBABE- puro-
mCherry-H2B vector.

MCF10AneoT and MCF10CA1 cells were obtained from Dr. Polin Lisa (Barbara Ann
Karmanos Cancer Institute, Detroit, USA) and were maintained in Dulbecco’s
Modified Eagle Medium : Nutrient Mixture F-12 (DMEM/F12) medium (Gibco)
supplemented with $5\%$ horse serum, 1.05 mM $\text{CaCl}_{\text{2}}$, 10 mM
Hepes, 0.5 mg ml-1 hydrocortisone, 10 $\mu$g ml-1 insulin and 20 ng ml-1 EGF.
Constitutive expression of mCherry-H2B was achieved by retroviral infection
with pBABE-puro-mCherry-H2B vector.

Transfections were performed using either calcium phosphate or FuGENE HD
Transfection Reagent (cat. no. E2311, Promega), according to the
manufacturer’s instructions. Lentiviral and retroviral infections were
performed as described if Ref. malinverno2017endocytic . Cells were seeded in
six-well plates in complete medium and cultured until a uniform monolayer had
formed. For all the experiments, the seeding density was $1.5$ × $10^{6}$
cells per well, with the exception of the experiments indicated as LC (low
confluency), for which the seeding density was $0.75$ × $10^{5}$ cells per
well.

RAB5A expression was induced, where indicated, 16 h before performing the
experiment by adding fresh complete media supplemented with $2.5$ g ml-1
doxycycline to cells. At the time of recording, fresh media containing EGF and
doxycycline was added. After cell induction, doxycycline was maintained in the
media for the total duration of the time-lapse experiment.

An Olympus ScanR inverted microscope with a 5× objective (experiments on
jamming monolayers of MCF10.DCIS.com, MCF10AneoT, and MCF10CA1a cells) or a
10× objective (all the other experiments) was used to acquire images with a
frame rate of $0.5$ frames/min (experiments on jamming monolayers), $0.1$
frames/min (experiment on MCF10A flocking monolayers) and $0.4$ frames/min
(experiment on MCF10A.DCIS.com flocking monolayers) minutes over a $24$ h
period. For each sample, 4 or 5 independent FOVs, much smaller than the entire
culture plates and far from the boundaries, are captured. Each FOV is imaged
both in fluorescence and in phase-contrast microscopy.

###  2.2 Nuclear tracking and segmentation

In order to extract both static and dynamic features of the monolayer at the
single nucleus level, a Matlab algorithm was developed for the automatic
nuclear segmentation and tracking. The algorithm operates on time-lapse
microscopy images of monolayers of cells with fluorescently tagged nuclei. The
main steps of the algorithm are detailed below.

####  2.2.1 Image background subtraction and registration

Images are first corrected for background intensity inhomogeneities and
microscope stage positioning errors.

Background subtraction is obtained through the implementation of the iterative
procedure described in Ref. wahlby2002 . We first perform a cubic
interpolation of the intensity in each image, by randomly choosing a set of
points within a mask which initially coincides with the whole image. After the
interpolation, the mask is redefined by excluding the points corresponding to
high intensity values, which are likely to belong to a nucleus rather than to
the background, and a second interpolation is performed based on the redefined
mask. After an adequate number of iterations (in our case ten), the mask
should include only regions belonging to background, evenly distributed across
the image. The final interpolation step provides thus a reliable estimate of
the background intensity distribution. This process is repeated over different
frames equally spaced in time, covering the whole duration of the experiment.
The average of all the so-obtained background images represents our best
estimate for the background intensity distribution, which is thus subtracted
to each image in the sequence.

Microscope stage translations between consecutive frames, due to positioning
errors while imaging different locations at each time point, can introduce a
bias in the reconstructed monolayer dynamics, by adding a spurious
contribution to the velocity within the field of view (FOV). In order to
address this issue, the Matlab intensity-based registration function
imregister is used to recover the displacement
$\Delta\vec{x}_{reg}\left(t_{i}\right)$ of the FOV center of mass between each
pair of consecutive frames $(i-1)$ and $i$, considering pure translations. We
then calculate the trajectory
$\vec{x}_{reg}\left(t_{i}\right)=\sum_{i^{\prime}=2}^{i}\Delta\vec{x}_{reg}\left(t_{i^{\prime}}\right)$
of the FOV center of mass. In the presence of both directed cell migration and
stage movement, $\vec{x}_{reg}$ can be written as
$\vec{x}_{reg}\left(t_{i}\right)=\vec{x}_{cm}\left(t_{i}\right)+\vec{x}_{noise}\left(t_{i}\right)$,
where $\vec{x}_{noise}\left(t_{i}\right)$ a the random, delta-correlated,
noise associated to stage movement, and $\vec{x}_{cm}\left(t_{i}\right)$ is
the genuine displacement of the center of mass of the cell monolayer, which we
estimate with a twentieth degree polynomial fit of $\vec{x}_{reg}$.
Subtraction of this term from $\vec{x}_{reg}$ enables isolating the noise
contribution $\vec{x}_{noise}$.

Figure 1:  a. Close-up view of a fluorescence microscopy image of a
MCF10.DCIS.com monolayer during a jamming experiment. The yellow bar
corresponds to $10$ µm. b. $L_{G}$ map resulting from the application of LoG
filtering to (a), after the Wiener filter application. Colorscale goes from
blue to red passing from light-blue and yellow. c. Binary map $L_{BW}$
obtained by $L_{G}$f via a thresholding operation, as described in the text.
Red points marks the center of mass of each connected domain of white pixels.
d. Euclidean distance transform $L_{ED}$ of the $L_{BW}$ map reported in (c).
Yellow points mark the result of watershed segmentation of $L_{ED}$. e.
Internal (red) and external (yellow) seed points obtained as shown in panel
(c) and (d) respectively are superimposed on the original image. The result of
seeded watershed segmentation made on the gradient of (a) is also shown in
green. f. Result of subpixel segmentation at the end of segmentation process.

####  2.2.2 Nuclear segmentation procedure

Once a registered stack of background-subtracted images $I\left(t_{i}\right)$
is obtained, we process each frame $i$ to identify single nuclei. Typical
nuclei segmentation methods identify nuclei either localizing their centers
(for example, looking for intensity maxima or centers of symmetry) or through
the identification of nuclear edges xing2016robust . In our case, the
identification of the nuclei from their centers is made extremely difficult by
the presence of intensity heterogeneity mirroring the chromatin nuclear
distribution (Fig. 1a); we thus prefer to rely on nuclear edges rather than on
nuclear centers. To identify nuclear edges, we apply a watershed transform to
the image spatial gradient $\nabla I\left(t_{i}\right)$, which locates nuclear
edges of $I\left(t_{i}\right)$ as ridge lines. However, direct watershed
segmentation is not very efficient in properly segmenting the nuclei, in
particular in the case of jamming monolayers, where the signal-to-noise ratio
is relatively low and partial superposition of two or more nuclei is
relatively frequent. In order to improve the quality of the result, we
introduce, before the application of watershed transform, a pre-processing
where the image gradient $\nabla I\left(t_{i}\right)$ is set to zero in
correspondence of suitable ”seed points” located both within and outside the
boundary of each nucleus. The procedure leading to the identification of the
seed points is described in detail in the following.

Noise in each frame $I\left(t_{i}\right)$ is first reduced by applying a
Wiener filter, an adaptive noise-removal filtering that preserves nuclei edges
lim1990two . A Laplacian-of-the-gaussian (LoG) filter is then applied to
enhance nuclear edges. In the obtained map $L_{G}$ (see Fig. 1.b), higher
values (towards yellow), corresponding to nuclear edges of the original image,
surround deep intensity wells (blue), within which the nuclei centers of mass
locate. Differences in the fluorescent intensity of different nuclei are
reduced by dividing $L_{G}$ by an intensity map obtained via bicubic
interpolation of the local minima of $L_{G}$. The corrected $L_{G}$ is then
binarized through the application of a suitable threshold $k_{th}$, setting to
$0$ ($1$) the pixels whose intensity is larger (lower) than $k_{th}$. The
threshold value is determined as the one maximizing the number of domains of
connected pixels equal to $1$ (Fig. 2.d). Indeed, as it can be appreciated by
comparing figures 2.a and 2.b, when, starting from the minimum value of
$L_{G}$, the value of $k_{th}$ is increased, the number of connected domains
initially increases as well. This is due to the fact that the correction based
on bicubic interpolation of the local minima of $L_{G}$ does not perfectly
level out the intensity differences between nuclei. On the other hand, when
$k_{th}$ reaches the typical intra-nuclear value of $L_{G}$, domains
corresponding to different nuclei starts to merge and their number tends to
decrease (Fig. 2.b-c).

Figure 2:  a-b-c. Examples of the details of $L_{BW}$ resulting from the
application of different threshold to $L_{G}$ in increasing order of $k_{th}$.
Binarized images are overplotted with different colors on the original detail
of $I$. Yellow bar correspond to 10 µm. d. Number of connected domains of
pixels recovered in $L_{BW}$ as a function of the threshold $k_{th}$. Colored
dots mark the same thresholds represented in (a), (b) and (c).

Once the proper threshold $k_{th}$ is imposed, the resulting binarized image
$L_{BW}$ (Fig. 1.c) is an ensemble of domains of connected pixels, each one of
which is a first approximation of the inner part of each nucleus.

On $L_{BW}$ additional minor operations are then performed to correct false
counting and spurious merging of neighbouring nuclei. False counting are in
most cases due to background noise and results in domains of few connected
pixels. A morphological ”open” operation (an erosion followed by a dilation)
has been implemented to remove these artifacts.

Three examples of mistakenly merged nuclei are marked in red in figure 3. In
order to correct this problem, the possible cases in which it may be occurred
are automatically identified considering the domains with an aspect ratio
larger than 3 (against a typical distribution of domains aspect ratios with
median around 1.4 and standard deviation of 0.5)111The domains’ aspect ratio
limit has been determined through experience, has no real nuclei has been
found in the analyzed data set with larger aspect ratios. It can be eventually
adjusted for samples with different domains’ aspect ratio distribution.. In
order to roughly localize in them the centers of the merged nuclei, for each
selected domain the two deeper local minima of $L_{G}$ are identified and two
disks of radius 2 pixels are centered on them, replacing the connected pixels
domain. The result of such correction is illustrated in figure 3, where bi-
nuclear domains are outlined in red, while the domains obtained upon
correction are marked in green222Note that for a proper segmentation the non-
perfect centering of the domain within the nucleus does not represent a
problem: it is enough that each nucleus contains one and only one entire
domain. The precise center of each nucleus is determined after nuclear
segmentation..

Figure 3:  a-c. Examples of successful separation of domains of connected
pixels belonging to different nuclei. Red shadows highlight domains identified
with the LoG filter procedure, which are corrected with the subsequent
separation step. Green shadows mark the domains determined at the end of the
separation procedure. Yellow segments mark 5 µm scale.

The geometrical centers of mass of the domains of connected pixels are chosen
as inner seeds for the watershed segmentation (red points in Fig.1.c).

The determination of the seeds external to the nuclei is simpler. The
connected domains in $L_{BW}$ do not only indicate the position of the
corresponding nuclei, but have also an area which is typically larger for
nuclei with wider projected area. Consequently, the space between two nuclei
in $I\left(t_{i}\right)$ is typically equidistant from the edges of the
corresponding connected domains in $L_{BW}$. For the identification of the
external seeds, we first determine the Euclidean distance transform $L_{ED}$
of $L_{BW}$ maurer2003linear , a map where each pixel value is the distance
from the closest domain (Fig.1.d). On $L_{ED}$ a watershed segmentation is
finally applied to determine the equidistance lines between noighbouring
domains in $L_{BW}$ (yellow points in Fig.1.d), which are then used as
external seeds.

In figure 1.e the recovered internal (red) and external (yellow) seeds are
superimposed on the original image $I\left(t_{i}\right)$. By setting $\nabla
I\left(t_{i}\right)$ to zero in correspondence of all seeds and evaluating the
watershed transform, the nuclei segmentation is finally obtained (Fig. 1.d,
green pixels). The recovered segmentation is pixel resolved; sub-pixel
resolution can be achieved by running the Matlab code subpixelEdges 2.13
trujillo2013accurate  on the obtained profiles (Fig. 1.e).

After nuclei segmentation, several static quantities are easily accessible. In
this work, we consider in particular the center of mass $\vec{x}_{j}$, the
projected area $A_{j}$ and the aspect ratio $\AR_{j}$ of the $j$-th nucleus.
The latter is defined as the ratio between the nuclear major and minor axis,
evaluated as the square root of the eigenvalues of the covariance matrix of
the segmented objects jolliffe1986principal .

####  2.2.3 Dynamic parameters

Time-resolved information on nuclear shape and mobility is recovered by
linking the positions $\vec{x}_{j}$ of each nucleus in subsequent frames. To
this aim we employed the publicly available Matlab implementation by D Blair
and E Dufresne blair2008matlab  of the linking algorithm of JC Grier and DG
Crocker crocker1996methods . Once nuclear trajectories are built, the time
evolution of the single-nucleus parameters $\vec{x}_{j}\left(t_{i}\right)$,
$A_{j}\left(t_{i}\right)$, $p_{j}\left(t_{i}\right)$ and
$\AR_{j}\left(t_{i}\right)$ can be also obtained. The instantaneous velocity
on the $j$-th nucleus across frames $i$ and $(i+1)$ is estimated as
$\vec{v}_{j}\left(t\right)=\left[\vec{x}_{j}\left(t_{i+1}\right)-\vec{x}_{j}\left(t_{i}\right)\right]/\delta
t$, where $\delta t$ is the time between two acquired frames.

Our definition of the instantaneous velocity in terms of the average velocity
over a time interval equal to $\delta t$ relies on the fact that, on this
timescale, the cellular motion can be assumed to take place with approximately
constant velocity. We checked this assumption by measuring the nuclear mean
square displacement and estimating the characteristic persistence time of the
ballistic-like motion observed for short time delays. As discussed in detail
in Appendix A, for all the considered data sets, the persistence time was
found to be larger than $\delta t$.

The velocity of the center of mass of the monolayer is evaluated as the
instantaneous mean velocity of the cells within the FOV

$\vec{v}_{CM}\left(t\right)=\left<\vec{v}_{j}\left(t\right)\right>_{j},$ (1)

where $\left<...\right>_{j}$ denotes the average over all the nuclei in the
field of view. The amplitude of velocity fluctuations in the monolayer is
evaluated as the root mean square velocity of the nuclei in the center of mass
reference frame:

$v_{RMS}\left(t\right)=\sqrt{\left<|\vec{v}_{j}\left(t\right)-\vec{v}_{CM}\left(t\right)|^{2}\right>_{j}}.$
(2)

Comparison between nuclei tracked automatically by the algorithm and manually
by an operator reveals that the described segmentation procedure is effective
in identifying about $80$-$90$% and $90$-$95$% of the nuclei present in each
FOV for the jamming and the flocking monolayers, respectively. This difference
in tracking efficiency reflects differences in noise level and spatial
resolution between different experiments (see previous section for details).
Despite the effort to minimize segmentation artifacts, such as multiple
segmentation of the same nucleus and spurious merging of multiple nuclei, such
segmentation errors occur, especially in those cases where the signal-to-noise
ratio is low or partial superposition of different nuclei is frequent. To
minimize the impact of segmentation errors on the analysis of nuclear
features, a ”quality check” has been implemented to select a subset of
reliably segmented nuclei. To this end, we evaluate the total instantaneous
intensity $J_{j}\left(t_{i}\right)$, which is obtained as the sum of the
intensity of all pixel within the segmented area of the $j$-th nucleus at
frame $i$. The instantaneous value $J_{j}\left(t_{i}\right)$ is then compared
with its median evaluated over the previous 10 frames. If the difference
between $J_{j}\left(t_{i}\right)$ and the median is larger than $10$%, the
segmentation of the $j$-th nucleus at frame $i$ is considered unreliable, and
the corresponding parameters are not included in the statistics. Trajectories
which, after the application of this quality filter, lose more than $20$% of
frames are entirely excluded. The number of nuclei that ”pass” the quality
check varies from sample to sample between 600 and 5000, depending on the FOVs
size and image quality.

###  2.3 Particle image velocimetry

Particle image velocimetry (PIV) of fluorescent microscopy images of confluent
monolayers is performed by using the Matlab PIVLab software
thielicke2021particle . We choose an interrogation area with size slightly
larger than the average inter-nuclear distance, corresponding to approximately
$14$ µm. Outliers in the reconstructed velocity field, whose components
$v_{x}$ and $v_{y}$ exceed fixed threshold values, are identified and replaced
with the median value of the velocity over neighboring grid points 333For each
component $v_{\xi}$ of the velocity, upper and lower bounds are identified as
$med(v_{\xi})-10\cdot\sigma_{v_{\xi}}$ and
$med(v_{\xi})+10\cdot\sigma_{v_{\xi}}$, respectively, where $med(v_{\xi})$ and
$\sigma_{v_{\xi}}$ are the median and the standard deviation of $v_{\xi}$
evaluated over the entire time lapse .. From the velocity field obtained from
PIV, $\vec{v}_{CM}\left(t\right)$ and $v_{RMS}\left(t\right)$ are recovered as
for PT averaging over grid points coordinates $\vec{x}$ in place of $j$.
Additionally, the local divergence of the velocity field
$\nabla\cdot\vec{v}\left(\vec{x},t\right)$ is computed from
$v_{x}\left(\vec{x},t\right)$ and $v_{y}\left(\vec{x},t\right)$ using the
Matlab divergence function.

For the same samples, the PIV algorithm has been tested on both fluorescent
and phase contrast time lapses. As shown in appendix B, both the PIV analysis
return velocity mean values and distributions in agreement with the ones
obtained from PT, thus validating the obtained dynamics.

##  3 Results and discussion

###  3.1 Jamming monolayers

A first set of experiments was performed on mature highly confluent monolayers
undergoing a jamming transition, a progressive slowing down towards a
dynamically arrested state angelini2011glass ; garcia2015physics ;
cerbino2021disentangling .

We considered three distinct MCF10A-derived cell lines (MCF10.DCIS.com,
MCF10AneoT, and MCF10CA1a cells), each one of those seeded at two different
densities (see section 2.1 for details).

####  3.1.1 Confluent monolayers evolve toward a high density kinetically
arrested state

In figure 4 we report the time evolution of the RMS velocity ${v}_{RMS}$ (fig.
4.a) and of the center of mass velocity (fig. 4.b) for the jamming monolayers,
as obtained from PIV. Each color refers to a different cell line (thick and
thin lines correspond to high and low seeding density, respectively). For each
cell line and seeding density, 5 independent FOVs are considered.

The general decreasing trend of both the considered indicators towards a
plateau follows the expected dynamic arrest characterizing jamming transition
and is accompanied by a progressive reduction of the mean cell area, evaluated
of the order of $30\%$ 444The cell area reduction as been evaluated through a
manual counting of the number of nuclei within a FOV of fixed size at the
beginning and at the end of the time lapses., as the number of cells grow.

The expected increasing density within the cell monolayer can be observed in
figure 5.a-d where four different snapshots of the same FOV of an
MCF10A.DCIS.com monolayer are reported at times $0$ h, $8$ h, $16$ h and $24$
h.

Figure 4: Time evolution of RMS velocity ${v}_{RMS}$ (a) and of the modulus of
the center of mass velocity (b) obtained from PIV on fluorescent time lapses.
Shadowed error bars are evaluated as the standard deviation of the mean
evaluated for each sample over the 5 different FOVs. Reported data refer to
jamming samples: DCIS HC (red thick line) and LC (dark red thin line), MCF10AT
HC (blue thick line) and LC (dark blue thin line), MCF10CA HC (green thick
line) and LC (dark green thin line).

####  3.1.2 Dynamic arrest correlates with a nuclear projected area decrease
and nuclear aspect ratio increase

Figure 5:  a-d. Four different snapshot of the same FOV region of DCIS HC
jamming experiment. Snapshots are taken at times $0$ h (a), $8$ h (b), $16$ h
(c) and $24$ h (d). e, f. Time evolution of the projected nuclear area (e) and
aspect ratio (f) for DCIS HC (red thick line) and LC (dark red thin line), for
MCF10AT HC (blue thick line) and LC (dark blue thin line) and for MCF10CA HC
(green thick line) and LC (dark green thin line). g. Scatter plot of single
cell instantaneous measured projected area versus the corresponding aspect
ratio. Data are obtained from a DCIS HC jamming monolayer and colors, in order
from dark red to yellow, refer to time intervals $0$-$8$h, $4$-$12$h,
$8$-$16$h, $12$-$20$h and $16$-$24$h, respectively. Reported points are
averages evaluated over equally spaced bins along the horizontal axis.

Nuclear average projected area, reported in figure 5.e, has the same
decreasing trend in time of cell area, progressively reducing due to cell
proliferation.

Besides these trends, we are here interested in understanding the interplay
between the dynamics of the monolayer and changes in nuclear shape at the
single cell level.

On top of the monotonic trends described above, each nucleus undergoes shape
fluctuations, as it can be seen in figure 6 where the contours of the same
nucleus are reported for different times.

Figure 6:  Example of resulting subpixel segmentation of a nucleus followed in
time. Images obtained from MCF10A RAB5A overexpressing monolayer in time
interval $0$-$5$ hours. Color scale indicates elapsing time from blue to
yellow. a-i. Different images centered on the same tracked nucleus followed in
time (timestep $30$ minutes) with the corresponding segmentation overplotted
on the original image. White scale bar correspond to $5$ µm. l. Overplot of
segmentations at different times (time step $10$ minutes) relative to the same
nucleus of (a)-(i). Radial axis units are micrometers. m. The same of (l), but
here each segmentation is corrected by a rigid rotation equal to the recovered
nuclear orientation in order to remove rotational contribution and only
visualize deformations.

Two parameters we have direct access to in order to characterize nuclear shape
fluctuations are the area and the aspect ratio. They appear to be strongly
correlated. As it can be appreciated from figure 5.f, where we report the time
evolution of the average nuclear aspect ratio, the aspect ratio displays a
monotonically increasing trend with time. This indicates that, while nuclear
projected area reduces as the density increases, nuclear shape anisotropy
tends to increase. In order to deepen this observation, we consider the
association between the instantaneous aspect ratio $\AR_{j}$ and the
instantaneous area $A_{j}$ at the single nucleus level. An example of the
resulting scatter plot, binned in equally spaced intervals of area, is
reported in figure 5.g. Different colors refer to 5 subsequent and overlapping
time intervals of $8$ hours covering the whole duration (24h) of the
experiment. A striking negative correlation is observed, proving that in
nuclei smaller projected areas are systematically more anisotropic. This
observation is in agreement with what observed in Ref. haase2016extracellular
. No great differences emerge for different times, suggesting that, at a first
glance, the correlation between $A$ and $\AR$ is independent on age and
monolayer dynamics.

Since the two parameters are strongly correlated, in the following we only
focus on the projected area, whose value can be determined with lower
uncertainty.

####  3.1.3 Nuclear deformation rates

Due to the non stationarity of the mean projected area on the experimental
time scales, deformation amplitudes are difficult to be uniquely retrieved, as
the obtained value are strongly dependent on the way the mean dynamics of the
projected area is subtracted. In analogy with the mean squared displacement,
we therefore introduce the mean square strain as

$MSS\left(\tau\right)=\left<\left<\Delta
a_{j}^{2}\left(\tau|t\right)\right>_{t}\right>_{j},$ (3)

where the nuclear strain $\Delta a_{j}\left(\tau|t\right)$ of the $j$-th
nucleus between time $t$ and $t+\tau$ is calculated as

$\Delta
a_{j}\left(\tau|t\right)=\frac{\left[A_{j}\left(t+\tau\right)-A_{j}\left(t\right)\right]}{\left<A_{j}\left(t\right)\right>_{t}}.$
(4)

Representative $MSS\left(\tau\right)$ curves obtained for the DCIS HC
monolayer are reported in figure 7.a. Each curve refers to one of the five,
partially overlapping time intervals the experimental window was divided into.
The linearity of $MSS\left(\tau\right)$ at low $\tau$ points to a diffusive-
like evolution of the area at short time scales. For larger time delays,
$MSS\left(\tau\right)$ becomes sublinear. The limited experimental time window
does not allow establishing whether an asymptotic plateau value $\sim 1$ is
eventually attained. In order to extract the key parameters characterizing
nuclear deformation, we use an exponential model function
$MSS\left(\tau\right)=\sigma_{w}+\dot{a}_{0}\tau_{c}\left(1-\\-e^{-\tau/\tau_{c}}\right)$
to fit the data. The diffusive-like growth of the area fluctuations is
captured by the model in the limit of small $\tau$:
$MSS\left(\tau\right)\sim\dot{a}_{0}\tau$, with a characteristic strain rate
$\dot{a}_{0}$. The short-time regime is followed by an exponential-like
relaxation trend towards a plateau value $\dot{a}_{0}\tau_{c}$. The model also
includes an offset $\sigma_{w}$ accounting for random, delta-correlated noise
associated to the determination of projected area.

Figure 7:  a. Mean square normalized area relative to jamming DCIS HC
monolayers. Different colors, in order from dark red to yellow, refer to time
intervals $0$-$8$h, $4$-$12$h, $8$-$16$h, $12$-$20$h and $16$-$24$h,
respectively. Reported data are scaled to the offset noise $\sigma_{w}$
obtained from the fit. Lines represent corresponding fits with the model
function described in the text, also scaled on $\sigma_{w}$. b. Distribution
of the divergence evaluated interpolating on nuclei centers of mass data from
PIV on fluorescent images. Data refer to the same sample and the same time
intervals of (a). c. Scatter plot of the instantaneous single nucleus strain
rate $\dot{a}_{n}$ as a function of the corresponding divergence. Reported
points are average over hundreds/thousands of strain rate-divergence couples
made over bins equally spaced along the abscissae axis. Error bars are
evaluated as standard deviations of the mean of the data contained in the
bins. Data refer to the same sample and the same time intervals of (a). d.
Nuclear characteristic strain rates as a function of time, obtained from the
fit of the corresponding MSS relative to the jamming experiments. Time
intervals are the same of (a). e. Divergence standard deviation as a function
of time. Time intervals are the same of (a). f. Reciprocal $e_{n}$ of the
slope of $\dot{a}_{n}$ versus $\nabla\cdot\vec{v}$ scatter plot. Time
intervals are the same of (a).

In figure 7.d we report the values of $\dot{a}_{0}$ obtained from the fitting
procedure for all the jamming experiments, at different time points. A common
feature that can be observed in all datasets, is the progressive decrease in
the nuclear strain rate over time, indicating that mechanical deformation of
nuclei decreases as kinetic arrest is approached (see also figure 4).

####  3.1.4 Nuclear deformation rates correlate with the local dynamics of the
monolayer

Figure 8:  a,b. Overplot of the divergence field over two FOV of an MCF10A
cell monolayer. Divergence filed is represented as a shaded colormap where the
passage from negative to positive values of divergence are marked by the color
progression from red to white. Arrows represent the velocity field. Blue
shades higlights two segmented nuclei subjected to a negative (a) and to a
positive (b) divergence.

In order to investigate in more detail the interplay between nuclear
deformation and cell motility, we consider the divergence of the velocity
field as a suitable parameter to capture local density fluctuations (i.e.
compressions or dilations) within the monolayer, as schematically illustrated
in figure 8. We evaluate the divergence field by computing for each time $t$
the divergence $\nabla\cdot\vec{v}\left(t,\vec{x}\right)$ of the velocity
field obtained from PIV analysis (see appendix B for details). We then
estimated the value $\nabla\cdot\vec{v}_{j}\left(t\right)$ of the divergence
in the position $\vec{x}_{j}\left(t\right)$ corresponding to the center of
mass of the tracked nuclei via cubic interpolation.

In figure 7.b we report the experimentally determined probability distribution
functions (PDFs) of $\nabla\cdot\vec{v}_{j}\left(t\right)$ for the DCIS HC
monolayer over the same time intervals considered also in panel a. We observe
that the PDF of the divergence tends to become narrower and narrower as the
monolayer ages. This can be also appreciated from figure 7.e, where the
standard deviation $\sigma_{\nabla}$ of the PDF of the divergence is reported
as a function of time. The overall trend is very similar to the one found for
$\dot{a}_{0}$.

The origin of the similarity between the average values reported in figures
7.d and 7.e can be investigated by considering the direct association between
the projected area of the $j$-th nucleus at a given time and the corresponding
instantaneous local value of the divergence
$\nabla\cdot\vec{v}_{j}\left(t\right)$. To this end, we computed the
instantaneous single nucleus area strain rate, defined as

$\dot{a}_{j}\left(t\right)=\frac{2}{{\color[rgb]{0,0,0}\delta
t}}\frac{A_{j}\left(t_{i+1}\right)-A_{j}\left(t_{i}\right)}{A_{j}\left(t_{i+1}\right)+A_{j}\left(t_{i}\right)}$
(5)

In order to probe the correlation between $\dot{a}_{j}\left(t\right)$ and
$\nabla\cdot\vec{v}_{j}\left(t\right)$, we evaluate the corresponding time
cross-correlation function $K_{\dot{a},\nabla\vec{v}}\left(\tau\right)$
averaged over the nuclei in the FOVs. For the precise definition of
$K_{\dot{a},\nabla\vec{v}}\left(\tau\right)$ and for consideration and
comparison with corresponding self correlations, see appendix C. In figure
18.c of appendix C we report the cross-correlation function for both negative
and positive values of time delay $\tau$. As it can be appreciated from the
figure, a marked cross-correlation peak, centered in $\tau=0$, is present. The
characteristic width of this peak varies between $5$ and $20$ minutes,
depending on the cell line and on the age. The fairly symmetric shape of the
cross-correlation peak with respect to the vertical axis indicates that there
is no systematic delay between divergence and nuclear deformation, at least
within the experimental temporal resolution ($2$ minutes).

The evidence of an instantaneous correlation between
$\dot{a}_{j}\left(t\right)$ and $\nabla\cdot\vec{v}_{j}\left(t\right)$ prompts
us to further investigate it to understand its character and origin. For each
one of the $8$ hours time intervals, we therefore generate a scatter plot of
$\dot{a}_{j}\left(t\right)$ versus $\nabla\cdot\vec{v}_{j}\left(t\right)$. To
make the plot more readable, we introduced a uniform binning of the horizontal
axis. For each bin, containing on average hundreds or thousands of points of
the scatter plot, we report in figure 7.c the average value $\dot{a}_{n}$ of
the strain rate as a function of the average value $\nabla\cdot\vec{v}$ of the
divergence. Remarkably, at least in the vicinity of the origin, the obtained
curves display a fairly linear behavior with zero intercept, thus denoting a
direct proportionality between the nucleus strain rate and the corresponding
divergence of the velocity field the nucleus is subjected to. In order to
rationalize the physical meaning of the proportionality constant $m$ between
$\nabla\cdot\vec{v}$ and $\dot{a}_{n}$, we note that, under the hypothesis of
weak spatial dependence of the monolayer cell number density $\rho_{c}$, the
continuity equation for the monolayer density can be written as

$\nabla\cdot\vec{v}=-\dot{\rho_{c}}=\dot{a}_{c},$ (6)

where $a_{c}$ is instantaneous cell area divided by its mean value and
$\dot{a}_{c}$ is the cell strain rate. The coefficient $m$ can be therefore
interpreted as the ratio between the nuclear area strain rate and the cell
area strain rate. The reciprocal value $e_{n}=m^{-1}$ is related to the
relative stiffness of the nucleus compared to the one of the entire cell: the
larger is $e_{n}$, the lower is the deformation per unit time of the nucleus
compared to the one of the cell. Estimates of $m$ for the different cell lines
and ages are obtained via a linear fit of the corresponding $\dot{a}_{n}$
versus $\nabla\cdot\vec{v}$ in the neighborhood of the origin. The
corresponding values of $e_{n}$ are reported in figure 7.f. Interestingly
enough, the relative stiffness has a time dependence which is rather different
from the one displayed by both the divergence and the strain rate. While in
all jamming monolayers both $\sigma_{\nabla}$ and $\dot{a}_{0}$ reach a stable
plateau value after about $12$ hours, the relative stiffness typically keep on
decreasing over time, with a roughly constant rate. As it can be appreciated
from figure 7.c, for large values of $\nabla\cdot\vec{v}$ a significant
deviation from linearity is observed, with $\dot{a}_{n}$ seemingly approaching
a saturation value for both negative and positive strain rates. The ratio
between nuclear strain rate and cell strain rate becomes thus lower and lower
as the imposed strain rates increase.

###  3.2 Flocking monolayers

####  3.2.1 In flocking monolayers, nuclei experience stronger deformations
and are stiffer

To further investigate the impact of cell motility on nuclear deformation, we
repeated the experiments using two different cell lines (MCF10A and
MCDF10A.DCIS.com) whose dynamical state is perturbed by inducing the
overexpression of the RAB5A protein, a master regulator of endocythosis
zerial2001rab . Upon RAB5A overexpression, mature, almost completely
kinetically arrested monolayers experience a dramatic reawakening of motility,
characterized by highly coordinated directed migration (flocking) and by the
presence of local cell rearrangments (fluidization). malinverno2017endocytic ;
giavazzi2018flocking ; giavazzi2017giant . For both considered cell lines, we
compare RAB5A overexpressing monolayers with the corresponding controls in
order to investigate the impact of the flocking transition on nuclear
deformation. As detailed in appendix D, PIV analysis confirms the striking
motility phenotype of RAB5A overepressing monolayer previously reported
malinverno2017endocytic ; palamidessi2019unjamming . In figure 9 we show the
average nuclear strain rate (fig. 9.a) and the standard deviation of the
divergence of the velocity field (9.b) for the two cell lines. As it can be
seen, both $\dot{a}_{0}$ and $\sigma_{\nabla}$ are systematically larger in
flocking monolayers, where the relative motion of cells is enhanced (leading
to larger values of divergence) and, as a consequence, nuclei are subjected to
larger stresses and deform more. Despite the significant difference in
$\dot{a}_{0}$ and $\sigma_{\nabla}$ between the two cell lines, a striking
correspondence between the two series of experiments can be observed in the
slopes of the linear regime of the scatter plots reported in figure 9.d and
9.e: as it can be clearly seen in figure 9.c, values of $e_{n}$ are very
similar for the corresponding conditions, with a systematic difference by a
factor $\sim 2$ between the RAB5A overexpressing monolayers and the control
ones. This suggests that the nuclei of RAB5A overexpressing cells deform less,
compared to the control ones, when subjected to the same local density
variation.

Figure 9:  Data relative to control (blue) and RAB5A overexpressing (red)
cells in experiments on flocking monolayers. Results relative to MCF10A (plain
circles) and MCF10A.DCIS.com (void triangles) are presented. a. Nuclear
characteristic strain rates $\dot{a}_{0}$ as a function of time, obtained from
the fit of the corresponding MSS relative to the jamming experiments. b.
Divergence standard deviation as a function of time. c. Reciprocal $e_{n}$ of
the slope of $\dot{a}_{n}$ versus $\nabla\cdot\vec{v}$ scatter plot, obtained
from linear fit of the linear part of the scatter plots. d,e. Scatter plot of
the instantaneous single nucleus strain rate $\dot{a}_{n}$ as a function of
the corresponding divergence. Data relative to MCF10A (d) and MCF10A.DCIS.com
(e) monolayers. Reported points are average over hundreds/thousands of strain
rate-divergence couples made over bins equally spaced along the abscissae
axis. Error bars are evaluated as standard deviations of the mean of the data
contained in the bins. Lines are linear fits with zero intercept made over the
linear part of the scatter plots. In the inset of (d), it is reported a zoom
of the scatter plot over the area where both samples exhibit a linear
behaviour.

The $\dot{a}_{n}$ versus $\nabla\cdot\vec{v}$ curves reported in panels 9.d,e
display some interesting features even beyond the linear regime. As far as
MCF10A cells are concerned, control monolayers display a clear deviation from
the linear regime when $|\dot{a}_{n}|$ is above $\sim 0.05$ h-1, while the
corresponding curve for RAB5A overexpressing cells remains linear for
$|\dot{a}_{n}|<0.1$ h-1. On the other hand, for MCF10A.DCIS.com cells, both
control and RAB5A overexpressing monolayers display a non-linear behaviour at
nuclear strain rates larger than $\sim 0.1$ h-1 and, apparently, the nuclear
strain rate seems to attain a plateau at $|\dot{a}_{n}|\sim 0.2$ h-1. The
occurrence of such a plateau suggests the presence of mechanisms preventing
the nuclei to be deformed above a certain threshold. Moreover, the marked
difference between the two cell lines is even more interesting considering
that MCF10A.DCIS.com is a tumoral cell line deriving from MCF10A
millerMCF10DCIS . Further investigation is required to clarify whether the
expanded linear regime (and the wider dynamic range of $\dot{a}_{n}$)
displayed by MCF10A.DCIS.com control monolayers compared to the corresponding
MCF10A ones could be related to the tumorogenic potential of the cells.

###  3.3 Model

Figure 10:  Schematic representation of a cell deforming under the action of
an in-plane compressive stress. On the right, 1D model of the cell as a series
of two elastic elements with different elastic moduli $E_{n}$ and $E_{cy}$,
corresponding to the nucleus and the cytoplasm, respectively. Corresponding
areas $A_{n}$ and $A_{cy}$ are represented as lengths of the springs in the 1D
model. The total area $A_{c}$ given by the sum of $A_{n}$ and $A_{cy}$.

Our in-plane observations do not allow discriminating whether the observed
fluctuations in nuclear and cellular projected area can be assumed to occur at
constant volume, with each change in the projected area corresponding to a
height variation of opposite sign. In principle, our observations are
compatible with an actual variation of cell and nuclear volume, mediated by
fluid exchange between nucleus and cytoplasm and extracellular medium. Large
volume fluctuations in multicellular aggregates, due to active water transport
in and out of the cell have been indeed reported. The complex interplay
between single-cell volume fluctuations and collective monolayer dynamics has
been investigated via a combination of velocimetry and single-cell
segmentation, similar to one exploited also in the present work
zehnder2015multicellular ; zehnder2015cell . Although the reported timescale
of these volume fluctuations (a few hours) zehnder2016langevin  is typically
much larger than the correlation time of the projected area fluctuations
observed in our experiments (a few minutes), we cannot rule out the
possibility that the observed behavior is actually due to a poroelastic-like
response of the cell darcy1856fontaines . Even without making assumptions on
the details of the underlying process, we can provide a description of the
mechanical properties of cells and nuclei in terms of effective quantities,
building on the experimentally accessible observables described in the
previous sections. Given the limited temporal resolution of our experiments
(2-10 minutes) we do not expect to be able to capture the full rheological
response of the cells. In particular, we can hardly probe the process of
viscoelastic relaxation upon stress application, which typically occurs on the
scale of few tens of seconds khalilgharibi2019stress . We are then most likely
sensitive to the elastic response of cell components. Indeed, the lack of
delay between $\dot{a}_{n}$ and $\nabla\cdot\vec{v}$ observed from their
cross-correlation (see figure 18) is consistent with this hypothesis.

We introduce a simple phenomenological model (sketched in figure 10) for the
mechanical response of the cell to in-plane homogeneous compressive or tensile
stresses, where an elastic nucleus and an elastic cytoplasm, with constant
elastic moduli $E_{n}$ and $E_{cy}$, respectively, are connected in series
Within this simple model, we can derive an explicit relationship between the
ratio $E_{n}/E_{cy}$ of the elastic moduli and relative stiffness $e_{n}$
introduced in section 3.1.4. For two elastic components in series, the total
area deformation $\Delta A_{c}$ is the sum of the deformations $\Delta A_{n}$
and $\Delta A_{cy}$ of the components, while the external stress $\sigma$
applied on the system coincides with the one applied on each element

$\sigma\propto\left(E_{n}a_{n}\right)=\left(E_{cy}a_{cy}\right).$ (7)

By introducing the ratio between the average projected nuclear area and the
average projected cell area
$\beta\doteq\left<A_{n}\right>/\left<A_{c}\right>$, it is possible to write
the total strain associated to the cell projected area $a_{c}=\Delta
A_{c}/A_{c}$ in terms of the nuclear $a_{n}$ and the cytoplasmic strain
$a_{cy}$, as

$a_{c}=\beta a_{n}+\left(1-\beta\right)a_{cy}.$ (8)

Combining this equality with the time derivative of equation 7, we get

$\frac{E_{n}}{E_{cy}}=\frac{e_{n}-\beta}{1-\beta}$ (9)

A time-resolved estimate of $\beta$ can be obtained directly from the images,
by combining the average nuclear area (Fig. 5) with the average cell area,
which can be simply estimated as the total area of the FOV divided by the
number of identified cells. In figure 11.a we report the resulting time
evolution of $\beta$ for the jamming monolayers, evaluated over the same time
intervals considered in panels 7.d-e-f. We note that $\beta$ takes slightly
different values for different samples, and it is approximately time
independent.

Figure 11:  Time evolution of the surface fraction $\beta$ covered by the
nuclei (a.) and of the elastic moduli ratio $E_{n}/E_{cy}$ (b.) for the set of
experiments on jamming monolayers. Data reported as a function of time in the
same time intervals of figures 7.d-f. Error bars evaluated by propagating the
corresponding errors.

Once $\beta$ is known, we can use equation 9 to obtain an estimate of the
ratio between the elastic moduli As it can be seen from figure 11.b, DCIS
moduli ratio uniformly decreases from values above unity to values less than
$1$. Their nuclear modulus is therefore initially larger than the one of
chytoplasm, but the relation reverses as far as cell packing increases.
Similarly, in MCF10AT initially both nucleus and cytoplasm have similar moduli
but with time cytoplasm becomes stiffer than nucleus. Finally, in MCF10CA
$E_{cy}$ is systematically larger than $E_{n}$. Possibility to observe nuclei
either with lower or higher elastic modulus is not surprising as already
observed in literature fischer2020effect  Robustness of the presented results
is granted by the small relative errors, obtained from the standard deviation
of the mean evaluated over $5$ different FOVs for each sample, and by the
similarity of values and trends among difference seeding concentrations of
same samples.

Through the same procedure, we estimate $\beta$ (fig. 12.a) and $E_{n}/E_{cy}$
(fig. 12.b) for the flocking monolayers. While in control samples the elastic
modulus of the nuclei is always close to the one of the cytoplasm, in RAB5A
overexpressing cells it is almost twice as large as the cytoplasmic one, for
both MCF10A and MCF10A.DCIS.com. In principle, these results would be
compatible with both nuclear stiffening and cytoplasmic softening, as a
consequence of RAB5A overexpression. However, AFM measurements of monolayer
rigidity (which are known to probe mainly nuclear stiffening) have already
documented in MCF10A cell line an increase by a factor two in the Young’s
modulus in RAB5A overexpressing confluent cells malinverno2017endocytic .

Figure 12:  Time evolution of the surface fraction $\beta$ covered by the
nuclei (a.) and of the elastic moduli ratio $E_{n}/E_{cy}$ (b.) for the set of
experiments on flocking monolayers. Data reported as a function of time in the
same time intervals of figures 9.a-c. Error bars are evaluated by propagating
the corresponding errors.

##  4 Conclusions

In this work, we have introduced and demonstrated an experimental non-invasive
procedure aimed at extracting information on the relative stiffness of nuclei
compared to the cytoplasm, for epithelial cell monolayers whose state is
monitored during time-lapse microscopy experiments. The procedure involves
segmentation and tracking of the fluorescently-tagged nuclei, combined with
PIV and PT analysis of the nuclear motion.

More specifically, space-resolved and time-resolved studies of the monolayer
dynamics and the nuclear deformation allows the study of their interplay, from
which we find that for small strain rates, a linear correlation holds between
the local density changes, estimated via the divergence of the velocity field,
and the nuclear strain rate. In this linear regime, consistent with an elastic
response,[In this linear elastic regime,] we extract information on the
effective relative stiffness of the nuclei, as we show with experiments
performed with both jamming and flocking monolayers. For jamming monolayers of
three different cell lines, we find nuclear moduli close to the ones of
cytoplasm; we also find that, as the monolayers jam over time and the dynamics
lowers and internal agitation within the monolayer decreases, the nucleus
softens compared to the cytoplasm. In agreement with this observation, a set
of experiments on flocking monolayers obtained by overexpression of the RAB5A
protein in two cell lines shows that the nuclear relative stiffness is higher
than in control samples, despite the fact that nuclei deform more. The
difference in the relative nuclear moduli is significant (about a factor of
two) and is surprisingly similar for the two cell lines.

For large strain rates, we observe nuclear stiffening, with strain rates
eventually approaching a plateau as a function of the divergence of the
velocity field. We interpret this result as a mechanoprotective response of
the nuclei to oppose large stresses that arise from the increased internal
agitation of the nuclei in a highly dynamic monolayer. Interestingly, our
results are compatible with the view that the value of the nuclear elastic
modulus for small strains is mainly determined by chromatin, whereas for
larger strains the strain-stiffening contribution of lamins dominates
stephens2017chromatin ; stephens2018separate .

The possibility of non-invasive monitoring of the relationship between
mechanical stimuli, nuclear deformations and nuclear relative stiffness in
cell monolayers may have immediate biological applications. An early version
of the methodology that we describe here was recently used in Ref.
frittoli2021tissue  to assess in real time the nuclear mechanical stress and
response and the associated role in causing a long-term transcriptional-
dependent phenotype as a consequence of the short-term adaptive response to
stress. We thus believe that our methodology represents a useful addition to
the portfolio of non-invasive tools to characterize the mechanobiological
response of tissues.

## Acknowledgements

This work has been supported by: the Associazione Italiana per la Ricerca sul
Cancro (AIRC) to F.G. and S.V. (MFAG#22083), and to G.S. (IG#18621 and 5XMille
#22759); the Italian Ministry of University and Scientific Research (MIUR) to
G.S. (PRIN: Progetti di Ricerca di Rilevante Interesse Nazionale – Bando
2017#2017HWTP2K).

## Author contribution statement

F.G. conceived research. R.C., G.S. and F.G. designed research. A.P. and E.F.
prepared the samples and performed the experiments. S.V. developed the image
processing pipeline, wrote the codes and analyzed all the data. S.V., R.C. and
F.G. wrote the manuscript. All Authors contributed discussing the results and
the manuscript.

##  Appendix A . Mean square displacement analysis

In this Appendix, we discuss the mean square displacements (MSDs) associated
to the trajectories obtained from nuclear tracking, as described in section
2.2. For each set of experiments, the MSD is calculated within each of the
five time windows considered in the main text. In figure 13.a-b-c, we display
averages over all five time windows as a function of the lag time $\tau$.
There, continuous and dashed black lines indicates the slope corresponding to
a ballistic scaling (quadratic in $\tau$), and a diffusive one (linear in
$\tau$), respectively. As it can be appreciated also visually, in all
experiments a crossover from a ballistic-like scaling at short time delays to
a diffusive-like scaling at larger time scales is observed. The characteristic
time $\tau_{p}$ marking the transition between the two regimes, which
corresponds to the persistence time during which the velocity of a given cell
is constant, can be estimated by fitting to the data the function

$MSD\left(\tau\right)=\sigma_{MSD}+2\tau_{p}v_{b}^{2}\left[\tau+\tau_{p}\left(e^{-\tau/\tau_{p}}-1\right)\right]$
(10)

where $v_{b}$ is the ballistic velocity and $\sigma_{MSD}$ is an offset to
account for experimental noise uhlenbeck1930theory . The persistence times
recovered from the fit are reported in figure 13.d-e-f. As it can be seen,
they are systematically larger than the time step between two acquired frames
in the corresponding experiments, represented by horizontal dotted lines. This
observation implies that, on the timescale $\delta t$ corresponding to the
delay between two consecutive images, cells are moving with good approximation
with constant velocity. This confirms the validity of the assumption made in
section 2.2.3 that the instantaneous nuclear velocity can be reliably
estimated as the average velocity between two consecutive frames.

Figure 13:  First row: nuclear mean square displacements (MSD) obtained from
experiments on jamming monolayers (a.), on MCF10A CTRL and RAB5A
overexpressing flocking monolayers (b.), and on DCIS.com CTRL and RAB5A
overexpressing flocking monolayers (c.). In the legend of panel (a.) HC and LC
correspond to high and low seeding density, respectively. The reported MSD are
averaged over the entire duration of the experiment. Corresponding best
fitting curves with a ballistic-to-diffusive transition model (see Eq. 10) are
reported as lines of the same colors. Black continuous and dashed lines
represent guidelines for ballistic and diffusive dynamics, respectively. Each
curve is plotted after subtraction of a constant offset $\sigma_{MSD}$,
obtained from the fitting procedure (see Eq. 10). In the second row, we
display the persistence times obtained for the jamming experiments (d.), the
MCF10A flocking experiment (e.), and the DCIS.com flocking experiment (f.),
within each one of the five considered time windows. In each panel, the dotted
line correspond to the value $\delta t$ of the delay time between two
consecutive acquired frames.

##  Appendix B . Comparison of PIV results from fluorescent and phase contrast
time lapses

In this Appendix, we compare the dynamics parameters obtained from PT and from
PIV made over fluorescent and phase contrast time lapses. We then deepen the
discrepancies between the velocities and velocity divergences recovered from
PT and from PIV analysis of fluorescent and phase contrast time lapses, in
order to motivate the choice of the fluorescent PIV analysis to evaluate the
divergence of the velocity field as an indicator of the local density
fluctuations within the monolayer.

In figure 14 are reported the velocity distributions along the two main axis
for a flocking MCF10A RAB5A overexpressing monolayer. The reported dynamics
refers to the $8$-$9$ hours time interval, corresponding to the cell motility
peak malinverno2017endocytic . Blue, orange and green lines refer to PIV
computed over phase contrast and fluorescent time lapses and to PT,
respectively. The recovered distributions are very similar to each other, thus
validating the obtained dynamics.

Figure 14:  Velocity components $v_{x}$ (a.) and $v_{y}$ (b.) distribution
within an MCF10A RAB5A overexpressing FOVs. Distributions are obtained from
PIV computed on phase contrast images (blue), from PIV computed over
fluorescent images (orange) and from PT (green). Distribution are recovered in
the time interval $8$-$9$ h.

Mean values obtained from PT and PIV are also very close, as shown in figure
15 where are reported for the six jamming experiments the time evolution of
$v_{RMS}$ as recovered from PIV made over fluorescent time lapses (fig. 15.a)
and from PT (fig. 15.b).

Figure 15:  Time evolution of RMS velocity ${v}_{RMS}$ obtained from PIV on
fluorescent time lapses (a.) and PT (b.). Shadowed error bars are avaluated as
the standard deviation of the mean evaluated for each sample over the 5
different FOVs. Reported data refer to jamming monolayers: DCIS HC (red thick
line) and LC (dark red thin line), for MCF10AT HC (blue thick line) and LC
(dark blue thin line) and for MCF10CA HC (green thick line) and LC (dark green
thin line).

Given the correspondence between the results obtained from PT and from PIV, it
is possible to choose which of the two methods to use to compare velocity
divergence with the nuclear deformations shown in section 3.1.4. A first
possibility would be to start from the single nuclei velocities
$\vec{v}_{j}\left(t_{i}\right)$ to build up frame by frame a velocity field
through a proper interpolation of the known velocities in the plane. This
choice has been discarded in order to avoid spurious correlations deriving
from the common origin from the localisation and segmentation of nuclei of the
correlating data (divergence of the velocity field and time derivative of the
area). We therefore chose to use a divergence obtained from PIV to ensure the
authenticity of any correlation. Additionally, even a 10% missed nuclei in PT
while allowing a fair enough agreement with PIV results in the velocity
distribution (Fig. 14), at the same time prevents a proper interpolation of
the velocity field. Test made revealed a discrepancy even by a factor 2 in the
width of the $\nabla\cdot\vec{v}$ distributions between PT and PIV, either
from fluorescent and phase contrast analysis.

A deeper analysis revealed that even between the velocity fields recovered
applying PIV analysis on fluorescent and phase contrast images have some
discrepancies. We evaluate the normalized space cross-correlation between the
velocity fields moduli $v_{Fluo}\left(\vec{x}\right)$ and
$v_{PhC}\left(\vec{x}\right)$ obtained from the PIV analysis of the two data
images respectively:

$K_{v_{Fluo},v_{PhC}}\left(\vec{r}\right)=\frac{\left<v_{Fluo}\left(\vec{x}\right)\cdot
v_{PhC}\left(\vec{x}+\vec{r}\right)\right>_{\vec{x}}}{\sqrt{\left<v_{Fluo}\left(\vec{x}\right)\right>_{\vec{x}}^{2}\left<v_{PhC}\left(\vec{x}\right)\right>_{\vec{x}}^{2}}}$
(11)

where $\left<...\right>_{\vec{x}}$ stands for the spatial averages over all
the positions where the velocity fields are evaluated. In figure 16.b it is
reported the azimuthal average of $K_{v_{Fluo},v_{PhC}}$ together with the
space self-correlations of $v_{Fluo}\left(\vec{x}\right)$ and
$v_{PhC}\left(\vec{x}\right)$. The spatial scale of decorrelation is similar
in the three cases, but the cross-correlation value in $r=0$ is lower by $20$%
than the self-correlations, thus testifying a weak but existent difference in
the velocity fields.

Figure 16:  a. Relative difference between the modulus of the velocity field
obtained from PIV from fluorescent and phase contrast images. b. Azimuthal
average of the space correlations $K_{v_{Fluo},v_{Fluo}}$ (orange),
$K_{v_{PhC},v_{PhC}}$ (blue) and $K_{v_{Fluo},v_{PhC}}$ (black). c. Relative
difference between the divergence of the velocity field obtained from
fluorescent and phase contrast images. d. Azimuthal average of the space
correlations $K_{\nabla_{Fluo},\nabla_{Fluo}}$ (orange),
$K_{\nabla_{PhC},\nabla_{PhC}}$ (blue) and $K_{\nabla_{Fluo},\nabla_{PhC}}$
(black).

Similarly, in figure 16.d we report the azimuthal average of the cross-
correlation relative to the divergences of the velocity fields defined as:

$K_{\nabla_{Fluo},\nabla_{PhC}}\left(\vec{r}\right)]=\frac{\left<\left[\nabla\cdot\vec{v}_{Fluo}\left(\vec{x}\right)\right]\left[\nabla\cdot\vec{v}_{PhC}\left(\vec{x}+\vec{r}\right)\right>_{\vec{x}}\right]}{\sqrt{\left<\nabla\cdot\vec{v}_{Fluo}\left(\vec{x}\right)\right>_{\vec{x}}^{2}\left<\nabla\cdot\vec{v}_{PhC}\left(\vec{x}\right)\right>_{\vec{x}}^{2}}}$
(12)

As in the previous case, self-correlations of
$\nabla\cdot\vec{v}_{Fluo}\left(\vec{x}\right)$ and
$\nabla\cdot\vec{v}_{PhC}\left(\vec{x}\right)$ are reported in the same plot.
For the divergences, the cross-correlation in $r=0$ is much lower than the
correspondent self-correlations, thus pointing out much larger differences
between the divergence fields evaluated with the two methods. Discrepancy
between the fields of $\nabla\cdot\vec{v}$ much more pronounced than the ones
in the corresponding $v$ can also be appreciated comparing figures 16.a and
16.c, where the relative difference between the fluorescent and the phase
contrast analysis are represented for the modulus of the velocity field and
the divergence respectively. While for $v$ the discrepancies only emerges in
the regions where there is a low density of nuclei (dark blue regions in fig.
16.a), very different appear to be $\nabla\cdot\vec{v}_{fluo}$ from
$\nabla\cdot\vec{v}_{phc}$. We note that such differences are only present
when spatial resolution is considered. Distributions of the divergence, frame
by frame, are very close in the two cases.

Given the differences, a choice has to be made between the two methods for the
comparison with the nuclear strain rates. To this aim, we report in figure 17
scatter plots of the single nuclei velocity components $v_{x,j}$ (fig. 17.a,c)
and $v_{y,j}$ (fig. 17.b,d) obtained from PT as a function of the
corresponding components of the cubic interpolation of the velocity field
$\vec{v}_{fluo}$ (fig. 17.a,b) and $\vec{v}_{phc}$ (fig. 17.c,d) in the center
of mass of the same nuclei. Looking at dispersions, one can find out that
correspondence is better between PT and PIV analysis of fluorescent time
lapses. This can be more easily sensed by looking at the width of the
distributions along their minor dimension. This can be made by compuing the
minor eigenvalue of the covariance matrix between the velocity component
obtained from PT and from PIV. The resulting widths are $\sigma_{x}=0.012$ and
$\sigma_{y}=0.012$ for the fluorescent time lapses case and $\sigma_{x}=0.016$
and $\sigma_{y}=0.017$ for the phase contrast one.

Larger discrepancies from PT in the latter case can be partially rationalized
considering that in phase contrast inner dynamics of the cells is also
visible. Moreover, visual inspection of the videos revealed that several
bright spots appear in phase contrast localized at the edges between adjacent
cells. They move along cell edges faster than the cells themselves, thus
affecting the velocity obtained from PIV.

In conclusion, because of this larger discrepancy between PT and phase
contrast, we choose to compare area strain rate with the divergence obtained
from PIV analysis of fluorescent time lapses.

Figure 17:  a-b. Correlation between the velocities components recovered from
PT and the ones obtained interpolating at the nuclei center coordinates the
velocity field obtained from PIV on fluorescent time lapses. The minor widths
of the distributions, evaluated as the lower eigenvalue of the covariance
matrix, are $\sigma_{x}=0.012$ and $\sigma_{y}=0.012$ for $v_{x}$ and $v_{y}$
respectively. c-d. Correlation between the velocities components recovered
from PT and the ones obtained interpolating at the nuclei center coordinates
the velocity field obtained from PIV on phase contrast time lapses. The minor
widths of the distributions are now $\sigma_{x}=0.016$ and $\sigma_{y}=0.017$
for $v_{x}$ and $v_{y}$ respectively. Scatter plots are made considering the
frames in the time interval $1$-$8$ hours relative to one of the field of view
of DCIS HC jamming experiments.

##  Appendix C . Self and Cross correlations between local divergence and
nuclear strain rate

In the main text, we report the normalized cross-correlation between nuclear
strain rate $\dot{a}_{j}\left(t\right)$ and the corresponding divergence of
the velocity field $\nabla\cdot\vec{v}_{j}$. Normalized cross-correlation for
generic variables $q_{j}\left(t\right)$ and $p_{j}\left(t\right)$ is defined
as

$K_{p,q}\left(\tau\right)=\frac{\left<\left<p_{j}\left(t+\tau\right)\cdot
q_{j}\left(t\right)\right>_{t}\right>_{j}}{\sqrt{\left<\left<p_{j}^{2}\left(t\right)\right>_{t}\right>_{j}^{2}\cdot\left<\left<q_{j}^{2}\left(t\right)\right>_{t}\right>_{j}^{2}}}$
(13)

Figure 18:  a. Normalized autocorrelation of the instantaneous nuclear strain
rate $\dot{a}_{n}$. Data relative to jamming experiment DCIS HC in the time
interval $0$-$24$h. b. Normalized autocorrelation of the divergence of the
velocity field. Data relative to jamming experiment DCIS HC in the time
interval $0$-$24$h. c. Normalized cross-correlation between the instantaneous
nuclear strain rate $\dot{a}_{n}$ and the divergence of the velocity field
evaluated in the center of mass of the same nucleus. Data relative to jamming
experiment DCIS HC in the time interval $0$-$24$h. d. Same data of c. but
averaging over positive and negative $\tau$. Black dashed line represent the
best exponential fit detailed in the text. e. Decay rate of
$K_{\dot{a},\nabla}$ as a function of the considered time interval for the
different experiments on jamming monolayers. Red, blue and green refer to
DCIS, MCF10AT and MCF10CA monolayers, respectively. Thin and thick lines
correspond to monolayers seeded at low and high concentration, respectively.

In figure 18 are reported the nuclear strain rate (fig. 18.a) and divergence
(fig. 18.b) self-correlations $K_{\dot{a},\dot{a}}$ and $K_{\nabla,\nabla}$
together with the cross-correlation $K_{\dot{a},\nabla}$ (fig. 18.c). In
figure 18.d the average of positive and negative $\tau$ of fig. 18.c is made.
A strong decorrelation at the first delay can be observed in
$K_{\dot{a},\dot{a}}$. This probably is a consequence of the high level of
noise present on the determination of $\dot{a}$, and not of a real lack of
time correlation of the strain rate. White noise origin of the strong time
decorrelation is confirmed by the fact that $\dot{a}$ and $\nabla\vec{v}$
instead correlates in time (18.c), with a characteristic time comparable with
the self correlation characteristic time $\tau_{K}$ of
$\nabla\vec{v}\left(t,\vec{x}\right)$ (Fig. 18.b). The fact that an apparently
delta-correlated signal - as the nuclear strain rate seems to be - correlates
in time with the divergence might seem weird. As we show in the following, it
is explained assuming that also nuclear strain rate has a finite time
correlation hidden behind a strong delta correlated white noise.

For a generic quantity $q\left(t\right)$ whose real measurement
$\tilde{q}\left(t\right)=q\left(t\right)+\sigma_{q}$ is affected by the delta
correlated white noise $\sigma_{q}$, the non normalized self correlation is:

$\left<\tilde{p}\left(t\right)\cdot\tilde{p}\left(t+\tau\right)\right>_{t}=\left<p\left(t\right)\cdot
p\left(t+\tau\right)\right>_{t}+\sigma_{p}^{2}\delta\left(\tau\right)$ (14)

For the nuclear strain rate, it is possible to recover an indicative value of
$\sigma_{p}$ from the $MSS$ as the square root of the offset $\sigma_{w}$. For
the reported DCIS HC monolayer, $\sqrt{\sigma_{w}}=0.07$. If we compare this
value with the square root of the non normalized self-correlation evaluated in
$\tau=0$, equal to $0.08$, we perceive how small is the hidden genuine
correlation compared to white noise.

If the large white noise impedes the observation of time correlation of
nuclear strain rate, it does not hinder information on the correlation with
independently measured variables like divergence. Indeed, in cross correlation
between noisy data $\tilde{p}$ and $\tilde{q}$ with uncorrelated errors
$\sigma_{p}$ and $\sigma_{q}$ it results:

$\left<\tilde{p}\left(t\right)\cdot\tilde{q}\left(t+\tau\right)\right>_{t}=\left<p\left(t\right)\cdot
q\left(t+\tau\right)\right>_{t}$ (15)

as other terms vanish because of the delta correlation of noises $\sigma_{p}$
and $\sigma_{q}$.

Oscillations at negative values of correlation at finite $\tau$ in
$K_{\dot{a},\dot{a}}$ and $K_{\nabla,\nabla}$, and consequently in
$K_{\dot{a},\nabla}$, can be rationalized considering that, as
kubo1966fluctuation :

$\lim_{\tau\to\infty}\frac{\left<a\left(\tau|t\right)^{2}\right>_{t}}{2\tau}=\int_{0}^{\infty}\left<\dot{a}\left(t\right)\dot{a}\left(t+\tau\right)\right>\,d\tau,$
(16)

in order to reach a plateau of $MSS$ for $\tau\rightarrow\infty$ the integral
of the time correlations must be zero and thus correlations must assume
negative values.

From the cross-correlation, we recovered the characteristic time $\tau_{k}$
through an exponential fit excluding point in $\tau=0$. Inclusion of $\tau=0$
results in a failure of both exponential and stretched exponential fits.
Apparently, a first very fast decorrelation arrived at time scales not
accessible in our experiments. An example of the fit is superposed to data in
figure 18.d and the resulting decorrelation rates $1/\tau_{k}$ are reported in
figure 18.e as a function of time for the experiments on jamming monolayers.

##  Appendix D . Dynamics of flocking monolayers

Figure 19: Time evolution of RMS velocity ${v}_{RMS}$ (a.) and of the modulus
of the center of mass velocity (b.) obtained from PIV on fluorescent time
lapses. Shadowed error bars are evaluated as the standard deviation of the
mean evaluated for each sample over different FOVs. Reported data refers to
flocking monolayers (blue for the control monolayers, red for the RAB5A
overexpressing ones): MCF10A (continuous line) and MCF10A.DCIS.com (dashed
line).

In figure 19 are reported the time evolution of $v_{RMS}$ (fig. 19.a) and of
the modulus of $\vec{v}_{CM}$ (fig. 19.b). As expected malinverno2017endocytic
; giavazzi2018flocking , both MCF10A and MCF10A.DCIS.com cell lines exhibit a
reawakening of motility as a consequence of RAB5A overexpression, which is
characterized by the emergence of a collective migration as a consequence of
mutual velocity aligment of neighboring cells - testify by the motility peak
in $\vec{v}_{CM}$ \- and by the acquisition of a fluid-like dynamical state in
the center of mass reference system, as evident from the dramatic difference
in $v_{RMS}$ between RAB5A overexpressing monolayers and control ones.

## References

  * (1) J.A. Park, J.H. Kim, D. Bi, J.A. Mitchel, N.T. Qazvini, K. Tantisira, C.Y. Park, M. McGill, S.H. Kim, B. Gweon et al., Nature materials 14, 1040 (2015)
  * (2) M. Sadati, N.T. Qazvini, R. Krishnan, C.Y. Park, J.J. Fredberg, Differentiation 86, 121 (2013)
  * (3) E. Blauth, H. Kubitschke, P. Gottheil, S. Grosser, J.A. Käs, Frontiers in Physics pp. NA–NA (2021)
  * (4) C. Malinverno, S. Corallino, F. Giavazzi, M. Bergert, Q. Li, M. Leoni, A. Disanza, E. Frittoli, A. Oldani, E. Martini et al., Nature materials 16, 587 (2017)
  * (5) A. Palamidessi, C. Malinverno, E. Frittoli, S. Corallino, E. Barbieri, S. Sigismund, G.V. Beznoussenko, E. Martini, M. Garre, I. Ferrara et al., Nature materials 18, 1252 (2019)
  * (6) J.H. Kim, A.F. Pegoraro, A. Das, S.A. Koehler, S.A. Ujwary, B. Lan, J.A. Mitchel, L. Atia, S. He, K. Wang et al., Biochemical and biophysical research communications 521, 706 (2020)
  * (7) L. Oswald, S. Grosser, D.M. Smith, J.A. Käs, Journal of physics D: Applied physics 50, 483001 (2017)
  * (8) O. Ilina, P.G. Gritsenko, S. Syga, J. Lippoldt, C.A. La Porta, O. Chepizhko, S. Grosser, M. Vullings, G.J. Bakker, J. Starruß et al., Nature cell biology 22, 1103 (2020)
  * (9) C. Alibert, B. Goud, J.B. Manneville, Biology of the Cell 109, 167 (2017)
  * (10) V. Swaminathan, K. Mythreye, E.T. O’Brien, A. Berchuck, G.C. Blobe, R. Superfine, Cancer research 71, 5075 (2011)
  * (11) N. Wang, J.D. Tytell, D.E. Ingber, Nature reviews Molecular cell biology 10, 75 (2009)
  * (12) C.M. Denais, R.M. Gilbert, P. Isermann, A.L. McGregor, M. Te Lindert, B. Weigelin, P.M. Davidson, P. Friedl, K. Wolf, J. Lammerding, Science 352, 353 (2016)
  * (13) M. Raab, M. Gentili, H. de Belly, H.R. Thiam, P. Vargas, A.J. Jimenez, F. Lautenschlaeger, R. Voituriez, A.M. Lennon-Duménil, N. Manel et al., Science 352, 359 (2016)
  * (14) Y. Xia, I.L. Ivanovska, K. Zhu, L. Smith, J. Irianto, C.R. Pfeifer, C.M. Alvey, J. Ji, D. Liu, S. Cho et al., Journal of Cell Biology 217, 3796 (2018)
  * (15) C.R. Pfeifer, Y. Xia, K. Zhu, D. Liu, J. Irianto, V.M.M. García, L.M.S. Millán, B. Niese, S. Harding, D. Deviri et al., Molecular biology of the cell 29, 1948 (2018)
  * (16) S. Grosser, J. Lippoldt, L. Oswald, M. Merkel, D.M. Sussman, F. Renner, P. Gottheil, E.W. Morawetz, T. Fuhs, X. Xie et al., Physical Review X 11, 011033 (2021)
  * (17) T. Harada, J. Swift, J. Irianto, J.W. Shin, K.R. Spinler, A. Athirasala, R. Diegmiller, P.D.P. Dingal, I.L. Ivanovska, D.E. Discher, Journal of Cell Biology 204, 669 (2014)
  * (18) M.T. Parreira, K. Lavrenyuk, J.M. Sanches, K.N. Dahl, Cytoskeleton (2021)
  * (19) T. Fischer, A. Hayn, C.T. Mierke, Frontiers in cell and developmental biology 8, 393 (2020)
  * (20) H. Liu, J. Wen, Y. Xiao, J. Liu, S. Hopyan, M. Radisic, C.A. Simmons, Y. Sun, ACS nano 8, 3821 (2014)
  * (21) B. Fabry, G.N. Maksym, J.P. Butler, M. Glogauer, D. Navajas, N.A. Taback, E.J. Millet, J.J. Fredberg, Physical Review E 68, 041914 (2003)
  * (22) F. Guilak, J.R. Tedrow, R. Burgkart, Biochemical and biophysical research communications 269, 781 (2000)
  * (23) E.M. Darling, D. Di Carlo, Annual review of biomedical engineering 17, 35 (2015)
  * (24) A.R. Harris, L. Peter, J. Bellis, B. Baum, A.J. Kabla, G.T. Charras, Proceedings of the National Academy of Sciences 109, 16449 (2012)
  * (25) M. Radmacher, R. Tillamnn, M. Fritz, H. Gaub, Science 257, 1900 (1992)
  * (26) J. Mitchison, M. Swann, J. exp. Biol 31, 443 (1954)
  * (27) A.C. Rowat, L.J. Foster, M.M. Nielsen, M. Weiss, J.H. Ipsen, Journal of The Royal Society Interface 2, 63 (2005)
  * (28) I. Titushkin, M. Cho, Biophysical journal 90, 2582 (2006)
  * (29) C. Lim, E. Zhou, S. Quek, Journal of biomechanics 39, 195 (2006)
  * (30) C. Wählby, J. Lindblad, M. Vondrus, E. Bengtsson, L. Björkesten, Analytical Cellular Pathology 24, 101 (2002)
  * (31) F. Xing, L. Yang, IEEE reviews in biomedical engineering 9, 234 (2016)
  * (32) J.S. Lim, _Two-dimensional signal and image processing_ (Englewood Cliffs, 1990), p.584, equations 9.44, 9.45, and 9.46
  * (33) C.R. Maurer, R. Qi, V. Raghavan, IEEE Transactions on Pattern Analysis and Machine Intelligence 25, 265 (2003)
  * (34) A. Trujillo-Pino, K. Krissian, M. Alemán-Flores, D. Santana-Cedrés, Image and Vision Computing 31, 72 (2013)
  * (35) I.T. Jolliffe, _Principal components analysis_ (Springer, 1986)
  * (36) D. Blair, E. Dufresne, Particle-tracking code available at http://physics. georgetown. edu/matlab (2008)
  * (37) J.C. Crocker, D.G. Grier, Journal of colloid and interface science 179, 298 (1996)
  * (38) W. Thielicke, R. Sonntag, Journal of Open Research Software 9 (2021)
  * (39) T.E. Angelini, E. Hannezo, X. Trepat, M. Marquez, J.J. Fredberg, D.A. Weitz, Proceedings of the National Academy of Sciences 108, 4714 (2011)
  * (40) S. Garcia, E. Hannezo, J. Elgeti, J.F. Joanny, P. Silberzan, N.S. Gov, Proceedings of the National Academy of Sciences 112, 15314 (2015)
  * (41) R. Cerbino, S. Villa, A. Palamidessi, E. Frittoli, G. Scita, F. Giavazzi, Soft Matter 17, 3550 (2021)
  * (42) K. Haase, J.K. Macadangdang, C.H. Edrington, C.M. Cuerrier, S. Hadjiantoniou, J.L. Harden, I.S. Skerjanc, A.E. Pelling, Scientific reports 6, 1 (2016)
  * (43) M. Zerial, H. McBride, Nature reviews Molecular cell biology 2, 107 (2001)
  * (44) F. Giavazzi, M. Paoluzzi, M. Macchi, D. Bi, G. Scita, M.L. Manning, R. Cerbino, M.C. Marchetti, Soft matter 14, 3471 (2018)
  * (45) F. Giavazzi, C. Malinverno, S. Corallino, F. Ginelli, G. Scita, R. Cerbino, Journal of Physics D: Applied Physics 50, 384003 (2017)
  * (46) F.R. Miller, S.J. Santner, L. Tait, P.J. Dawson, JNCI: Journal of the National Cancer Institute 92, 1185a (2000)
  * (47) S.M. Zehnder, M.K. Wiatt, J.M. Uruena, A.C. Dunn, W.G. Sawyer, T.E. Angelini, Physical Review E 92, 032729 (2015)
  * (48) S.M. Zehnder, M. Suaris, M.M. Bellaire, T.E. Angelini, Biophysical journal 108, 247 (2015)
  * (49) S.M. Zehnder, F.M. Zegers, T.E. Angelini, Journal of biomechanics 49, 1286 (2016)
  * (50) H. Darcy, _Les fontaines publiques de la ville de Dijon: exposition et application…_ (Victor Dalmont, 1856)
  * (51) N. Khalilgharibi, J. Fouchard, N. Asadipour, R. Barrientos, M. Duda, A. Bonfanti, A. Yonis, A. Harris, P. Mosaffa, Y. Fujita et al., Nature physics 15, 839 (2019)
  * (52) A.D. Stephens, E.J. Banigan, S.A. Adam, R.D. Goldman, J.F. Marko, Molecular biology of the cell 28, 1984 (2017)
  * (53) A.D. Stephens, E.J. Banigan, J.F. Marko, Nucleus 9, 119 (2018)
  * (54) E. Frittoli, A. Palamidessi, F. Iannelli, F. Zanardi, S. Villa, L. Barzaghi, H. Ando, V. Cancila, G. Beznuskenko, G. Della Chiara et al., Under Review (2021)
  * (55) G.E. Uhlenbeck, L.S. Ornstein, Physical review 36, 823 (1930)
  * (56) R. Kubo, Reports on progress in physics 29, 255 (1966)